Error message when trying to compute EOF with covariance matrix.

5 views (last 30 days)
I am trying to compute EOFs from a code I used before, but it does not seem to work here since my dataset has a too good spatial resolution (1/4 degree).
I want to compute the firsts 2 EOFs of mean June SLP northward of 70N for 43 years.
latitude has size 81 (4x20 = 80)
longitude has size 1440 (4x360 = 1440)
The data I want to compute the EOFs from has size 1440x81x43 (latitude x longitude x time)
When I want to compute the covariance matrix, I can't. I get the following error message:
Error using *
Requested 116640x116640 (101.4GB) array exceeds maximum array size preference (62.6GB). This might cause
MATLAB to become unresponsive.
Error in cov (line 155)
c = (xc' * xc) ./ denom;
Related documentation
%% Compute EOFs %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%I CAN'T ATTACH MY DATA SINCE IT IS 31.1mb AND THE LIMIT IS 5 MB
load('mslpJune_dt.mat'); %load detrended data
mslpJune_dt_reshaped = reshape(mslpJune_dt,[1440*81 43]);
mslpJune_dt_reshaped_trans = transpose(mslpJune_dt_reshaped);
%I CAN'T COMPUTE COV, ERROR MESSAGE FROM MATLAB THAT IT IS TOO BIG
C = cov(mslpJune_dt_reshaped_trans);
[eigenvectors,eigenvalues] = eigs(double(C),2);
%Principal components
pcs = mslpJune_dt_reshaped*eigenvectors;
%Scaling pcs
pcs1_norm(:,1)=pcs(:,1)/std(pcs(:,1));
pcs2_norm(:,2)=pcs(:,2)/std(pcs(:,2));
Does anyone know how to fix this issue and compute the EOFs?
Thank you

Accepted Answer

Christine Tobler
Christine Tobler on 6 May 2022
The problem is that the covariance matrix becomes very large here. Luckily, it's not necessary to compute this matrix explicitly: It's enough to compute a matrix xc such that C = xc' * xc. Then, instead of computing eig(C) it's equivalent to compute svd(xc). I explain a bit more in the code below:
mslpJune_dt_reshaped = randn(1440*81, 43);
mslpJune_dt_reshaped_trans = transpose(mslpJune_dt_reshaped);
% COV would return a matrix that's too large to store. Let's take the lines
% of code that COV does internally (see "edit cov") and see where that gets us:
x = mslpJune_dt_reshaped_trans;
[m,n] = size(x);
denom = m - 1;
xc = x - sum(x,1)./m; % Remove mean
%c = (xc' * xc) ./ denom; % Still can't afford this last one
xc = xc ./ sqrt(denom); % now, C = xc' * xc
% Compute the SVD of xc:
[U, S, V] = svd(xc, 'econ'); % Could also use svds(xc, 2) to get only the two largest singular values
% We know that xc = U*S*V', therefore C = xc'*xc = V*S*U'*U*S*V'
% and because U is orthogonal, U'*U cancels out and we have
% C = V*S*S*V', and therefore we can treat V as the matrix of eigenvectors
% and the square of the singular values on the diagonal of S as the
% eigenvalues:
eigenvectors = V;
eigenvalues = S*S;
% Verify these are really the eigenvalues and eigenvectors of C = xc'*xc:
norm(xc'*(xc*eigenvectors) - eigenvectors*eigenvalues)
ans = 3.2419e-11
It seems like what you're computing may be similar to pca - perhaps that function could be useful.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by