Problem is related with my other question

So, generate exponentially correlated row $π‘₯_{0,0},π‘₯_{1,0}...π‘₯_{𝑛,0}$ and starting from π‘₯0,0 there will be a correlated column $π‘₯_{0,1},π‘₯_{0,2}...π‘₯_{0,n}$.

During generation we used a correlation $𝜌=𝑒π‘₯𝑝(βˆ’Ξ”πœ/πœ‰)$. And $π‘₯_{𝑗,0}=\sqrt{1βˆ’πœŒ^2}Rand[NormalDistribution[πœ‡,𝜎]]+𝜌 π‘₯_{π‘—βˆ’1,0}$

So...During 2 dimensions it should be $𝜌=𝑒π‘₯𝑝(βˆ’Ξ”πœπ‘₯/πœ‰π‘₯βˆ’Ξ”πœπ‘¦/πœ‰π‘¦)$. And in isotropic space $𝜌=𝑒π‘₯𝑝(βˆ’2Ξ”πœ/πœ‰)$.

But how to generate $π‘₯_{1,1}$ correctly with this $𝜌$? This one with square mean looks absolutely wrong: $π‘₯_{1,1}=\sqrt{1βˆ’πœŒ2}Rand[NormalDistribution[πœ‡,𝜎]]+𝜌\sqrt{π‘₯_{1,0} π‘₯_{0,1}}$

Tried to generate with arithmetic mean: $π‘₯_{1,1}=\sqrt{1βˆ’πœŒ2}Rand[NormalDistribution[πœ‡,𝜎]]+𝜌(π‘₯_{1,0}+ π‘₯_{0,1})/2$ But it look not correlated in Y direction... Any way to make it better?

my code:

ΞΌ = 0; Οƒ = 1; ΞΎ = 10; Ο„ = 1; l = {{0, 0, 
   RandomVariate[NormalDistribution[ΞΌ, Οƒ], 
     1][[1]]}}; size = 100; ρ = E^(-(
   Abs[Ο„]/ΞΎ)); ρ2 = E^(-2 Abs[Ο„]/ΞΎ);
For[j = 1, j < size, j++,
  Subscript[x, 0] = {l[[1, 3]]}; 
   j] = ρ2 Subscript[x, j - 1][[1]] + 
    Sqrt[1 - ρ2^2]
      RandomVariate[NormalDistribution[ΞΌ, Οƒ], 1]; 
  AppendTo[l, {j, 0, Subscript[x, j][[1]]}]];

For[i = 1, i < size, i++,
   0] = ρ2 l[[1 + size (i - 1), 3]] + 
    Sqrt[1 - ρ2^2]
      RandomVariate[NormalDistribution[ΞΌ, Οƒ], 1]; 
  AppendTo[l, {0, i, Subscript[x, 0][[1]]}];
  For[j = 1, j < size, j++,
    j] = ρ2 (Subscript[x, j - 1][[1]] + l[[size i + j, 3]])/2 + 
     Sqrt[1 - ρ2^2]
       RandomVariate[NormalDistribution[ΞΌ, Οƒ], 1]; 
   AppendTo[l, {j, i, Subscript[x, j][[1]]}]]];
  • 2
    $\begingroup$ The standard way to generate correlated noise, is to: 1) Compute the autocorrelation function; 2) Fourier Transform this to get the Power spectrum; 3) Generate independent random numbers; 4) scale by the sqrt of the Power spectrum; 5) Fourier transform; 6) Verify that you have the autocorrelation function you wanted. (This will work in any number of dimensions) $\endgroup$ – mikado Mar 9 at 22:24
  • 1
    $\begingroup$ @mikado: Not sure where I heard this many years ago: "We must love standards because we have so many of them." For me the standard way for this problem is to generate random samples directly from a multivariate normal with the desired covariance structure. $\endgroup$ – JimB Mar 9 at 23:03
  • $\begingroup$ @JimB if I want to generate (say) a 1000x1000 image of correlated noise, my method involves generating 10^6 random samples, weighting each one with a simple function and 2D FFT (quick). Yours involves generating a 10^6 x 10^6 covariance matrix. Mathematica would struggle to generate the resulting random vector. $\endgroup$ – mikado Mar 9 at 23:37
  • $\begingroup$ @mikado Yep, there would be a bit of a struggle. My weak/feeble defense is that I didn't see the OP state any size requirements. $\endgroup$ – JimB Mar 9 at 23:42

As noted by @mikado generating the random sample directly is more than a bit sluggish. If you have something on the order of 10,000 points (a 100 x 100 array), then the following code is simple and straightforward but takes about 15 minutes:

n = 100;
ρ = 0.1;
x = Flatten[Table[{i, j}, {i, n}, {j, n}], 1];
Σ = Table[If[i == j, 1, Exp[-ρ Norm[x[[i]] - x[[j]]]]], {i, Length[x]}, {j, Length[x]}];
z = Transpose[{x[[All, 1]], x[[All, 2]], 
  Flatten[RandomVariate[MultinormalDistribution[ConstantArray[0, Length[x]], Ξ£], 1]]}];]

Contour plot Scatter of sample points


Your Answer

By clicking β€œPost Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.