关于2012年数学建模a题葡萄酒评价的分析

这道题目内容还是很丰富的，第一问两组品酒员评价结果有无显著性差异.

这就好比概率论上给你两组分别来自正态总体,让你评价有无差异类似，但是很明显现在并无法确定这两组数据是正态分布的（事实上程序跑出来不是正态分布).

这就用到成对数据检验概念，简单来说就是两组数据作差，可以认为是呈正态分布，然后

 

针对这题就是先求每一组对于不同酒的打分（组内取平均数即可），然后得到两组数据这两组数据做差，（所有差值平均-0）/标准差/sqrt(n)得到的数值与查表得到的数值检验，如果小就接受。

为了严谨也可以对得到的差值进行正态性检验。

matlab有函数可以，这是我搜到比较好的介绍网站http://www.cnblogs.com/djcsch2001/archive/2012/02/05/2339199.html。

当然也可以进行符号与 秩的检验，反正对于程序员来说就是使用不同的函数罢了，总之不要直接把两组数据当做正态分布。

第二问是我们使用了PCA+kmeans，步骤就是先降维，再根据葡萄酒理化指标降维得到的各主成分的得分来聚类，聚类玩好坏的评定根据葡萄酒的质量来.

第三问重点在于皮尔逊相关系数矩阵的求解,这也给求解多元线性回归提供了思路

  (1)若你确定一定有关系且是线性的，那就先画散点图，说明有线性关系，然后直接regress回归（这是建立在比较好的时候）,顺便还可以画残差图.

 (2)变量较多，且不一定相关，就先建立相关系数矩阵，找到相关性大于0.5，然后用自变量和找到的回归.

 (3)实在效果太差，可以考虑用stepwise逐步回归，这种做法是一定有结果，缺点是很有可能就是本来不相关的硬被拉起来找关系。

第四问其实比较扯淡，我们做法就是先拿出20组数据进行拟合，剩下几组数据根据拟合的系数算出数值，这组数值和原来的得分进行显著性检验，接下去采用第问方法。

总体来说：这里面部分算法接触过，所以写的还算顺利，问题在于竟然用spss可以不用写程序，用matlab三天编了近25k的代码，真是苦逼，一个附录弄了20页。

这里就展示由pca得到的数据进行kmeans..

clc;
clear;
time=0;
k=input('请输入分为几类:');
x=xlsread('problem two three red.xls','A1:B27');
for i=1:size(x,2)
x(:,i)=(x(:,i)-mean(x(:,i)))/std(x(:,i));
end
z=x(1:k,1:2);
z1=zeros(k,2);
while time<=1000
count=zeros(k,1);%%%每一个聚类中心有多少个属于它的点
allsum=zeros(k,2);%%%每一个属于该聚类中心的点到该一个聚类中心的横纵距离之和
num=[];%%%%记录下属于该聚类中心点的编号
temp=[];%%%%记录下每一个点到聚类中心的距离，找到最小的那个
sum(1,1)=0;
for i=1:size(x,1)
for j=1:k
temp(j,1)=sqrt((z(j,1)-x(i,1)).^2+(z(j,2)-x(i,2)).^2);%第i个点到第j个聚类中心的距离
temp(j,2)=j;
end
temp=sortrows(temp,1);%对于temp按照第一列降序排序
c=temp(1,2);%%%找到最小距离的聚类中心的编号
if temp(k,1)>sum(1,1)%%%到最远的聚类中心比现在的大
sum(1,1)=temp(k,1);
sum(1,2)=temp(k,2);%%%记录下是第几个点
end
count(c)=count(c)+1;%%%属于这个中心的点的数量加1
num(c,count(c))=i;
allsum(c,1)=allsum(c,1)+x(i,1);%%%求属于这个聚类中心所有点的横坐标之和
allsum(c,2)=allsum(c,2)+x(i,2);%%%求属于这个聚类中心所有点的纵坐标之和
end
for i=1:k
if count(i)~=0
z1(i,1)=allsum(i,1)./count(i);%%%每一个聚类中心所有点横坐标之和/属于该聚类中心的点数为新的横坐标
z1(i,2)=allsum(i,2)./count(i);
else%%%如果存在空聚类，则为较远的那个点
z1(i,1)=x(sum(1,2),1);
z1(i,2)=x(sum(1,2),2);
end
end
% z1(:,1)=allsum(:,1)./count(:);
% z1(:,2)=allsum(:,2)./count(:);
if (z==z1)%%%如果点不再变化
break;
else
z=z1;
end
time=time+1;
end
% plot(x(:,1),x(:,2),'r*');
% hold on;
% plot(z1(:,1),z1(:,2),'bo');
num(num==0)=NaN;
colorstring={'y*','m*','c*','r*','g*','b*','k*','r*','bo','cx'};
for i=1:k
disp(['第',num2str(i),'类为:',num2str(num(i,:))]);
end
H=[];
for i=1:k
for j=1:size(num(i,:),2)
if ~isnan(num(i,j))
hold on;
plot(x(num(i,j),1),x(num(i,j),2),colorstring{i});
H=[H,plot(z1(i,1),z1(i,2),'bo')];
end
end
end
legend(H([1]),'聚类中心');