Matlab统计工具箱1

Matlab统计工具箱

Getting Started Toolbox functions and GUIs

Organizing Data Data arrays and groups

Descriptive Statistics Data summaries

Statistical Visualization Data patterns and trends

Probability Distributions Modeling data frequency

Random Number Generation Topics in random number generation

Hypothesis Tests Inferences from data

Analysis of Variance Modeling data variance

Regression Analysis Continuous data models

Multivariate Methods Visualization and reduction

Cluster Analysis Identifying data categories

Classification Categorical data models

Markov Models Stochastic data models

Design of Experiments Systematic data collection

Statistical Process Control Production monitoring

方差分析

ANOVA Operations

方差分析模型

方差分析是分析试验（或观测）数据的一种统计方法。在工农业生产和科学研究中，经常要分析各种因素及因素之间的交互作用对研究对象某些指标值的影响。这时需要用到方差分析。利用方差分析，我们能推断哪些因素对所考察指标的影响是显著的，哪些是不显著的。

方差分析包括了：单因素方差分析、双因素方差分析和多因素方差分析。

单因素方差分析

单因素方差分析问题：某一因素（A）对结果（）的影响分析。一般情况下，假设因素A有个不同的水平。单因素方差分析模型如下：

其中表示在第组（水平）第个观察值，为总的（的 ）均值，为第组所的均值且满足，表示第组第个样本的误差，为第组的数据个数且。方差分析的目的是要基于数据分析确定这组数据的均值有没有显著性差异，即假设检验问题：

v.s. 不全相等。 （1）

原假设表明因素A对的值没有影响。

格式

p = anova1(X)
p = anova1(X,group)
p = anova1(X,group,displayopt)
[p,table] = anova1(...)
[p,table,stats] = anova1(...)

实例 1

Create X with columns that are constants plus random normal disturbances with mean zero and standard deviation one:

X = meshgrid(1:5)

X =

1 2 3 4 5

1 2 3 4 5

1 2 3 4 5

1 2 3 4 5

1 2 3 4 5

X = X + normrnd(0,1,5,5)

X =

1.3550 2.0662 2.4688 5.9447 5.4897

2.0693 1.7611 1.4864 4.8826 6.3222

2.1919 0.7276 3.1905 4.8768 4.6841

2.7620 1.8179 3.9506 4.4678 4.9291

-0.3626 1.1685 3.5742 2.1945 5.9465

Perform one-way ANOVA:

p = anova1(X)

p =

7.9370e-006

the very small p value indicates that differences between column means are highly significant. The probability of this outcome under the null hypothesis (that samples drawn from the same population would have means differing by the amounts seen in X) is equal to the p value.

实例 2

The following example is from a study of the strength of structural beams in Hogg. The vector strength measures deflections of beams in thousandths of an inch under 3,000 pounds of force. The vector alloy identifies each beam as steel (铁'st'), alloy 1 (合金1'al1'), or alloy 2 (合金2'al2'). (Although alloy is sorted in this example, grouping variables do not need to be sorted.) The null hypothesis is that steel beams are equal in strength to beams made of the two more expensive alloys.

strength = [82 86 79 83 84 85 86 87...

74 82 78 75 76 77 ...

79 79 77 78 82 79];

alloy = {'st','st','st','st','st','st','st','st',...

'al1','al1','al1','al1','al1','al1',...

'al2','al2','al2','al2','al2','al2'};

p = anova1(strength,alloy)

p =

1.5264e-004

双因素方差分析模型

双因素方差分析问题。一般地，不妨假设因素A有个水平，因素B有个水平。在双因素（多因素方差）分析中，不仅要考虑每个因素的影响，有时还要考虑因素之间的交互作用。下面是考虑有交互作用的双因素方差分析的模型如下：

其中表示A因素第组B因素第组第个样本的观察值，为总的均值，为A因素第组的均值，为B因素第组的均值，表示A因素第组B因素第组交互作用的均值，满足。方差分析的目的是要基于数据分析确定两因素是否存在交互作用，各因素不同组的均值有没有显著性差异，也就是假设检验问题：

v.s. 不全相等， （2）

v.s. 不全相等， （3）

v.s. （4）

格式

p = anova2(X,reps)
p = anova2(X,reps,displayopt)
[p,table] = anova2(...)
[p,table,stats] = anova2(...)

实例

