Path: news.mathworks.com!not-for-mail
From: "Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid>
Newsgroups: comp.soft-sys.matlab
Subject: Re: correlation matrix estimation
Date: Sat, 17 May 2008 04:08:02 +0000 (UTC)
Organization: The MathWorks, Inc.
Lines: 81
Message-ID: <g0lln2$nqc$1@fred.mathworks.com>
References: <e963b995-e578-4c43-b848-711ec847fdf8@c65g2000hsa.googlegroups.com>
Reply-To: "Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid>
NNTP-Posting-Host: webapp-03-blr.mathworks.com
Content-Type: text/plain; charset="ISO-8859-1"
Content-Transfer-Encoding: 8bit
X-Trace: fred.mathworks.com 1210997282 24396 172.30.248.38 (17 May 2008 04:08:02 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Sat, 17 May 2008 04:08:02 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 1187260
Xref: news.mathworks.com comp.soft-sys.matlab:468962


msmscarlatti@googlemail.com wrote in message <e963b995-e578-4c43-
b848-711ec847fdf8@c65g2000hsa.googlegroups.com>...
> I have written some code which calculates an EWMA estimate for the
> correlation matrix of two time series. Suppose that X and Y are our
> time series vectors of length N (say N daily observations of two
> market interest rates) with X_1 the oldest observation and X_N the
> most recent one.
> 
> First the code transforms X and Y by subtracting off the average, i.e.
> 
> X --> X - E(X)
> Y --> Y - E(Y)
> 
> Then it concatenates X and Y forming a (Nx2) matrix, let's call it Z.
> Then it multiplies the ith row by lambda*(N-i). The covariance matrix
> is calculated by
> 
> Cov = (1-lambda) * Z' * Z
> 
> Finally, the correlation matrix is computed in the usual way by
> 
> Corr = S*Cov*S
> 
> where S is a diagonal matrix whose elements are given by 1/
> sqrt(Cov_ii).
> 
> My code seems to work fine for small N, but when I try it for very
> large N (several hundreds) I get nonsensical results. Have I missed
> something obvious?
> 
> I include my MATLAB code below.
> 
> function [CovMat,CorrMat]=EWMACovariance(X,Y,lambda)
> % COVARIANCE Computes covariances and correlations
> n=size(X,1);
> m=mean(X); % Compute the means
> X=X-m(ones(n,1),:); % Subtract the means
> m = mean(Y);
> Y=Y-m(ones(n,1),:);
> lambdaVector = zeros(n,1);
> for i = 1:n
>    lambdaVector(i) = lambda^(n-i);
> end
> X = [X Y];
> lambdaMatrix = diag(lambdaVector);
> X = lambdaMatrix*X;
> CovMat=(1-lambda)*X'*X; % Compute the covariance
> if nargout==2 % Compute the correlation, if requested
>    s=sqrt(diag(CovMat));
>    CorrMat=CovMat./(s*s');
> end
----------------
  There are a couple of aspects to your code that I question.  First, it is my 
understanding that in computing an exponentially-weighted moving average 
(EWMA) covariance, the sum of the applied weights should be precisely 1.  
This does not appear to be the case in your code.  In every row [X(i),Y(i)] you 
are first multiplying each of the two terms by lambda^(n-i) and then taking 
the product of those two terms, which gives you lambda^(2*(n-i))*X(i)*Y(i), as 
you calculate X'*X.  Finally you multiply by just a single (1-lambda).  If we 
sum these weights for all n products we get

 (1-lambda)*(1 + lambda^2 + lambda^4 + ... + lambda^(2*(n-1)))

which is very far from 1.  It looks to me as though you should either be 
initially multiplying the rows of [X,Y] by sqrt(lambda^(n-i)) or else waiting 
until you have computed the individual products X(i)*Y(i) before multiplying 
by lambda^(n-i).  (The latter seems preferable.)  However, even when this 
correction is made, the sum is still not precisely 1.  Instead of multiplying by 
(1-lambda) at the last step, you should multiply by (1-lambda)/(1-lambda^n) 
to get a sum of weights of exactly 1.

  My second worry is that, to be consistent with the EWMA "decay-with-time" 
philosophy, the means you have computed should be done the same way, 
namely with decaying weights that add up to 1.  If a member of a time series 
in computing covariance is to be given a low weight because it is old, then 
why shouldn't it have the same low weight in determining a mean value?  In 
your code you used the ordinary matlab 'mean' function which applies equal 
weights in obtaining means.

Roger Stafford