Autocorrelated random field generation

Question

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]]}; 
  Subscript[x, 
   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++,
  Subscript[x, 
   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++,
   Subscript[x, 
    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]]}]]];
ListContourPlot[l]
ListPlot3D[l]

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) — mikado, Mar 09 '19 at 22:24
@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. — JimB, Mar 09 '19 at 23:03
@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. — mikado, Mar 09 '19 at 23:37
@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. — JimB, Mar 09 '19 at 23:42

score 3 · Answer 1 · answered Mar 10 '19 at 01:43

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]]}];]
ListContourPlot[z]
ListPointPlot3D[z]

Autocorrelated random field generation

1 Answers1

Linked