The data below come from a study of popcorn brands and popper type (Hogg 1987). The columns of the matrix popcorn are brands (Gourmet, National, and Generic). The rows are popper type (Oil and Air.) The study popped a batch of each brand three times with each popper. The values are the yield in cups of popped popcorn.

load popcorn

popcorn

popcorn =

5.5000 4.5000 3.5000

5.5000 4.5000 4.0000

6.0000 4.0000 3.0000

6.5000 5.0000 4.0000

7.0000 5.5000 5.0000

7.0000 5.0000 4.5000

p = anova2(popcorn,3)

p =

0.0000 0.0001 0.7462

回归分析

Linear Regression

回归分析模型

模型假设和相关问题

设变量与多个变量(p>1)有线性相关关系，相应的观察数据为，。我们讨论的多元线性回归模型为：

其中和是未知参数。

格式

b = regress(y,X)
[b,bint] = regress(y,X)
[b,bint,r] = regress(y,X)
[b,bint,r,rint] = regress(y,X)
[b,bint,r,rint,stats] = regress(y,X)
[...] = regress(y,X,alpha)

实例

导入数据cars; 确定weight 和 horsepower 作为自变量（预测变量）, mileage 为响应变量:

load carsmall

x1 = Weight;

x2 = Horsepower; % Contains NaN data

y = MPG;

计算有交互作用项的线性回归模型的回归系数:

X = [ones(size(x1)) x1 x2 x1.*x2];

b = regress(y,X) % Removes NaN data

b =

60.7104

-0.0102

-0.1882

0.0000

Plot the data and the model:

scatter3(x1,x2,y,'filled')

hold on

x1fit = min(x1):100:max(x1);

x2fit = min(x2):10:max(x2);

[X1FIT,X2FIT] = meshgrid(x1fit,x2fit);

YFIT = b(1) + b(2)*X1FIT + b(3)*X2FIT + b(4)*X1FIT.*X2FIT;

mesh(X1FIT,X2FIT,YFIT)

xlabel('Weight')

ylabel('Horsepower')

zlabel('MPG')

view(50,10)

主成分分析

Principal Component Analysis

模型

在研究某个问题时涉及到个变量（指标）（记为列向量，且有均值向量和方差协方差矩阵）。通常较大且这些指标之间存在一定程度的相关性。因此我们考虑用以下线性变换将这个变量转换成新的个两两不相关的变量：

记，且限定，则上述的线性变换中的变换矩阵的转置。由于这组新的变量中任意两个不相关，即对任意有。而且变量中每一个变量的方差依次是余下所有变量中最大的，因此有。方差越大，包含的信息就越多。我们称为的第主成分。理论上可以证明第个主成分对应于的方差协方差矩阵的第个特征值，与特征值相应的标准化特征向量记为。且有的特征根是主成分的方差，即，这里有。因此主成分的方差协方差矩阵为对角矩阵diag()。

格式

[COEFF,SCORE] = princomp(X)
[COEFF,SCORE,latent] = princomp(X)
[COEFF,SCORE,latent,tsquare] = princomp(X)
[...] = princomp(X,'econ')

实例

Compute principal components for the ingredients data in the Hald data set, and the variance accounted for by each component.

load hald;

[pc,score,latent,tsquare] = princomp(ingredients);

pc,latent

pc =

0.0678 -0.6460 0.5673 -0.5062

0.6785 -0.0200 -0.5440 -0.4933

-0.0290 0.7553 0.4036 -0.5156

-0.7309 -0.1085 -0.4684 -0.4844

latent =

517.7969

67.4964

12.4054

0.2372

The following command and plot show that two components account for 98% of the variance:

cumsum(latent)./sum(latent)

ans =

0.86597

0.97886

0.9996

1

biplot(pc(:,1:2),'Scores',score(:,1:2),'VarLabels',...

{'X1' 'X2' 'X3' 'X4'})

因子分析

Factor Analysis

模型假设

类似于在主成分分析中的情形，在研究某个问题时涉及到个容易观测的变量（记为列向量）。其均值向量，方差协方差矩阵）。我们假定存在个潜在相互独立的的公共因子（变量）（记为列向量）和个相互独立的特殊因子变量（记为列向量）满足以下的模型：

其中，，记为。相互独立。模型用矩阵向量形式可表示为

其中

被称为因子载荷矩阵。

格式

