Problem with Mathematica's KarhunenLoeveExpansion's output

Question

I have a matrix called realizationMat, it contains 101 measurements of 90 realizations of a stochastic process. I will post instead an abridged 11*5 smaller example

meshgrid[x_List, y_List] := {ConstantArray[x, Length[x]], 
Transpose@ConstantArray[y, Length[y]]}     
{xx, yy} = meshgrid[Range[0, 1, .1], Range[0, 1, .1]];
kernel = Exp[(-1)*Abs[xx - yy]]; 
{eval, evec} = Eigensystem[kernel];
realizationNumber = 5;
longitudeOfEigA = Length[evec];
realizationNumber = 5;
evecRand = 0*ConstantArray[1, {longitudeOfEigA}];
realizationMat = 0*ConstantArray[1, {realizationNumber}];
For[i = 0, i < realizationNumber + 1, i++,
For[j = 0, j < longitudeOfEigA, j++;
evecRand[[j]] = 
eval[[j]]^0.5*RandomVariate[NormalDistribution[]]*
evec[[j]]; (*Print[i];*)realization = Total[evecRand, {1}]; 
realizationMat[[i]] = realization]];
VA1 = realizationMat[[1]];
VA2 = realizationMat[[2]];
VA3 = realizationMat[[3]];
VA4 = realizationMat[[4]];
VA5 = realizationMat[[5]];
KLVariablesB = KarhunenLoeveDecomposition[{VA1, VA2, VA3, VA4, VA5}];

In the documentation is written: rows of the transformation matrix m are the eigenvectors of the covariance matrix formed from the arrays ai. Also, The transformed arrays bi are uncorrelated, are given in order of decreasing variance, and have the same total variance as ai.

However, if I write ListPlot[{KLVariablesB[[2, 1]]}, Joined -> True], it doesn't look like an eigenvector at all, see figure 1

. Also, if I instead plot ListPlot[{KLVariablesB[[1, 1]]}, Joined -> True], it looks way more as an eigenvector, see figure 2, yet these aren't the eigenvectors of the original covariance matrix.

Can someone please tell me what is wrong with my code?

Best regards.

I am also attaching an image of the original eigenvectors of the kernel

Answer 1

Running your code as given

meshgrid[x_List, y_List] := {ConstantArray[x, Length[x]], 
Transpose@ConstantArray[y, Length[y]]}     
{xx, yy} = meshgrid[Range[0, 1, .1], Range[0, 1, .1]];
kernel = Exp[(-1)*Abs[xx - yy]]; 
{eval, evec} = Eigensystem[kernel];
realizationNumber = 5;
longitudeOfEigA = Length[evec];
realizationNumber = 5;
evecRand = 0*ConstantArray[1, {longitudeOfEigA}];
realizationMat = 0*ConstantArray[1, {realizationNumber}];
For[i = 0, i < realizationNumber + 1, i++,
For[j = 0, j < longitudeOfEigA, j++;
evecRand[[j]] = 
eval[[j]]^0.5*RandomVariate[NormalDistribution[]]*
evec[[j]]; (*Print[i];*)realization = Total[evecRand, {1}]; 
realizationMat[[i]] = realization]];
VA1 = realizationMat[[1]];
VA2 = realizationMat[[2]];
VA3 = realizationMat[[3]];
VA4 = realizationMat[[4]];
VA5 = realizationMat[[5]];
KLVariablesB = KarhunenLoeveDecomposition[{VA1, VA2, VA3, VA4, VA5}];

lp = ListLinePlot[{KLVariablesB[[2, 1]]}, PlotMarkers -> Automatic];

Manipulate[Show[Plot[a Sin[b + c x], {x, 0, 3 Pi}], lp,
  PlotRange -> {{0, 5}, {-6, 6}}],
 {{a, 6}, 0.5, 8}, {{b, 0.21}, 0, 4}, {{c, 3.1}, 1, 10}]

It appears you don't have enough resolution. A sine wave can be matched up to the output.

Note, there is a Stack Exchange site dedicated to Mathematica here: https://mathematica.stackexchange.com/