%自己编的EOF程序
clear; clc
data=csvread('sst.csv'); % input data, each row is a observation, each col
% represents a time
threshold = 95.0; %input var

% first
[r,c] = size(data);
nobs = r; %观测站点数
ntimes = c;%测量的时间维

% each row is a observation
% each col is a time

% eof transformation
A = data*data';
[eof0,lamda0,explain0]=pcacov(A); %pca 变换， eof0即de空间场

pca0 = eof0'*data;

% select the eof vars based on the threshold

% base on the input explain to determin the nuber select
dd = find(explain0>threshold,1); % get the index，返回满足条件的索引值



rmse = sqrt(sum(lamda0(dd+1:end))/(r*c)); % compute the error,output
%output eof(space) and pca(time)
eof = eof0(:,1:dd); %every col is a spatial pattern ,output var
pca = pca0(1:dd,:);  %every row is time pattern, output var

% eof*pca can get the original data

% verify the results;
newdata = eof*pca ; % get the original data use the eof*pca

diff = data - newdata ; % compute the difference,
rmse2 = sqrt(sum(diff(:).^2)/(r*c)); % error ,the same as the rmse


实现步骤
1读取CSV文件，获取数据的行列数
2矩阵与矩阵转置的成绩（得到协方差矩阵）
3进行PCA变换，得到主成分分析的结果
功能：运用协方差矩阵进行主成分分析
   格式：PC=pcacov(X)
             [PC,latent,explained]=pcacov(X)
   说明：[PC,latent,explained]=pcacov(X)通过协方差矩阵X进行主成分分析，返回主成分(PC)、协方差矩阵X的特征值(latent)和每个特征向量表征在观测量总方差中所占的百分数(explained)。
4得到PCA原数据，主成分eof转置矩阵与原矩阵相乘
5通过threshold（95）作为判断值，得到满足百分数大于95条件数据的索引
6 通过索引值，计算满足条件的空间维和时间维上面的矩阵(eof空间维，pca是时间维)
7通过上面得到的eof和pca矩阵进行乘法运算得到模拟的新数据
8最后将新数据与旧数据相减，并计算误差。