lambda = factoran(X,m)
[lambda,psi] = factoran(X,m)
[lambda,psi,T] = factoran(X,m)
[lambda,psi,T,stats] = factoran(X,m)
[lambda,psi,T,stats,F] = factoran(X,m)
[...] = factoran(...,param1,val1,param2,val2,...)

实例1

导入carbig数据, 用两因子拟合默认的模型.

load carbig

X = [Acceleration Displacement Horsepower MPG Weight];

X = X(all(~isnan(X),2),:);

[Lambda,Psi,T,stats,F] = factoran(X,2,...

'scores','regression');

inv(T'*T) % Estimated correlation matrix of F, == eye(2)

Lambda*Lambda'+diag(Psi) % Estimated correlation matrix

Lambda*inv(T) % Unrotate the loadings

F*T' % Unrotate the factor scores

biplot(Lambda,... % Create biplot of two factors

'LineWidth',2,...

'MarkerSize',20)

实例 2

Although the estimates are the same, the use of a covariance matrix rather than raw data doesn't let you request scores or significance level:

[Lambda,Psi,T] = factoran(cov(X),2,'xtype','cov')

[Lambda,Psi,T] = factoran(corrcoef(X),2,'xtype','cov')

实例 3

斜交旋转Use promax rotation:

[Lambda,Psi,T,stats,F] = factoran(X,2,'rotate','promax',...

'powerpm',4);

inv(T'*T) % Est'd corr of F,

% no longer eye(2)

Lambda*inv(T'*T)*Lambda'+diag(Psi) % Est'd corr of X

Plot the unrotated variables with oblique axes superimposed.

invT = inv(T)

Lambda0 = Lambda*invT

line([-invT(1,1) invT(1,1) NaN -invT(2,1) invT(2,1)], ...

[-invT(1,2) invT(1,2) NaN -invT(2,2) invT(2,2)], ...

'Color','r','linewidth',2)

hold on

biplot(Lambda0,...

'LineWidth',2,...

'MarkerSize',20)

xlabel('Loadings for unrotated Factor 1')

ylabel('Loadings for unrotated Factor 2')

聚类分析

Hierarchical Clustering

系统聚类法是聚类分析中应用最为广泛的一种方法.它的基本原理是：首先将一定数量的样品(或指标)各自看成一类，然后根据样品(或指标)的亲疏程度，将亲疏程度最高的两类合并，然后重复进行，直到所有的样品都合成一类.衡量亲疏程度的指标有两类：距离、相似系数.

1.常用距离

1）欧氏距离

假设有两个维样本和，则它们的欧氏距离为

2）标准化欧氏距离

假设有两个维样本和，则它们的标准化欧氏距离为

其中：表示个样本的方差矩阵，，表示第列的方差.

3）马氏距离

假设共有个指标，第个指标共测得个数据(要求)：

，

于是，我们得到阶的数据矩阵，每一行是一个样本数据.阶数据矩阵的阶协方差矩阵记做.

两个维样本和的马氏距离如下：

马氏距离考虑了各个指标量纲的标准化，是对其它几种距离的改进.马氏距离不仅排除了量纲的影响，而且合理考虑了指标的相关性.

4）布洛克距离

两个维样本和的布洛克距离如下：

5）闵可夫斯基距离

两个维样本和的闵可夫斯基距离如下：

注：时是布洛克距离；时是欧氏距离.

6）余弦距离

这是受相似性几何原理启发而产生的一种标准，在识别图像和文字时，常用夹角余弦为标准.

7）相似距离

2.MATLAB中常用的计算距离的函数

假设我们有阶数据矩阵，每一行是一个样本数据. 在MATLAB中计算样本点之间距离的内部函数为

y=pdist(x) 计算样本点之间的欧氏距离

y=pdist(x,'seuclid') 计算样本点之间的标准化欧氏距离

y=pdist(x,'mahal') 计算样本点之间的马氏距离

y=pdist(x,'cityblock') 计算样本点之间的布洛克距离

y=pdist(x,'minkowski') 计算样本点之间的闵可夫斯基距离

y=pdist(x,'minkowski',p) 计算样本点之间的参数为p的闵可夫斯基距离

y=pdist(x,'cosine') 计算样本点之间的余弦距离

y=pdist(x,'correlation') 计算样本点之间的相似距离

另外，内部函数yy=squareform(y)表示将样本点之间的距离用矩阵的形式输出.

3.常用的聚类方法

常用的聚类方法主要有以下几种：最短距离法、最长距离法、中间距离法、重心法、平方和递增法等等.

4. 创建系统聚类树

假设已经得到样本点之间的距离y，可以用linkage函数创建系统聚类树，格式为z=linkage(y).

其中：z为一个包含聚类树信息的(m-1)×3的矩阵.例如：

z=

2.000 5.000 0.2

3.000 4.000 1.28

则z的第一行表示第2、第5样本点连接为一个类，它们距离为0.2；则z的第二行表示第3、第4样本点连接为一个类，它们距离为1.28.

在MATLAB中创建系统聚类树的函数为

z=linkage(y) 表示用最短距离法创建系统聚类树

z=linkage(y,'complete') 表示用最长距离法创建系统聚类树

z=linkage(y,'average') 表示用平均距离法创建系统聚类树

z=linkage(y,'centroid') 表示用重心距离法创建系统聚类树

z=linkage(y,'ward') 表示用平方和递增法创建系统聚类树

格式

T = cluster(Z,'cutoff',c)
T = cluster(Z,'cutoff',c,'depth',d)
T = cluster(Z,'cutoff',c,'criterion',criterion)
T = cluster(Z,'maxclust',n)

实例

比较Fisher iris数据的聚类结果和类别数据:

load fisheriris

d = pdist(meas);

Z = linkage(d);

c = cluster(Z,'maxclust',3:5);

crosstab(c(:,1),species)

ans =

0 0 2

0 50 48

50 0 0

crosstab(c(:,2),species)

ans =

0 0 1

0 50 47

0 0 2

50 0 0

crosstab(c(:,3),species)

ans =

0 4 0

0 46 47

0 0 1

0 0 2

50 0 0

判别分析

Discriminant Analysis

判别分析概说

在已知研究对象分成若干类型，并已取得各种类型的一批已知样品的观测数据，在此基础上根据某些准则建立判别式，然后对未知类型的样品进行判别分类。

距离判别法—首先根据已知分类的数据，分别计算各类的重心，计算新个体到每类的距离，确定最短的距离（欧氏距离、马氏距离）

Fisher判别法—利用已知类别个体的指标构造判别式（同类差别较小、不同类差别较大），按照判别式的值判断新个体的类别

Bayes判别法—计算新给样品属于各总体的条件概率，比较概率的大小，然后将新样品判归为来自概率最大的总体

判别分析是利用原有的分类信息，得到体现这种分类的函数关系式（称之为判别函数，一般是与分类相关的若干个指标的线性关系式），然后利用该函数去判断未知样品属于哪一类。

Matlab中对于给定的数据，用classify函数进行线性判别分析，用mahal函数计算马氏距离。

格式

class = classify(sample,training,group)
class = classify(sample,training,group,'type')
class = classify(sample,training,group,'type',prior)
[class,err] = classify(...)
[class,err,POSTERIOR] = classify(...)
[class,err,POSTERIOR,logp] = classify(...)
[class,err,POSTERIOR,logp,coeff] = classify(...)

实例

load fisheriris

SL = meas(51:end,1);

SW = meas(51:end,2);

group = species(51:end);

h1 = gscatter(SL,SW,group,'rb','v^',[],'off');

set(h1,'LineWidth',2)

legend('Fisher versicolor','Fisher virginica',...

'Location','NW')

[X,Y] = meshgrid(linspace(4.5,8),linspace(2,4));

X = X(:); Y = Y(:);

[C,err,P,logp,coeff] = classify([X Y],[SL SW],...

group,'Quadratic');

Visualize the classification:

hold on;

gscatter(X,Y,C,'rb','.',1,'off');

K = coeff(1,2).const;

L = coeff(1,2).linear;

Q = coeff(1,2).quadratic;

% Function to compute K + L*v + v'*Q*v for multiple vectors

% v=[x;y]. Accepts x and y as scalars or column vectors.

f = @(x,y) K + [x y]*L + sum(([x y]*Q) .* [x y], 2);

h2 = ezplot(f,[4.5 8 2 4]);

set(h2,'Color','m','LineWidth',2)

axis([4.5 8 2 4])

xlabel('Sepal Length')

ylabel('Sepal Width')

title('{\bf Classification with Fisher Training Data}')

本文来源：https://www.2haoxitong.net/k/doc/59c58fff58fafab068dc0269.html