数学建模MATLAB之分析法（一）

本篇分别有主成分分析，层次分析,聚类分析,因子分析的简单介绍以及在MATLAB（2018a）上使用的简单介绍 ( 附源码 )

主成分分析

主成分分析-百度百科
简单讲:主成分分析就是取出几个较少的变量尽可能多地反映原来变量的信息
(只看主要对结果有影的几个变量)
主成分分析适用于数学建模中筛选出少量变量对结果造成主要影响并用函数关系表征的情况

%% 主成分分析 (降维)
clc, clear

load example_1.txt  %数据导入（数据要求：前几列为自变量，最后一列为因变量）
data = example_1;

[m, n] = size(data); 
num = 3;   % 选取的主成分的个数

mu = mean(data); 
sigma = std(data);  %标准差
%z-score标准化方法适用于属性A的最大值和最小值未知的情况，或有超出取值范围的离群数据的情况
%标准化的新数据=（原数据-均值）/标准差
std_data = zscore(data);
b = std_data(: , 1:end-1);     % 四个变量x1, x2, x3, x4
r = cov(b);                                % 变量的协方差矩阵

% 运用协方差矩阵进行PCA
[PC, latent, explained] = pcacov(r);  %返回主成分(PC)、协方差矩阵X的特征值 (latent)和每个特征向量表征在观测量总方差中所占的百分数(explained)
% 新的主成分z1 = PC(1,1)*x1 + PC(2,1)*x2 + PC(3,1)*x3 + PC(4,1)*x4  , z2 = ...
f = repmat(sign(sum(PC)), size(PC, 1), 1);            %sum(PC)表示对矩阵PC的列求和
PC = PC .* f;


%1.普通的最小二乘法回归
regress_args_b = [ones(m, 1), b] \ std_data(:, end);   %标准化数据的回归方程系数
bzh = mu ./ sigma;
% 原始数据的常数项
ch10 = mu(end) - bzh(1:end-1) * regress_args_b(2:end) * sigma(end);
fr_1 = regress_args_b(2:end); fr_1 = fr_1';
% 原始数据的自变量的系数
ch1 = fr_1 ./ sigma(1:end-1) * sigma(end);
% 此时模型为 y = ch10 + ch1[1]*x1 + ch1[2] * x2 + ch1[3] * x3 + ch1[4] * x4
% 计算均方误差
check1 = sqrt(sum( (data(:, 1:end-1) * ch1' + ch10 - data(:, end)) .^2 ) / (m - n))


%2.主成分回归模型
pca_val = b * PC(:, 1:num);
%主成分数据的回归方程系数
regress_args_pca = [ones(m, 1), pca_val] \ std_data(:, end);
beta = PC(:, 1:num) * regress_args_pca(2:num+1);   %标准化数据的回归方程系数
% 原始数据的常数项
ch20 = mu(end) - bzh(1:end-1) * beta * sigma(end);
fr_2 = beta';
% 原始数据的自变量的系数
ch2 = fr_2 ./ sigma(1:end-1) * sigma(end);
% 此时模型为 y = ch20 + ch2[1]*x1 + ch2[2] * x2 + ch2[3] * x3 + ch2[4] * x4
% 计算均方误差
check2 = sqrt(sum( (data(:, 1:end-1) * ch2' + ch20 - data(:, end)) .^2 ) / (m - num - 1))

层次分析

层次分析-百度百科
简单讲：就是自己将两两影响因素作重要性对比得出判断矩阵，然后得出各个影响因素的权重
层次分析适用于数学建模多个变量对结果不同影响程度并用函数关系表征的情况

%层次分析法（AHP)
disp('请输入判断矩阵A（n阶）');
A = input('A=');
[n,n] = size(A);
x = ones(n,100);
y = ones(n,100);
m = zeros(1,100);
m(1) = max(x(:,1));
y(:,1) = x(:,1);
x(:,2) = A*y(:,1);
m(2) = max(x(:,2));
y(:,2) = x(:,2)/m(2);
p=0.0001; i=2; k=abs(m(2)-m(1));
while k>p
    i=i+1;
    x(:,i) = A*y(:,i-1);
    m(i) = max(x(:,i));
    y(:,i) = x(:,i)/m(i);
    k=abs(m(i)-m(i-1));
end
a = sum(y(:,i));
w = y(:,i)/a;
t = m(i);
disp(w);
%一致性检验
CI = (t-n)/(n-1);
RI = [0 0 0.52 0.89 1.12 1.36 1.41 1.46 1.49 1.52 1.54 1.56 1.58 1.59];
CR = CI/RI(n);
if CR<0.10
    disp('此矩阵一致性可以接受！');
    disp('CI=');disp(CI);
    disp('CR=');disp(CR);
end

聚类分析

聚类分析-百度百科
简单讲：就是多个影响因素按不同联系程度分类

%聚类分析
X;   %X为导入数据的名称(类型：数值矩阵)
X=X';    %将矩阵X转置，按情况看是否需要此语句
Y=pdist(X);     %X为要聚类的数组
SF=squareform(Y);
Z=linkage(Y,'average');
dendrogram(Z,0);
T=cluster(Z,'maxclust',3);

因子分析

因子分析-百度百科
可看作是主成分分析的推广，涉及的计算与主成分分析很类似

%因子分析
clc,clear
load ssgs.txt   %把原始数据保存在纯文本文件ssgs.txt中
n=size(ssgs,1);
x=ssgs(:,[1:4]); y=ssgs(:,5); %分别提出自变量x1...x4和因变量x的值
x=zscore(x); %数据标准化
r=corrcoef(x)  %求相关系数矩阵
[vec1,val,con1]=pcacov(r)  %进行主成分分析的相关计算
f1=repmat(sign(sum(vec1)),size(vec1,1),1);
vec2=vec1.*f1;     %特征向量正负号转换
f2=repmat(sqrt(val)',size(vec2,1),1); 
a=vec2.*f2   %求初等载荷矩阵
num=input('请选择主因子的个数：');  %交互式选择主因子的个数
am=a(:,[1:num]);  %提出num个主因子的载荷矩阵
[bm,t]=rotatefactors(am,'method', 'varimax') %am旋转变换,bm为旋转后的载荷阵
bt=[bm,a(:,[num+1:end])];  %旋转后全部因子的载荷矩阵,前两个旋转，后面不旋转
con2=sum(bt.^2)       %计算因子贡献
check=[con1,con2'/sum(con2)*100]%该语句是领会旋转意义,con1是未旋转前的贡献率
rate=con2(1:num)/sum(con2) %计算因子贡献率
coef=inv(r)*bm          %计算得分函数的系数
score=x*coef           %计算各个因子的得分
weight=rate/sum(rate)  %计算得分的权重
Tscore=score*weight'   %对各因子的得分进行加权求和，即求各企业综合得分
[STscore,ind]=sort(Tscore,'descend')      %对企业进行排序
display=[score(ind,:)';STscore';ind'] %显示排序结果
[ccoef,p]=corrcoef([Tscore,y])    %计算F与资产负债的相关系数
[d,dt,e,et,stats]=regress(Tscore,[ones(n,1),y]);%计算F与资产负债的方程
d,stats  %显示回归系数，和相关统计量的值