多元统计分析2012

多元统计分析multivariate

data

analysisOutlines 介绍多元统计特征，包括多元数据之间的距离描述、多元数据图形表达（可视化）、多元正态分布、统计量及显著性检验。

结合主成分分析和因子分析，探讨线性维度减小方法。以多维标度分析为基础，讨论非线性维度减小方法。

介绍判别分布、聚类分析等多元数据分类和聚类方法。X=3+randn(100,1)Y=5+6*X+2*randn(100,1)Z=8-3*X+4*X.^2+3*randn(100,1)X1=X(1:50)X2=X(51:100)Y2=Y(51:100)Y1=Y(1:50)Z1=Z(1:50)Z2=Z(51:100)B=[X1

X2

Y1

Y2

Z1

Z2]save

B.mat

B多元数据距离测量

度量空间距离公理论。数学上，距离的真实测量，称为度量（metric）。它遵守三个公理：设dxy表示对象x和对象y之间的距离，则公理1：dxy=dyx。公理2：dxy≥0，当且仅当x=y时，dxy=0。

公理3：dxy≤dxo+doy，o为任意第三个对象，即满足三角形任意两边长度和大于

第三边长度。多元数据之间的距离:观测之间的距离（1）欧几里得距离（Euclidian

distance）但是在异质性的环境中，样本之间并不能直接可比。

标准化欧几里得距离（standardized

Euclideandistance）非欧几里得距离（non-Euclidean

distance）

欧几里得距离是一种具体的距离，因为它与我们距离的物理概念一致。非欧几里得距离则有不同方式，用于满足数学度量空间的测量。流行的非欧几里得距离包括：Bray-Curtis差异（dissimilarity）,L1距离（即曼哈顿距离），Jaccard指数（用于数据元素的有-无测量）。非欧几里得距离是异质环境数据之间差异测量的基础。Bray-Curtis差异有-无数据（presence-absence

data）的差异测量两总体或两样本之间距离总体与观测之间的距离基于比率的距离多元数据可视化指标变量制图

采用2、3个指标变量（index

variable），以投影图形方式描述4个以上的变量,称为投影寻踪

(projection

pursuit)。事实上，主成分分析等多元统计分析就是生成指标变量，以指标变量描述整个变量。指标变量制图优点是在二维或三维上，可绘制更多变量。其潜在缺陷是损失或模糊对象之间的可能差别。

MATLAB中，提供一组函数scatter3、gscatter、plot3和rotate3d等，用于绘制3D散点图。x=2+rand(20,1)y=3-rand(20,1)z=5+rand(20,1)subplot(1,2,1);plot3(x,y,z,"*b");

grid

onsubplot(1,2,2);scatter3(x,y,z,"*r")多元散点图

多元散点图（draftsman’s

plot）显示包含一数值图，是一个变量与其它变量之间对应的二维制图。它包括一组图形矩阵，也称散点图矩阵（scatterplot

matrix）。每个图形足够小，以在一个空间可以同时查看所有图形。这种制图优点是在二维空间绘制图形，适合显示变量对之间关系。缺陷是难以表达一个变量与其它多个变量之间关系，不适合显示多变量之间关系。

MATLAB中，采用函数plotmatrix、gplotmatrix绘制散点图矩阵x

=

randn(50,3);y

=

x*[-1

2

1;2

0

1;1

-2

3]";plotmatrix(y,"*b")单个数据点表示

这种图形显示每个对象的多元数据，每个变量是一个具有大小和方向的符号。不同变量有不同的

符号特征，这些符号组合为人脸或星形玫瑰图。

这种方法是1979年Fienberg提出的。它的优点是可以同时显示全部变量情形，缺陷是图形显示可能

依赖对象获取的顺序，而且只适合于小数据量的

多元数据（样本数小于10，变量数小于15）。

以MATLAB数据carsmall.mat为例，采用glyphplot函数，直接绘出多元数据的星形玫瑰图或Chernoff脸部图load

carsmall

X

=

[Acceleration

Cylinders

DisplacementHorsepower

MPG

Weight];glyphplot(X,"grid",[2,5])

Chernoff脸谱图用于说明多元数据的趋势，广泛用于聚类、判别和时间序列分析。选择face参数可绘制出脸谱图。不同部位有不同参数值，它们决定脸型，如描述耳可以有耳的水平面位置和半径；眼眉有高度、坡度和长度；眼睛有高度、大小、分隔距离、偏心率等；脸有宽度、半脸高度偏心率；嘴有中心位置、曲率、长度和张开程度鼻有宽度和长度。采用不同变量表达脸谱参数，比较脸谱，显示出样本数据类型:glyphplot(X,

"Glyph","face",

"grid",[2,5])变量剖面图

变量剖面图是以横轴为变量组，逐个表示变量；变量数值是纵轴，连接对象或以顺序条形图表示每个对象。该图优点是每个对象连接成一个折线。缺陷是不同量纲的变量不具有可比性，容易误导读者。

MATLAB中，还提供一个函数slice，可以选择坐标位置，完成3维实体的切片表示。首先，创建一个3D栅格:[x,y,z]

=

meshgrid(-3:.1:3,-3:.1:3,-3:.1:3);[n,d]=size(x(:));建立一个3D分布的数据：a

=

(2*pi)^(3/2);si=

(x.^2

+

y.^2

+

z.^2);prob

=

exp((-.5)*si)/a;选择x=0,

