Karhunen Loeve Expansion in Matlab

  • #1

Main Question or Discussion Point

Hello everyone.

I am trying to generate the KL expansion of a stochastic process. I use a Monte Carlo sampling method to generate the process, which involves two random variables and I compare it with it's theoretical mean for 50 values of time and they look quite similar. Then, I calculate the autocovariance matrix, from where I extract the eigenvalues and eigenvectors which are used to calculate the new variable Y which is the stochastic process S in a base where the random variables are not correlated. However, if I calculate the mean of Y it doesn't look like the the mean of S at all.

Can someone please tell me what am I doing wrong?
Thank you very much for reading.

clear
close all
V = unifrnd(1,2,1,10000);
A = betarnd(2,2,1,10000);
t=50;
for i=1:t
S(i,:)=V*i+0.5*A*i^2;
theoreticalmeanS(i)=3/2*i+1/4*i^2;
meanS(i)=mean(S(i));
end
for i=1:t
for j=1:t
AutocorrelationS(i,j)=(S(i,:)*S(j,:)')/10000;
AutocovarianceS(i,j)=(S(i,:)*S(j,:)')/10000-meanS(i)*meanS(j);
end
end

[eigenvectors, eigenvalues]=eig(AutocovarianceS);
eigenvalues_column=eig(AutocovarianceS);
Y=eigenvalues*S;


Note: I am following the book of Alberto Leon Garcia (probability statistics and random processes for electrical engineering), pages 324-325.
 

Answers and Replies

  • #2
JMz
329
78
Not sure if I can help here (I am extremely familiar with the core of MATLAB, but not with its stats toolbox, which you appear to be using), but let me ask a question: For the data sizes you're dealing with, is there a reason not to use SVD to give you the PCs and their weights directly, instead of going through the covariance intermediate and finding eigenvalues?

My reason for asking is partly just that, in theory and in my experience, SVD is more stable and well-behaved than EVD (and there are deep reasons for that). For instance, with round-off, you might have a complex eigenvalue for the covariance matrix; SVD is simply not capable of that.
 
  • #3
Stephen Tashi
Science Advisor
7,176
1,315
for i=1:t
S(i,:)=V*i+0.5*A*i^2;
theoreticalmeanS(i)=3/2*i+1/4*i^2;
meanS(i)=mean(S(i));
end
What does ##S(i,k)## represent? The ##k## value of the ##i##th realization of the process? - or the ##i##th value of the ##k##th realization of the process? (If the process takes place in time, which index represents time?)
 

Related Threads on Karhunen Loeve Expansion in Matlab

Replies
3
Views
1K
Replies
1
Views
1K
Replies
2
Views
1K
  • Last Post
Replies
1
Views
2K
  • Last Post
Replies
3
Views
2K
  • Last Post
Replies
1
Views
2K
  • Last Post
Replies
1
Views
2K
Replies
8
Views
686
Replies
1
Views
2K
Top