I need to generate two correlated gaussian random fields. As far as I know, this question and this other provides the means to generate a single autocorrelated process. However, I am clueless about how to include a second random field, correlated with the first. The only way of doing this that I have found is on this blog and isn't a implemented in Mathematica, but in R. I am looking for a Mathematica code which does the same thing as the second blog or presents an alternative algorithm that can generate two processes as the answer to the first or second links but correlated as per the third?

  • $\begingroup$ What exactly do you mean by the two fields being correlated? For example, if you generate one random field, you can then generate a second field just by adding noise to the first field. Or is the correlation structure you seek more complicated than that? $\endgroup$
    – JimB
    Aug 5, 2021 at 19:56

1 Answer 1


I'm not sure what you mean by "correlated" in "correlated gaussian random fields" but you can generate two independent fields and then compare the first with the sum of the first and second. The individual points are then correlated.

Stealing from the excellent answer you posted, one can construct two independent fields:

GaussianRandomField[size : (_Integer?Positive) : 256, 
  dim : (_Integer?Positive) : 2, Pk_ : Function[k, k^-3]] := 
  Module[{Pkn, fftIndgen, noise, amplitude, s2}, 
  Pkn = Compile[{{vec, _Real, 1}}, 
    With[{nrm = Norm[vec]}, If[nrm == 0, 0, Sqrt[Pk[nrm]]]], 
    CompilationOptions -> {"InlineExternalDefinitions" -> True}];
  s2 = Quotient[size, 2];
  fftIndgen = ArrayPad[Range[0, s2], {0, s2 - 1}, "ReflectedNegation"];
  noise = Fourier[RandomVariate[NormalDistribution[], ConstantArray[size, dim]]];
  amplitude = Outer[Pkn[{##}] &, Sequence @@ ConstantArray[N@fftIndgen, dim]];

BlockRandom[SeedRandom[42, Method -> "Legacy"]; (*for reproducibility*)
  grf1 = Chop[GaussianRandomField[]]];

BlockRandom[SeedRandom[12345, Method -> "Legacy"]; (*for reproducibility*)
  grf2 = Chop[GaussianRandomField[]]];

GraphicsRow[{MatrixPlot[grf1, PlotLabel -> Style["1", Bold, 18]],
  MatrixPlot[grf2, PlotLabel -> Style["2", Bold, 18]]}]

Two random fields

Now sum the two fields:

MatrixPlot[grf1 + grf2, PlotLabel -> Style["1+2", Bold, 18]]

Sum of two random fields

The individual points of 1 and 1+2 are positively correlated:

points = Transpose[{Flatten[grf1], Flatten[grf1 + grf2]}];
Correlation[points][[1, 2]]
(* 0.712464 *)
ListPlot[points, PlotLabel -> Style["1 vs 1+2", Bold, 18]]

Scatter plot of individual points

  • $\begingroup$ Thanks for the idea. I will study it. I am working with this other article in which a covariance kernel for two processes is generated and then, the realizations of said process generated using Karhunen-Loève can be split in two. I need to do that but using 2D fields instead of 1D gaussian processes. I hope I have explained myself. $\endgroup$ Aug 6, 2021 at 7:11
  • $\begingroup$ Thanks for the accept but if someone else posts something that's more relevant to your particular objective, don't hesitate to accept that answer. $\endgroup$
    – JimB
    Aug 6, 2021 at 20:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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