y=0,

z=0的平面位置进行切片：slice(x,y,z,prob,0,0,0)xlabel("X"),ylabel("Y"),zlabel("Z")Andrews曲线X=[7

15

5

15

2

9

5;

7

7

3

24

25

5

3;

5

25

7

6

15

9

2]t

=

linspace(-pi,pi)f1=7/sqrt(2)+15*sin(t)+5*cos(t)+15*sin(2*t)+2*t)+9*sin(3*t)+5*cos(3*t)f2=7/sqrt(2)+7*sin(t)+3*cos(t)+24*sin(2*t)+25*t)+5*sin(3*t)+3*cos(3*t)f3=5/sqrt(2)+25*sin(t)+7*cos(t)+6*sin(2*t)+15*t)+9*sin(3*t)+2*cos(3*t)plot(t,f1,".",t,f2,"*",t,f3,"o")legend("F1","F2","F3");xlabel("t")平行坐标制图

MATLAB中，平行坐标制图函数

parallelcoords也可绘制n×p的多元数据矩阵

X。在笛卡儿坐标系中，坐标轴是正交的，如此可以表达3D。扩展笛卡儿坐标系，以一组平行坐标轴，则能够可视化多元数据。

以MATLAB中的carsmall.mat数据为例，取其中部分变量和观测结果，parallelcoords函数的运行结果…load

carsmall

X

=

[Acceleration

Cylinders

DisplacementHorsepower

MPG

Weight];X1=X(30:40,:)parallelcoords(X1)合成数据表达与饼图

合成数据是一种特定数据，表示一个样本

可能由多个成分组合而成。MATLAB中，pie、pie3函数用于表达合成数据。X

=

[19.3

22.1

51.6;

34.2

70.3

82.4;

61.4

8290.8;

50.5

54.9

59.1]x

=

sum(X)A=zeros(size(x))[c,offset]

=

max(x)A(offset)

=

1h

=

pie(x,A);

colormap

summer

得出一个饼图，其最大比例部分与其它部分分离。以下给它增加注释：textObjs

=

findobj(h,"Type","text");oldStr

=

get(textObjs,{"String"});val

=

get(textObjs,{"Extent"});oldExt

=

cat(1,val{:});Names

=

