50 views (last 30 days)

Hi, I'd like to calculate a 'running' correlation coefficient of two variables. t seems that my problem is that the corrcoef(x,y) function writes a (2,2)-Matrix but I only wanna save the actual coefficient, i.e. (1,2) or (2,1), into my resuls. How can I do this? Any help is very much appreciated! > > > I've tried the following:

Here some fake data:

x = [1 3 4 8 9 12 13 14 40 30];

y = [2 3 4 7 8 9 13 14 41 32];

window=2; % window size for the correlation coefficient

the cyclist
on 3 Apr 2015

Roger Stafford
on 3 Apr 2015

Edited: Roger Stafford
on 3 Apr 2015

By "running correlation coefficient" I assume you mean you want the correlation computed for each possible length of x and y from 1:1, 1:2, ..., 1:N for some large maximum number N. If N is very large it becomes inefficient to repeatedly compute using 'corrcoef'. Fortunately it is possible to express the correlation in terms of five running sums which can greatly reduce the total computation required. You can see this form in the Wikipedia site

http://en.wikipedia.org/wiki/Correlation_and_dependence

near the end of the first section.

The following code carries out that computation in terms of the Wikipedia formula.

sx = x(1); sy = y(1);

sx2 = x(1)^2; sxy = x(1)*y(1); sy2 = y(1)^2;

c(1) = 1; % Avoid a NaN on the first step

for n = 2:N

sx = sx + x(n);

sy = sy + y(n);

sx2 = sx2 + x(n)^2;

sxy = sxy + x(n)*y(n);

sy2 = sy2 + y(n)^2;

c(n) = (n*sxy-sx*sy)/sqrt((n*sx2-sx^2)*(n*sy2-sy^2));

end

Note that for each step in the loop there are only about twenty operations performed for updating the correlation value c(n), whereas with a large value for N, each step would involve some multiple of n of such operations.

You should not necessarily take the above code literally, but its method may serve to give you a procedure for performing the computation more efficiently for your "running" type situation.

Roger Stafford
on 3 Apr 2015

Just change the way the five "running" sums are generated.

t = 10; % <-- You choose the width of the window

sx = cumsum(x); sx = sx(t+1;end)-sx(1:end-t);

sy = cumsum(y); sy = sy(t+1;end)-sy(1:end-t);

sx2 = cumsum(x.^2); sx2 = sx2(t+1;end)-sx2(1:end-t);

sxy = cumsum(x.*y); sxy = sxy(t+1;end)-sxy(1:end-t);

sy2 = cumsum(y.^2); sy2 = sy2(t+1;end)-sy2(1:end-t);

c = (t*sxy-sx.*sy)./sqrt((t*sx2-sx.^2).*(t*sy2-sy.^2));

Sign in to comment.

Sign in to answer this question.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.