1. First, a question about the sense of EOF's: (Note: SCORE=PC, COEFF=EOF)

  1. How many values (numbers) are in your input data array? 240 x 144 (time x lon)
  2. How many values (numbers) are needed to build each term on the left? 240+144 (PCi(i,240) x EOFi(i,144)', i=1,2,...144)
  3. If 5 EOFs capture most of the data's variance, how much smaller (in the above sense) is the EOFxPC representation compared to the full data set? 240x144 - (240+144)x5 ??


2. Read in your field1 (let's call it x again).


Perform and display an EOF analysis of your first field (i.e. SST).
chen_EOF1.jpg
chen_EOF2.jpg


SSTAOr.jpgSSTARe.jpgTruncation1.jpgTruncation2.jpg

The following is the main part of my code for computing EOF and PC:
%% EOF=COEFF; PC=SCORE
sst = nc_varget('data.nc','sst');
sstm = mean(sst);
sstm = repmat(sstm,240,1);
ssta = sst-sstm;
 
[COEFF,SCORE,latent,tsquare] = princomp(ssta); % n x m; SCORE=PC, COEFF=EOF!!!
recon_ssta = SCORE*transpose(COEFF);
 
explained = latent/sum(latent(:)) .*100;
explained(1:10)
There is another way to do it; here is a simple demonstration:
X = [2 6 1 5 2;
     9 4 0 5 4];
X(1,:) = X(1,:)-mean(X(1,:)); X(2,:)=X(2,:)-mean(X(2,:));
 
% co-variance matrix
C=X*X'/5;
 
[EOF,E] = eig(C); % EOF: eigenvectors; E:eigenvalues
PC = EOF*X;
% reverse the order
E = fliplr(flipud(E));
lambda = diag(E); % retain eigenvalues only; exactly the same as variable 'latent' in the foregoing code.
EOF = fliplr(EOF);
PC = flipud(PC);
 
%% check
EOF*EOF'   % = I
PC*PC'/5   % = lambda
EOF*PC     % = X