{"X:

";"Y:

";"Z:

"};newStr

=

strcat(Names,oldStr);set(textObjs,{"String"},newStr)val1

=

get(textObjs,

{"Extent"});newExt

=

cat(1,

val1{:});

offset

=

sign(oldExt(:,1)).*(newExt(:,3)-oldExt(:,3))/2;pos

=

get(textObjs,

{"Position"});

textPos

=cat(1,

pos{:});textPos(:,1)

=textPos(:,1)+offset;set(textObjs,{"Position"},num2cell(textPos,读者可以自己尝试pie3函数绘制3维饼图。

此外，MATLAB还提供bar，bar3，bar3h，barh，comet，comet3，contour，contour3，contourf，cylinder，errorbar，feather，fill，

fill3，lmesh，meshc，

meshz，pie，pie3，plot，plot3polar，quiver，quiver3，ribbon，scatter，scatter3，semilogx，semilogy，

stairs，stemstem3，surf，surfc，surfl，surfnorm和waterf函数，它们适合于不同制图情形，读者可以根据需要，借助它们合理地完成数据可视化。

2009年版的MATLAB统计工具库，还有

scatterhist等图形函数。多元数据分布

多元分布不是新概念，只是单变量概率分布的扩展。对应单变量的二项分布、正态分布和

beta分布，多元分布中有多项式分布（multinomial

distribution）、多元正态分布（multivariate

normal

distribution）和Diri布。多元正态分布：也称为多元高斯分布（multivariate

Gaussian

distributio或矩阵正态分布

多元正态分布是多元统计学中最受欢迎的分布，常常通过数学变换将数据拟合为正态模型。注

意，这种变换依然有争议，因为有时变换是由

于方便和简化计算，结果却可能误导或难以反

回原始变量。尽管如此，多数推断性多元统计

学利用多元正态分布。

从单变量正态分布，到两变量正态分布的扩展过程，是理解多元正态分布的基础。单变量X的正态分布函数：

MATLAB中的多元正态分布：MATLAB中，函数

mvnrnd产生多元正态分布的随机矩阵。mvnpdf用于计算多元正态概率密度函数。2009年版的

MATLAB中，提供gmdistribution函数，它采用期望最大化（ExpectationMaximization，EM）算法拟合数据，由高斯混合模型生成一个多元正

态分布的混合成分集。

以下是一个MATLAB的多元正态概率密度函数实例mu

=

[0

0];sigma

=

[.25

.3;

.3

1];x1

=

-3:.2:3;

x2

=

-3:.2:3;[X1,X2]

=

meshgrid(x1,x2);F

=

mvnpdf([X1(:)

X2(:)],mu,sigma);F

=

reshape(F,length(x2),length(x1));surf(x1,x2,F);caxis([min(F(:))-.5*range(F(:)),max(F(:)axis([-3

3

-3

3

0

.4])xlabel("x1");ylabel("x2");zlabel("概率密

2009年版的MATLAB中，提供函数mvncdf计算多元正态分布的累积分布。Dirichlet分布多元数据的统计量分布1．Wishart分布

在多元统计中，Wishart分布协方差矩阵分析具有重要意义。Wishart分布名称来自JohnWishart，是单变量χ2分布的多元扩展或gamma分布的非整数自由度情形。

MATLAB中，函数wishrnd用于产生Wishart随机矩阵。注意，函数中的sigma参数即为协方差矩阵Σ。此外，MATLAB中，还有一个逆Wishart分布（inverse

Wishart

distribution）函数iwishrnd用于生成一个逆Wishart分布的随机数矩阵。2．Hotelling

T2分布

Hotelling

T2统计因Harold

Hotelling而得名是Student’s

t统计量在多元假设检验中的扩展。

MATLAB中，函数mvtrnd产生一个多元t分布的随机矩阵。2009年版的MATLAB中，函数

mvtpdf和mvtcdf用于多元t分布的概率密度与累积概率分布。3．Wilks分布

统计学中，Wilks

λ分布是F-分布的扩展。因

Samuel

S.Wilks而得名，是用于多元假设检验的概率分布，尤其是似然比率检验。MATLAB中，

Wilksλ是用于多元方差分析（MANOVA）函数的检验统计量，用于检验关联变量组合的可识别群的均值是否有差异。例如，我们可以将一个班级的学生分为城市和农村来源的两组，考虑8门专业课程（8个相关变量的组合），采用

Wilksλ比较这两组学生的成绩组合均值是否有显著差异。1．多元两样本均值比较：Hotelling-T^2检2．Bonferroni校正

多元试验与观测产生大量数据，不同变量类型、数据类型、观测数量和分析目的，决定不同统计量的采用。比较单变量与多变量显著性检验，两者有区别。设某单变量参数的零假设在0.05显著性水平上是真实，我们可以得到95%置信概率水平上结论是正确的。然而，对双变量，获得正确结论的概率是0.95×0.95=0.90。依此类推，p个变量的情形将是0.95p，即至少出错一次的概率是1-0.95p。这样，就产生了一个难以相信的结果：当p越大，正确结论的置信度就越低。这显然不符合实际情形，因此，意大利数学家Carlo

EmilioBonferroni提出Bonferroni校正。

多元显著性检验的Bonferroni校正是简单地给检验值乘以一个p值（如果有p次独立检验，则每次乘以1/p），将校正后的数值与0.05显著性水平比较，确定是否显著。如果校正值仍小于0.05，则拒绝零假设。

Bonferroni校正也存在一些有争议的问题，有些应用校正情形需要慎重。3．多元多样本均值比较Bryan

F.J.Manly归纳有4个统计量可以描述：（1）Wilks的Λ（lambda）检验对应的统计量.（2）Roy最大特征根检验（3）Pillai迹统计量。（4）Lawes-Hotelling迹统计量.PP.318-3194．多元两样本与多样本方差比较MATLAB中多元统计量比较分析1．MATLAB中多元方差分析m多an元ov方a是差一分种析“（门mu户lt”iv检ar验ia（tegantaelwyasyitseosft）va。ri如an果ce多，元MAFN-O检VA验）是。显著的，则可以考虑个别单变量分析。manova也可以检测组合差异。当变量之间存在多重共线性时，单纯的单变量分析难以发现这种差异。但是manova对异常值非常敏感，尤其是样本数量比较小时。manova假定因变量之间存在线性关系，也不能给出主要因子与重复因子之间的交互影响效应。变ma量no均va满用足于方同差时分分析析相多关个假相设互。依其赖基的本多假元设情是形，。待它分假析定数：据①集所X有满足：每个变量数据集/组（group）对应的总体服从正态分布；每个总体的方差、协方差矩阵相同；所有的观测相互独立。②因变量之间是线性或多重共线性。MATLAB中，函数manova1用于实现单因素多元方差分析。组MA均TL值AB为中基，础还，有生一成个谱函系数图ma.novacluster。它在MANOVA之后，以2．MATLAB中均值多重比较检验

MATLAB中，multcompare函数用于完成多重均值比较。3．多元等方差检验多元等方差检验，也称为多元方差齐性（homogeneity

of

variance）检验，即检验k群方差是否相同。它包括两种检验方法，Bartlett检验和

Levene检验，检验的零假设都是k群具有相同的方差。MATLAB中，采用vartestn函数实现这两种检验。（1）Bartlett检验（2）Levene检验减小多元数据维度：线性方法

多元分析的基本问题是Richard

Bellman提出的维度灾难（curse

of

dimensionality）。在高维空间的观测（采样）比低维空间的更稀疏和更小代表性；高

维空间中，服从某分布的典型数据点被采样的频率更小；数据集中的某些维度及相应的测量特征在建

模上也没有意义，即存在冗余。因此，多元数据分

析中，常包含有维度减少（dimension

reduction）方法。维度减少是以更小维度空间，更少的变量表

示多元数据。其目的是探索高维数据（多元数据），找出数据模式或结构；以二维或三维散点图可视化

多元数据。

经典线性维度减少方法包括：多元特征选择和多元特征转换。其中多元特征转换方法可分为主成分分析（principal

component

analysis）和因子分析

(factor

analysis)。

主成分分析面对的问题是，对大数据量测量变量，在数据集中含有多少冗余属性信息？该问题陈述包括三个特定问题：

①相关矩阵包含的线性关系模式是什么？设属性信息矩阵X，该问题的回答等价于评价X中多重共线性（multicol-linearity）程正是主成分分析目标。

②能否给相关矩阵一个最为简洁的描述？该问题的回答等价于，消除矩阵X中的冗余信息，正是因子分析的定位。

③一个给定的相关矩阵中，是否潜藏着某种合理维度？该问题的回答等价于假设检验，是因子分析要解决的问题。

主成分分析与因子分析目的和术语有大量重叠，两者都是维度减少技术，都是以更少新变量替

代原有多个观测量。有些人甚至认为两者没有

本质区别，因子分析中的一些算法也包含在主

成分分析中。

但是，两者目的和基本模型还是有不同。概括地说，如果需要更少的维度去可视化数据，简化和近似数据，应采用主成分分析；如果需要一个解释性模型，以理解数据之间相关性，则需要因子分析。

可以认为，因子分析是主成分分析的扩充。主成分分析通过坐标旋转，解释原始p个变量的方差分布。因子分析解释原始p个变量与公共因子之间相关性。主成分分析不依赖任何特定统计模型，而因子分析则以模型为基础。1．主成分分析概念与原理

1901年，Karl

Pearson首先提出主成分分析（PCA）概念与描述，直到1933年Hotelling才完成PCA计算方法。计算机技术出现后，才有可能分析更多变量的主成分问题。

作为最简单的多元分析方法之一，PCA通过检测变量之间线性独立性，以严格定量分析为基础，简化或合并相似特征的变量，或找出众多变量中对系统行为起控制作用的变量。用一组少量的变量，替代原来众多变量集，减少信息冗余。以PCA方法产生的新变量集，是原有数据变量集的子集，称为主成分（principalcomponents，PC）。主成分可以表达原有数据集特征，每个主成分是原有数据集中变量的线性组合。每个主

成分之间相互正交，这样就没有冗余问题。这正是多

元分析所涉及的如何避免重复属性信息的重复计算问

尽管线性分析中，以多变量为基础，建立正交基方法非常多（详细资料，请参阅线性代数相关内容），但是主成分分析的特殊是在于：

①线性空间中，第一主成分是一个坐标轴。将每个观测（样本）投影在该坐标轴，其结果形成一个新变量。在所有可能选择的坐标轴上，该变量方差最大。

②第二主成分是垂直于第一主成分的坐标轴。将每个观测（样本）投影在该坐标轴上，产生另一个新变量，该变量的方差在所有可能选择的第二坐标轴中是最大的。

理论上，主成分包含的变量集可以和原始数据中所有的变量数目一样多，但是实际上最初主成分之间方差和应当是不小于原始数据总方差的80%。

主成分常用于环境系统变化驱动力分析，如环境变化、土地利用与覆盖变化的主要驱动力分析。

应当提醒，主成分分析不是万能技术。主成分分析期待主成分的方差集中在少数成分上，这样可以采用少数具有较大方差的Z值，表达变量X，将X转换为更少变量的Z成分集。如果原始数据的变量X之间没有任何相关性，是相互独立的，则主成分分析没有任何意义。只有在原始数据中变量之间存在非常高的相关性，主成分分析才可能将原始数据中的众多变量（甚至变量数p＞20），缩减为2～3个主成分。变量之间强相关，意味着原始数据中存在非常大的冗余。

从这种意义上，目前土地评价、生态评价、环境评价工作中，凭个人喜好和经验选择变量（因子），不考虑变量之间相关关系，直接加权评分，得到评价结论的方法，其过程与结论都是值得怀疑的。如果变量之

间有强相关性，意味着工作重复和忽略了真正的重要

变量，尽管这些变量数值可能微不足道。个人喜好是

难以被人重复的，因此，评价工作中就将出现不同观

察角度，选择不同变量和权重。从统计学角度，这些

评价方法不严谨,也没有说服力。也许主成分分析考虑因子之间相关性，能让类似的评价工作更经得起讨论。3．主成分分析计算过程

主成分分析过程：①标准化变量X，使其均值为0、方差是1，以消除量纲影响。②计算标准化变量的协方差矩阵C。对标准化变量，矩阵C是相关矩阵。③找出矩阵C的特征根λi及其对应特征向量ai=[a1i,a2i,…，a第i个主成分系数是ai的元素，λi是它的方差。④舍弃任何小方差对应的主成分。例如，计算了10个变量的数据集，如果前三个主成分方差和占90%比例，则舍弃其它7个成分。4．主成分计算的两种方法这里进一步介绍主成分计算的两种方法：

①使用样本协方差矩阵计算主成分。以协方差矩阵为基础的主成分分析优点是，统计推断方法相对容易，易于理解。

②使用样本相关系数矩阵计算主成分。当原始维度上方差差异明显。例如，一些变量的方差相对其它变量的方差非常大，可以采用相关系数矩阵分析主成分。5．主成分几何意义主成分几何解释是以相关矩阵C的特征函数为基础。Daniel

A.Griffith以Z1和Z2为原始图形轴，主要特征向量表示一条穿越散布z-得分点的直线，z-得分点到该直线的垂直距离是最小的。所谓的主成分得分（principalcomponent

scores）是合成变量（主成分）的度量。主特征值给定了直线的相对长度。第二特征向量必须垂直第一特征向量，确保两个主成分之间相对独立，其相对长度由第二特征值给定。变量Z1与Z2之间相关系数是特征向量与原始z-得分坐标之间形成的夹角，该系数用于创建合成变量集F。假如变量x1与x2之间是完全随机（无相交），则相应相关系数是0(ρ=0),Z1与Z2之间也是随机的。这种情形下，z-得分值点散布是一个圆环，特征值等于1，特征向量是z-得分数轴。如果x1与x2之间是部分相关，则相应相关系数是0＜ρ＜1，z-得分值点散布是一个椭圆环，特征值等于椭圆半轴，特征向量与z-得分数轴有一个夹角。如果x1与x2完全相关，相关系数为ρ=1。z-得分值点散布椭圆变成直线，相对长度等于2，第二特征向量的相对长度是0。在三维情形下，z-得分值点散布是特征向量为轴的椭球或球。MATLAB中主成分分析实现

MATLAB中，princomp函数计算数据集的主成分，基本要求是需要一个原始的数据集，以供分析。假如缺乏原始数据，但是有样本数据之间协方差或相关系数，也可以采用pcacov函数完成主成分分析。X=3+randn(30,1)Y=5+6*X+2*randn(30,1)Z=8-3*X+4*X.^2+3*randn(30,1)X1=X(1:15)X2=X(16:30)Y2=Y(16:30)Y1=Y(1:15)Z1=Z(1:15)Z2=Z(16:30)B=[X1

X2

Y1

Y2

Z1

Z2]save

B.mat

Bboxplot(B)stdr=std(B)meanr=mean(B)meanr=repmat(meanr,15,1)sr=(B-meanr)./repmat(stdr,15,1)[coefs,

scores,

variance,t2]=princomp(srI=coefs"*coefsplot(scores(:,1),scores(:,2),"+")percent_explained=100*variance/sum(varipareto(percent_explained)[st2,index]=sort(t2,"descend")extreme=index(1)extreme=index(2)corrmatrix=corrcoef(sr)corrmatrix=flipud(corrmatrix)imagesc(corrmatrix),colormap(gray)Colorbar[coefs,scores,variances,t2]=princomp(sr

主成分的表达:采用图形,显示出不同主成分上的变量分布及其相关性,见PP338-340.因子分析 环境问题中，一个区域污染可能是大气污染，也可能是水污染和土壤污染，因此环境监测变量可以归于大气因子、水因子和土壤因子。就象竞技体育的十项全能运动一样，有的是力量型，有的是速度型项目。因子分析是以变量之间组合（相互依赖关系）因子为基础，以因子代替变量，达到减少维度的目的。

归纳起来，因子分析用于拟合多元数据模型，估计变量之间某种相互依赖关系。因子分析模型中，原始数据中测量变量依赖于更少的未测量因子，因为每个因子可以影响更多的变量，它们是公共因子（common

factor）。假定每个变量与公共因子之间有线性关联，相关系数称载荷（loadings），每个测量变量也包括一个独立随机变异成分，称个性变差（specificvariance），因为对每个变量是特定的。2．因子分析步骤与方法因子分析包括3个步骤：（1）确定临时因子载荷aij：临时因子F*之间可能有相关，或个性变差之间存在相关。aij的数值越接近于零，则xi与Fj之间越无关。aij的数值绝对值越大，则xi与F间相关关系越强。（2）因子旋转：方差最大旋转（varimax

rotation）方法是非常好的标准因子分析方法，由Kaiser提出并加以改进。它基于一个假设，可借助j因子的载荷平方方差解释j因子，即，如果方差大，则趋向于0或1。方差

最大旋转就是最大化所有因子方差和，最后找出一组

不相关的因子。在方差分析中，放弃因子之间相互独

立思想，使因子载荷尽可能简单。

（3）计算因子分值：因子分值用于校验变量与因子之间相关程度，采用类似主成分分析的方法计算它们，称为主成分因子分析（principal

components

factoranalysis）。3．MATLAB中因子分析实现

MATLAB中，采用函数factoran实现极大似然公共因子分析。

[lambda,psi,T,stats,F]=factoran(sr,2,"sco

ssion")具体解释见PP343-344减小多元数据维度：非线性方法

非线性维度缩减是将高维空间映射到低维空间，常用方法包括多维标度分析（multidimensional

scaling）,有些中文资料中也译为“多维标度”。多维标度分析概念

可视化数据的重要目标之一是直观感觉数据点之间相对距离的远近，散点图就是一种良好的

工具。有些统计分析情形中，数据根本就不是

以点的形式存在，而是以两两相似或差异的方

式出现。没有任何数据点，就不能绘制散点图。退一步说，即使有些数据是以点的形式出现，

不是两两点之间距离关系表示，散点图也可能

没有意义。因为有些数据之间两点邻近的度量

不是采用欧几里得距离，而是必须使用非欧几

里得距离，如城市街区距离。对这种情形的多

元变量，很难用一般图形刻画出观测（样本）之间的相似和差异关系。因此，减少维度，又

尽量不损失信息量的多维标度分析成为一种必

要手段。

多维标度分析（multidimensional

scaling，MDS是适应上述情形的一组方法集，它可视化点对

之间在多种距离上是如何相近或差异，在更低维度空间上表示数据。MDS也不必要求采用原始数据，只需要两两距离或相异的矩阵。有些中文多元统计分析著作中，多维标度分析被称为多维标度法。

多维标度分析最初可能是用于数学心理学领域。

1952年，Torgerson首先提出多维标度分析的数学模型。目前它已经广泛应用在交通、生态、

地学与社会经济领域。多维标度分析是通过模

型等方法，如采用距离、相似系数、稀疏系数

等,描述对象之间亲缘关系变量，把高维对象转换为低维情形，在较低维度空间标度、描述更

高维度的对象，以揭示数据中隐含结构。

与聚类分析比较，多维标度分析弥补了聚类分析上的不足。聚类分析归类相似样本，以谱系图反映样本之间亲缘关系，将高维样本纳入在一维谱系中，这种简化可能导致数据失真，使聚类结果具有一定主观性。多维标度分析则将高维对象简化为较低维度，却最大限度地保持对象之间原有相关关系。与主成分分析比较，两者都是降低样本空间维度。多维标度分析更趋向以图形方法，显示对象之间关系。这种图形可以是一维、二维、三维或更高的维度。一个多维标度环境问题实例

假设以天津n条街道/道路为对象，考虑它们之间环境质量差异。任意两条街道（r,s）之间差异

(dissimilarity)（有些中文教程，也将dissimilar译为不相似。）的描述量为,表示r与s之间没有差异，表示r与s之间环境质量完全不同。两两街道比较的结果，的取值有n(n-1)/2种可能。如果不同的判断角度，则有i角度的差异描述结果。

多维标度分析的基本思想，是寻找到一个更小维度空间，通常是欧几里得空间。以空间中的每个点表示一个对象（街道），以点之间的距离匹配初始的差异。两条相似的街道由多维空间中两个距离相近的点表示，两条相异的街道则由多维空间两个距离较远的点表示。多维标度分析类型（1）经典MDS：

（2）度量最小二乘标度分析（Metric

leastsquares

scaling）：（3）非度量MDS：4．邻近式邻上近的（接pr近ox。im对it象ie之s）间的邻字近面的理测解量，有就两是者在类空型间：、相时似间或其它方（similarity，）与差异（dissimilarity）。距差离异、可城以市采街用区距距离离来、表M达ah，al包an括ob欧is几距里离得、距Mi离nk、ow加sk权i度欧量几距里离得等。MATLAB中经典多维标度分析的实现

MATLAB中，采用cmdscale函数计算经典多维标度问题。cmdscale函数采用相关点距离矩阵为输入，创建点构图。理论上，这些点是在2维或3维的，以它们之间的欧几里得距离为基础，建立原始距离矩阵。这样，cmdscale创建点散点图可以可视表达原始距离。其句法是：

[Y，e]=cmdscale(D)D=pdist(B,"euclidean")[Y,e]

=

cmdscale(D)[e

e/max(abs(e))]dim

=

sum(e

>

eps^(3/4))maxerr2

=

max(abs(pdist(B)

-

pdist(Y(:,1:见PP348-349非度量多维标度分析PP.349-353Procrustes分析Procrustes是古希腊传说中Attica的强盗,他将被他抓到的人放在一张铁床上,身高比床长的人,被其砍去长出的部分；身高比床短的人,被其强行拉长。理解了Procrustes，就能知道Procrustes分析是比较形状的统计过程。Procrustes分析，有些中文翻译为“普鲁克分析”，起源于生物学研究中的动物形态

测定，但是很快就广泛用于考古、天文、土木工程、地理、网络设计和物理化学。相信它在环境污染模

式和环境模拟上也将有应用价值。

Procrustes分析常用于统计分析的预处理步骤，它需要应用一个标记数据集（landmark

dataset）,其中有意义的重要总体特征被标记为具体的几何位置。借助保留形状的欧几里得空间变换，最小化在位置、旋转和尺度上的变异，Procrustes分析计算出两个数据集中的标记最好拟合重叠情形。MATLAB中，函数procrustes

MATLAB中，函数procrustes执行procrustes分析，其基本句法为：d

=

procrustes(X,Y)判别分析1分类与判别分析概述1．分类的分数类据是模已型知,数称据为类分型类，器将（新cl数as据si归fi类er的s）过。程一。个我分们类将器类的型建响造应，首先需要训练数据，该数据的类型是已知的。在一定置信水平上，分类器将新的待分类数据指定不同的类型。与聚类（cluster）不一样，分类是事先已知类型，而聚类是事先不知道类型。（分d类is方cr法im可in归an纳t

a为na两ly种si：s）①。参它数以方训法练，数如据判为别基分础析,拟合参数模型，分类新数据。②非参数方法，如分类树（classification

trees）。类似于非线性回归中讨论的回归树方法。2．判别分析（判d别is函c-数ri分mi析na（ntdiasncarliymsinsa）nt。f判un别ct分io析n

a是na多ly元si统s）计，分简析称中判，别用分于析判别某观测或样本所属类型的一种统计分析方法。注意，不管采用什么样的判别分析方法，有两个基本假设是必需的：（1）对所有类型，类内的协方差矩阵是相等的。（2）对显著性检

验，类内数据服从多元正态分布。另外，还要求误差独立。

面对的问题是用某种方法把研究对象分成若干类型的情形下，判断和确认某新观测的类型归属。

判别分析首先指出最佳判别变量，以提供两个类群的最佳区别；然后,确定一个新观测（样本）与已知类型的接近程度，即建立判别函数。

判别函数是判别变量的线性组合，判别函数的值称为判别得分。判别函数目标是使判别组间平方和与组内平方和之比值为最大，再利用判别得分将观察值归类。

最后，指定一些判别规则，判定样本或观测的类型归属。判别规则可以是随机的，如Fisher准则和

Bayesian准则；也可以是确定的,只考虑判别函数值（得分）的大小。

判别分析有多种类型，按判别的组数分为两组判别分析（即样本只能属于A或B组）和多组判别分析；按模型分为线性判别和非线性判别；按处理过程分为逐步判别和序贯判别。

与主成分分析、因子分析比较，判别分析不需要标准化数据，判别分析结果不受单个变量尺度的影响。

判别分析更类似于回归分析，探索独立变量的线性组合，获取独立变量方差的最好解释。类似于回归问题，判别分析主要是针对区间变量和比率变量。当对数回归与判别分析的假设都满足时，它们趋于得到同样结果。与聚类分析比较，在判别分析中，至少有一个已经明确知道类别的“训练样本”。

利用这个数据，建立判别准则，判别未知类别观测值的归属。聚类分析中，一般人们事先并不知道或一定要明确应该分成几类，完全根据数据来确定。

环境统计中，判别分析面对环境质量分级或污染性质分类等具体情形，目的有二：

①评价区域环境对象分类或分区的合理性。多数环境问题是区域性问题，以多元分析为基础，指明区域之间的环境差异。因此，判别分析应当回答的问题是，p个预测变量最有力的线性组合关系是什么？以此为基础，才能建构环境特征归类或区域的最大分割。

②对新采样的环境样本，或面临新的环境区域，与以往类型进行比较，判断样本或区域属性，属于哪种类型。即提出适当计算规则，有效地归纳新样本到既有类型中。在区域问题处理上，判别分析也可以用于遥感图象数据解释与计算处理。2判别分析方法与准则

判别分析方法很多，包括Mahalanobis距离判别法、Fisher判别法、Bayesian判别法、典型判别函数、逐步判别函数分析、观测

的Jackknife分类和Logistic回归等。1．使用Mahalanobis距离的判别分析

最简单的判别分析办法，就是计算到哪个类型中心距离最近，按最近距离判别对象属于哪一类。通常使用的距离是Mahalanobis距离。

判别函数(discriminant

function)是用来比较个中心距离的数学函数。这种根据距离远近进行判别的方法，原理简单，直观易懂。利用

Mahalanobis距离进行判别分析，是以m个样本平均向量为基础，估计真实n个平均向量的类

型。计算个体到每个类中心的Mahalanobis距离判断每个个体更趋于归属哪个类。2．典型判别函数

典型判别函数（canonical

discriminantfunctions）思想来源于1936年Fisher的研究它建立在以方差分析为基础的F-统计量上，找出p个变量的线性组合，使类间差异最大，类内差异最小。找出线性组合实际是一个特征值问题，与类间的线性组合无关。3．Fisher判别分析(P359-361)4．Bayesian判别分析

Bayesian判别的统计思想是，假定对研究对象已有一定认识，常用先验概率分布来描述这种

认识；然后我们取得一个样本，用样本修正已

有知识（先验概率分布），得到后验概率分布，各种统计推断都通过后验概率分布来进行。将

这种思想用于判别分析，就得到Bayesian判别。5．逐步判别分析MATLAB中判别分析的实现

判别分析应用训练数据，估计自变量的判别函数参数。判别函数决定自变量不同类型空间的边界。1．函数classify

MATLAB中，函数classify用于判别分析。其语法规则是：class=classify(sample,

training,

group,

ty[class,

err,

posterior,

logp]=classify(…)

以MATLAB自带的数据fisheriris.mat为例，讨论判别分析过程。

（1）先选择训练数据，采用散点图观察数据类型：从图可见，该数据包括两种类型，而且没有很明显的区分界限。load

fisheririsSL

=

meas(51:end,1);SW

=

meas(51:end,2);group

=

species(51:end);h1

=

gscatter(SL,SW,group);set(h1,"LineWidth",2)legend("组1","组2","Location","NW")

（2）采用函数classify进行判别分类：为了表达分类在同样尺度上，可以采用meshgrid函数绘制格网。[X,Y]

=

meshgrid(linspace(4,8),linspace(X

=

X(:);

Y

=

Y(:);

[class,err,posterior,logp]

=

classify([XSW],group,"quadratic")（3）可视化分类结果（图）：hold

ongscatter(X,Y,class,"wb",".",1,"off")axis([4.5

8

2

4])xlabel("萼片长度（SL）")ylabel("萼片宽度（SW）")当然，也可以直接列举出分类结果。如：load

discrimsample

=

ratings(idx,:)training

=

ratings(1:200,:)g

=

group(1:200)class

=

classify(sample,training,g)first5

=

class(1:5)first5

=

2

2

2

2

2函数knnclassify和函数svmclassify

MATLAB中，函数svmclassify和函数knnclassi也是判别分类函数。其中函数svmclassify是以支持向量机（a

support

vector

machine）为基础的分类器，函数knnclassify是以近邻方法为基础的数据分类。knnclassify函数基本句法是：

class

=

knnclassify(Sample,

Training,

Groudistance,

rule)svmclassify函数基本句法是：

group

=

svmclassify(SVMStruct,

Sample,"PropertyName",

PropertyValue...)以下是一个简单例子：sample

=

[.9

.8;.1

.3;.2

.6]training=[0

0;.5

.5;1

1]group

=

[1;2;3]class

=

knnclassify(sample,

training,

group)class

=

3

1

2更大一些的数据集，我们可以随机产生一组随机正态分布训练数据：training

=

[mvnrnd([

1

1],

eye(2),

100);

mvnrnd([2*eye(2),

100)]指定类型，并绘制训练数据类型分布（图）：group

=

[repmat(1,100,1);

repmat(2,100,1)]gscatter(training(:,1),training

(:,2),group,"rlegend("训练组1","训练组2");hold

on再生成一个连续均匀分布的样本分类（图）：sample

=

unifrnd(-5,

5,

100,

2)c

=

knnclassify(sample,

training,

group)gscatter(sample(:,1),sample(:,2),c,"mc");

ho

legend("训练组1","训练组2","组1中的数据","中的数据")hold

off现在指定k=10，即每个中心有10个近邻（图）：gscatter(training(:,1),training(:,2),group,hold

onc3

=

knnclassify(sample,

training,

group,

10)gscatter(sample(:,1),sample(:,2),c3,"mc","o

legend("训练组1","训练组2","组1中的数据","中的数据")hold

offsvmclassify函数%

Load

the

data

and

select

features

for

classificationload

fisheririsdata

=

[meas(:,1),

meas(:,2)];%

Extract

the

Setosa

classgroups

=

ismember(species,"setosa");%

Randomly

select

training

and

test

sets[train,

test]

=

crossvalind("holdOut",groups);cp

=

classperf(groups);%

Use

a

linear

support

vector

machine

classifiersvmStruct

=

svmtrain(data(train,:),groups(train),"showpclasses

=

svmclassify(svmStruct,data(test,:),"showplot"%

See

how

well

the

classifier

performedclassperf(cp,classes,test);cp.CorrectRate聚类分析聚类分析(cluster

analysis)，也称为分割分析（segmentation

analysis）或分类分析（taxonomyanalysis）,根据数据特征将数据集分割为类型，是一种以对象组间距离最大，组内对象距离最小为基础，将对象分为多个类群的方法。聚类分析是

一种与判别分析功能相反的多元统计分析技术。

本质上，聚类分析属于非监督分类，事先并不知

道包含多少种类型。判别分析则是监督分类，事

先已知类型。明显地，相对判别分析而言，聚类

分析困难一些，因为判别分析事先已知类型情形，聚类分析类型事先是未知的。聚类分析非常有用，它可以用于探索分析多种类型的数据集，完成数

据挖掘，在工程、生物、医疗、语言、人类学、

心理学和市场等许多领域有广泛应用。

聚类分析算法实际是找出一些变量关系的距离函数以样本之间相似性，类群之间相离性为基础，分类样本，完成聚类。因此聚类分析中，需要选择聚类标准，确定将形成多少类。多数聚类算法只是简单地以最大相似性连接样本对象对，聚合成类群。首先计算对象间的相似性，然后以相似性计算结果为基础，生成一个层次树，即一种表示亲缘关系的谱系图（dendrogram）。有时还应当考虑变量的绝对值，绝对值大不一定表示重要。有些环境数值尽管绝对值小，但是其含量对揭示环境状态，却显得非常重要。因此需要采用自动缩放技术（autoscaling）,通过校正将数据变换均值为0，方差是1的数据。其次，为避免负相关，可使用对数比

率转换方法。聚类分析类型（a1g）gl谱om系er聚at类iv方e法hi（erhairecrhaircahlimceatlhcoldu，stAeHrMi）ng。）也：有又人称称凝之聚为层系次统方聚法类。它以对象之间距离计算为基础，合并邻近对象，最终产生表示亲缘关系的谱系图。聚（类2）（分di割vi聚si类ve（hpiaertairtcihoincianlgmceltuhsotde，riDnHgM）：或正k-好均与值前聚者类相。反这，种也方称法为有分许裂多法，基本算法是根据经验选择一些类型中心(比如两个类)，然后以对象到类型中心距离为基础，确认哪些对象该进或该出某些类型。一旦分割完成，按类

型平均值重新计算类型中心。再次合并或分割子类，直到得到满足限定条件

的分类结果。其计算基础是近邻方法（the

nearest-neighbor

method）。的谱类系群聚；类对和异k-常均值敏聚感类；都统有计一特定征的未主知观等性问，题甚。至出现空类；趋于特定分布形状以（概3）率基密于度模为型中聚心类的（有m限od混el合-模ba型se。d

c使lu用st概e率ri模ng型）进方行法聚是类以的概思率想模开型始为于基上础，如世纪六十年代，1965年Edwards等提出借助有限混合模型，进行聚类分析问题。有限混合模型方法是，假定以加权组分密度模拟概率密度函数，当使用有限混合进行聚类分析，聚类问题转换为估计假设混合模型的参数，如概率密度

估计。每个组分密度对应一个类群，采用后验概率确认聚类结果。常用的估

计有限概率密度参数方法是EM算法（Expectation-Maximization

algorithm，

EM）,它以极大似然估计为基础。高斯混合模型聚类（Gaussian

mixture

models）：是通过表示观测变量的概率密度函数为多元正态密度的混合体。聚类过程包括五个步骤择（（1）选模取式变陈量述）：、这特是征聚提类取的（初确始定步观骤测，样包本括数选目择、类选群择数测目量、尺特度征、选完成数据转换）。类（群2）之模间式的邻距近离度或量邻：近多程数度聚，类但方是法不要同求的度距量离观度测量样，本可与能其导可致能不归同属的数据分割结果。性（的3），分观组测：样分本组只是能一是个属将于数或据不分属割于为