精通MATLAB科学计算（第3版）(王正林)13-3r

第章概率统计计算

概率论与数理统计是一门研究随机现象并找出其统计规律的学科，广泛应用于社会、

经济、科学等各个领域。MATLAB中的统计工具箱(StatisticsToolbox)包含了200多个

用于概率统计方面的功能函数，且具有简单的接口操作。

通过本章学习，读者能利用统计工具箱中相关函数进行统计量计算、参数估计、假设

检验、方差分析和回归分析等计算，而且还能使用统计学中的一些特殊的绘图工具，直观

地认识样本的结构和分布特征。

13.1MATLAB统计工具箱介绍

MATLAB中的统计工具箱是一套建立在MATLAB数值计算环境下的统计分析工具，

能够支持范围广泛的统计计算任务。

统计工具箱主要包括如下两类工具。

(1)概率统计的数值计算函数

该函数可以通过命令行或用户的应用程序调用，多数为MATLAB中的M文件。

(2)交互式图形工具函数

该函数是能够通过图形用户界面(GUI)来使用的交互式图形工具。

表13-1给出了统计工具箱中函数的主要分类。

表13-1统计工具箱中函数的主要分类

函数分类具体描述

参数估计

概率分布累积分布函数

概率密度函数

逆累积分布函数

概率分布随机数产生器

分配函数的矩

实验设计构建统计模型

主成分分析对数据进行降维处理

假设检验假设的统计检验

描述统计学针对数据抽样的描述统计

第13章概率统计计算

统计作图统计的作图

统计的过程控制统计的过程控制

续表

函数分类具体描述

聚集分析有相似的特性的项目列入群体中

线性模型对数据进行线性模型拟合

非线性回归非线性回归模型的拟合

文件输入，输出对操作系统文件读写数据

数据示范的数据

13.2随机变量的数字特征

随机变量的数字特征是与随机变量有关的某些数值，它们虽然不能完整地描述随机变

量，但能描述随机变量在某些方面的重要特征。本节主要讲述用MATLAB求解随机变量

的期望(均值X方差、矩、协方差矩阵、相关系数等数字特征。该部分所用到的命令函

数属于统计工具的描述统计学(DescriptiveStatistics)部分。

13.2.1期望

设离散性随机变量的分布律为：

?{X=x*}=p*,&=1,2,…

则定义X的期望为：

co

E(X)=£xkPk

k=l

而对于来自总体X的一个样本，设其样本值为X=(X1,X2,…,X"),则定义样本均值为:

可以证明F依概率收敛于X的均值。因此通常用亍来近似代替X的均值。MATLAB

中提供了函数mean。来计算这种较为简单的均值，其调用格式如下。

1.m=mean(X)

若X为向量形式，则mean(X)为向量中各元素的平均值；

若X为矩阵，则mean(X)为一个行向量，其中的每一个元素是矩阵X各列元素的平均

值。

◄◄◄◄307

精通MATLAB科学计算(第2忡-----------------------------

2.m=mean(X,dim)

函数返回X在某一维上的均值，其中维数由标量dim设置。

例如X为一个矩阵时，mean(X,2)是一个列向量,包含矩阵X每一行的均值。dim的

默认值为lo

【例13-1]样本均值计算实例lo下面给出一组样本值来自于某一个未知分布，计

算该分布的样本均值。

[19・1,13・5,12.6,10,14.6,24,13.7,15,20,13.5]

解：通过在MATLAB中编写程序exl301.m来实现。

%exl301.m用mean()计算样本均值

x=[19.1,13.5,12.6,10,14.6,24,13.7,15,20,13.5];

m=mean(x)

命令窗口中输出结果为：

>>exl301

m=

15.6000

由以上计算结果可知，该分布的样本均值为15.6。

【例13・2】样本均值计算实例2O如下为来自于总体的三个样本的样本值，分别计

算其均值。

天：2.13.55.01.53.03.24.22.53.01.8

X2：3.23.04.54.21.41.55.04.53.53.0

X3：2.42.64.23.53.23.34.04.23.62.8

解：用矩阵的均值形式计算，在MATLAB中编写程序exl302.m来实现。

%exl302.m分别计算三个样本的均值

Xl=[2.13.55.01.53.03.24.22.53.01.8];备第一个样本

X2=[3.23.04.54.21.41.55.04.53.53.0];%第二个样本

X3=[2.42.64.23.53.23.34.04.23.62.8];%第三个样本

X=[X1;X2;X3];务将三个样本作为X的各行

ml=mean(X,2)%计算X各行的均值，即得各样本均值

,,8将三个样本作为X的各列

X=[X1'ZX2ZX3]；

m2=mean(X)省得各列均值

命令窗口输出为：

>exl302

ml=

2.9800

308►►►►

---------------------------------------------------------第章概率统计计算

3.3800

3.3800

m2=

2.98003.38003.3800

由以上计算结果可知，样本XI的均值为2.98,样本A2的均值为3.38,样本X3的均

值为3.38。

【例13-3】样本均值计算实例3。设随机变量X、丫的联合分布律如表13-2所示，

分别计算E(X),E(Y),E(R3)以及E(X+Y)O

表13-2x、丫联合分布律

1234p(y=j)

01/100001/10

104/102/101/107/10

20002/102/10

0g)1/104/102/103/101

解：用均值定义计算，通过在MATLAB中编写程序exl303.m来实现，具体代码如下

所示。

%exl303.m分别计算二维离散分布的均值

X=[0,l,2];

Y=[l,2,3,4];

Px=[0.1,0.7,0.2];%x各点对应的概率

Py=[0.1,0.4,0.2z0.3];%Y各点对应概率

mX=sum(X.*Px)%E(X)sum为求和函数。

mY=sum(Y.*Py)

z=X+3;

mX3=sum(z.*Px)%E(X+3)

mXY=mX+mY%E(X+Y)

命令窗口输出为：

>>exl303

mX=

1.1000%E(x)

mY=

2.7000%E(Y)

mX3=

4.1000%E(X+3)

mXY=

3.8000%E(X+Y)

在MATLAB中还给出了求其他形式均值的一些函数，其命令形式与mean相同。各

命令具体含义详见表13-3。

◄◄◄◄309

精通MATLAB科学计算（第2蚪

表13-3几种均值函数

函数名功能

geomean求样本几何均值

harmmcan求样本调和平均值

nanmean可忽略样本中的非数值输入，并计算相应均值

trimmean除去输入样本中一些数值变化太大的值，然后再计算样本均值

如需了解表13-3中函数的详细用法，可到MATLAB的Help中的相关部分检索查询

函数的详细使用说明。

310►►►►

第13章概率统计计算

13.2.2方差、标准差、矩

方差是用来刻画随机变量x取值分散程度的一个量。其一般用下式表达：

D(X)=E{[x-E(x)]1}

在应用上还引入与随机变量X具有相同量纲的量向方，记为b(X)，称为标准差或

均方差。

X的左阶中心矩应为：

E{[X-E(X)]k],k=2,3,-

可知方差即为二阶中心矩。

对于一个样本来说，样本方差通常有两种形式表达，如下所示。

无偏估计式：

有偏估计式：

)2

〃/=!

样本标准差也对应有如下两种形式：

s=*3a

或

S=g=I眩(XL42

V«,=1

样本的Z阶中心矩为：

MATLAB中用方差函数var()来计算样本方差，用标准差函数std()来计算样本的标准

差，用矩函数moment。计算样本的各阶中心矩，具体用法如下。

1.var()方差函数

k=var(X)若X为向量则返回向量的样本方差值，若X为矩阵则返回矩阵各列向

◄◄◄311

精通MATLAB科学计算(第2忡-----------------------------

量方差组成的行向量。其采用无偏式计算方差。

V=var(X,1)函数采用有偏估计式计算X的方差，即前置因子为var(X,O)等同

于var(X),其采用无偏式计算方差，前置因子为

K=var(X,w)函数返回X以w，为权的方差。对于矩阵X,w的元素个数必须等于X

的行数；对于向量X,w的元素个数与X元素个数相同。

V=var(X,flag,dim)函数返回X在特定维上的方差，dim指定维数，flag指定选择的

计算式。即flag=O,选择无偏式计算；flag=l,选择有偏式计算。

此外，nanvar()是忽略非数(NaNs)的方差计算函数，调用格式与var()相同。

2.std()标准差函数

s=std(X)函数返回向量(矩阵)X的标准差(前置因子\/n-\卜

s=std(X,flag)flag=O,前置因子为\!n-\;flag=l,前置因子为l/n0

s=std(X,flag,dim)函数返回X在特定维上的标准差，dim指定维数，flag指定选

择的计算式。

nanstd()函数用法同std(),忽略非数的标准差计算函数。

【例13-4】样本方差、标准差计算实例。分别计算向量*=[1,-1,1,2]的方差和标准

差以及加权后的方差和标准差，其中加权向量为»v=[l,2,3,4]o(加权方差和标准差值指的

是计算样本方差或标准差时，各样本值要乘以相应的权值，然后再除以总的权值。)

解：通过在MATLAB中编写程序exl304.m来实现，具体代码如下所示。

%exl304.m计算x的方差和标准差

x=[-l-112];3输入向量X

w=[l234];告权向量W

当用各种命令形式计算方差

vl=var(x)

v2=var(x,0)

v3=var(x,1)

v4=var(x,w)

v5=var(x,1,w)

多用各种命令形式计算向量标准差

sl=std(x)

s2=std(x,1)

s3=std(x,w)

在命令窗口中运行该文件，得如下结果：

>>exl304

vl=

2.2500%var(X)

v2=

312►►►>

第13章概率统计计算

2.2500%var(X,0)

v3=

1.6875%var(X,1)

v4=

1.5600%var(X,w)

v5=

0000%var(X,1,w)

si=

1.5000%std(X)

s2=

1.2990%std(Xz1)

s3=

1.2490%std.(X,w)

3.moment。矩函数

〃7=moment(X,order)函数返回向量(矩阵)X的左阶中心矩。order规定中心矩的阶数。

???=moment(Xorder,dim)函数返回dim维上的X的中心矩。

【例13-5】中心矩计算实例。随机产生正态分布的一组随机数，计算其3阶中心矩。

解：通过在MATLAB中编写程序exl305.m来实现。

%exl305.m计算正态随机数的3阶矩

X=randn([54])%产生正态分布的随机数矩阵，行数为5,列数为4

ml=moment(X,3)%计算矩阵X各列的3阶矩

m2=moment(X,3,2)%计算矩阵X各行的3阶矩，并返回2维上的中心矩

»exl305

X=

-0.43261.1909-0.18670.1139

-1.66561.18920.72581.0668

0.1253-0.0376-0.58830.0593

0.28770.32732.1832-0.0956

-1.14650.1746-0.1364-0.8323

ml=

-0.11350.03620.87800.0586

m2=

0.1983

-1.7088

-0.0227

0.7167

-0.0004

13.2.3协方差、相关系数

随机变量X、y的协方差和相关系数的定义式如下：

◄◄◄313

精通MATLAB科学计算(第2忡-----------------------------

cov(x,y)=£{[x-E(x)][y-E(y)]}

、cov(x,y)

cocft(x,y)=iJ

对于n维随机变量，通常用协方差矩阵描述它的2阶中心矩。如对于二维随机变量

(XD,定义协方差矩阵形式为：

(C\\C12]

。22J

其中：

2

cix=E{[x-E(x)]]

cl2=E{[x-E(x)][y-E(y)]}

c2i=E{[y-E(y)][x-E(x)]}

2

c22=E{[y-E(y)]]

其相应的样本协方差形式此处不再给出。与样本方差形式类似，也具有两种形式。

MATLAB中用cov()和corrcoef。计算协方差和相关系数。具体用法如下。

1.COV。计算协方差

C=cov(X)X为向量时，函数返回此向量的方差。X为矩阵时，矩阵的每一行表示一

组观察值，每一列代表一个变量。函数返回此矩阵的协方差矩阵，其中协方差矩阵的对角

元素是X矩阵的列向量的方差值。

c=cov(x,F)返回x、y的协方差矩阵，其中x、y行数和列数相同。

C=cov(X,l),C=cov(X,y,l)计算协方差矩阵时前置系数取1/rtocov(X,0)与cov(X)相

同，都是取前置系数为1/n-l,此用法可参考var函数。

2.corrcoe%)计算相关系数

R=corrcoef(X)返回矩阵X的相关系数矩阵，其各点值对应于相关矩阵的各点值除以

相应的标准差。

7?=corrcoef(x,y)返回x、y的相关系数矩阵。若x、y分别为列向量，则该命令等

同于7?=corrcoef([,ry])o

[7?,P]=concoef(...)返回的P矩阵是不相关假设检验的p值。

[R,P,RLO,RUP]=corrco砥…)对于每一个R值返回的95%置信区间为[RLO,RUP]。

[...]=corrcoef(...,'param1',vail,'param2',val2,...)设定特殊的可选项。

corrcoef函数功能强大，本节仅讨论它的前两种常用用法。

【例13-6】协方差和相关系数计算实例。计算向量x和向量y的协方差和相关系数。

314►►►►

第13章概率统计计算

向量X与向量J，分别如下所示：

x=[1,222,1,2]

y=[2,3,2,2,3,2]

解：通过在MATLAB中编写程序cxl306.m来实现，具体代码如下所示。

%exl306.m计算向量之间的协方差和相关系数

x=[1,2,2,2,1,2];

y=[2,3,2,2,3,2];

cx=cov(x)金计算向量X的协方差，可知其等于X方差

vx=var(x)3x方差

cxy=cov(x,y)用x、y协方差

cor=corrcoef(x,y)%x.y相关系数

>>exl306号运行exl306

ex=

0.2667

vx=

0.2667

exy=

0.2667-0.0667

-0.06670.2667

cor=

1.0000-0.2500

-0.25001.0000

【例13-7】随机变量数字特征综合计算实例。根据表13-4中所示的数据计算各变

量的均值、方差以及它们之间的协方差矩阵和相关系数。

表134几组变量数据

VOLHPMPGSP

894965.496

925555.997

927049.0105

925346.596

897046.2105

925545.497

948043.4107

897339.3103

896639.6100

917838.9106

解：通过在MATLAB中编写程序exl307.m来实现，具体代码如下所示。

%exl307.m计算各变量均值、方差以及它们之间的协方差矩阵和相关系数矩阵

◄◄◄315

精通MATLAB科学计算(第2版i------------------------

VOL=[89929292899294898991]，

HP=[49557053705580736678]1;

MPG=[65.455.949.046.546.245.443.439.339.638.9]1;

SP=[96971059610597107103100106]1;

%将变量数据按列存到一个矩阵中

X=[VOLZHPZMPG,SP];

m=mean(X)%计算X均值，得各变量样本的均值

v=var(X)3计算X方差，得各个变量样本值之间方差

cX=cov(X)%计算变量之间的相关矩阵

cv=diag(cX)%取相关矩阵对角元素

corX=corrcoef(X)S计算变量之间的相关系数

输出计算结果为：

>>exl307

m=

90.900064.900046.9600101.2000

v=

3.2111123.211168.647119.9556

cX=

3.21112.4333-0.78221.0222

2.4333123.2111-64.860047.8000

-0.7822-64.860068.6471-20.4356

1.022247.8000-20.435619.9556

cv=

3.2111

123.2111

68.6471

19.9556

corX=

1.00000.1223-0.05270.1277

0.12231.0000-0.70520.9640

-0.0527-0.70521.0000-0.5521

0.12770.9640-0.55211.0000

13.2.4偏斜度和峰度

偏斜度和峰度是用来描述随机变量分布的形状与对称形式，或正态分布型的偏离程度

的统计量，其分别定义如下。

随机变量工的偏斜度：

E(x-E(x))3

'一卬(x)产2

随机变量X的峰度：

,£(x-£(x))4

卬(x)「

在MATLAB中分别用函数skewness。、kurtosis。来实现上述计算。

316►►►>

第13章概率统计计算

1.skewness函数

调用格式如下：

y=skewness(X),X为向量时，函数返回此向量的倾斜度；X为矩阵时，矩阵的每一

行表示一组观察值，每一列代表一个变量，函数返回此矩阵的每一列的倾斜度。

y=skewness(X,flag),flag用于设定是否要修正偏差。若flag=O，则计算倾斜度时修正

系统偏差；flag=l(默认值)时不修正偏差。

y=skewness(X,flag,dim),在矩阵X的dim维上计算倾斜度。通常对于二维矩阵，

dim=l,函数计算各列的倾斜度；dim=2,函数计算各行的倾斜度。

2.kurtosis函数

调用格式如下：

k=kurtosis(X),函数返回X的峰度。

k=kurtosis(X,flag),flag意义同skewness。

k=kurtosis(X,flag,dim),参数意义同skewness©

通过上面各式可以计算出向量(矩阵)X的偏斜度和峰度。

式中flag默认值为I表示不进行偏移修正;flag可以赋值为0，表示进行偏移修正。dim

指定在矩阵X的特定维上求相应统计量。

下面是一段在MATLAB命令窗口中输入的计算偏斜度和峰度的例子：

»x=[20.1,12.4,22,2,18.0,24.0,25.2,20.0,21.2];%输入向量x

»y=[10.2,9.4,8.5;12.1,10.3,11.0;12.1,13.4,11.4];毛输入矩阵y

>>skewness(x)%求x偏斜度

ans=

-1.4381

>>skewness(xz1)%flag=1

ans=

-1.4381

>>skewness(xz0)%flag=0

ans=

-1.7433

>>kurtosis(x)*求x峰度

ans=

4.0773

>>skewness(y)金求y偏斜度

ans=

-0.70710.5636-0.6559

>>skewness(yz1,2)%flag=l,dim=2

-0.0719

0.2641

0.3480

◄◄◄317

精通MATLAB科学计算(第2版i---------------------------------------------

»kurtosis(y)超求y峰度，flag默认1,dim默认1

ans=

1.50001.50001.5000

13.2.5其他数字特征

除了上面介绍的几类重要的数字特征外，MATLAB的统计工具箱中还有很多对数据

性质进行描述的函数，下面通过一些表格对这些函数进行简单介绍。表13-5所示给出数据

比较和简单计算类函数；表13-6给出了几种特殊的统计学函数。

表13-5数据比较和简单计算类函数

函数名功能

max求样本中的最大值元素

nanmax忽略样本中的非数求最大值元素

tnin求样本中的最小值元素

nanmin忽略样本中的非数求最小值元素

median求随机变量的中值

续表

函数名功能

nanmedian忽略样本中的非数求中值

sum求样本元素的累和

nansum忽略样本中非数求和

cumtrapz梯形累和函数

cumsum求此样本值之前的元素和

mad计算样本绝对偏差的均值或中值，由参数设定

sort将样本值按从小到大排列

sortrows将样本的行作为一个整体进行排序

range求样本值中最大值和最小值之差

prctile求样本的在不同分位数上的值

iqr求四分位差值，即样本在75%分位点与25%分位点上取值之差

表13-6一些特殊函数

函数名功能

jackknife计算jackknife样本的特定统计量

bootstrp通过对数据进行bootstrap重采样计算相应的统计量

bootci计算bootstrap统计量的置信区间

tabulate以矩阵形式给出样本值的频数表

13.3特殊分布的概率计算

318►►►>

---------------------------------------------------------第章概率统计计算

概率统计类课程中曾给出一些具有特殊分布的随机变量，例如均匀分布、二项分布、

正态分布和泊松分布等。MATLAB中也有专门的对这些特殊分布进行操作处理的函数。

MATLABT的统计工具箱支持20多种特殊概率分布。对于每种分布，该工具箱都提

供相应的用于计算其常用统计量的5种基本函数。它们分别是：概率密度函数、累积分布

函数、逆累积函数、均值和方差函数以及随机数发生器，下面分别进行讲述。

13.3.1概率密度函数

表13-7给出了MATLAB-7支持的20多种分布以及它们名称的字母缩写。MATLAB

可以识别这些字母所代表的分布，并进行相应的计算和操作。

表13-7MATLAB支持的分布类型

分布名称字母缩写分布名称字母缩写

Beta分布beta负二项分布nbin

二项分布bino非中心F分布ncf

,2分布chi2非中心/分布net

指数分布exp非中心/分布ncx2

续表

分布名称字母缩写分布名称字母缩写

极值分布Ev正态分布norm

产分布f或F泊松分布poiss

7分布gam瑞利分布rayl

几何分布geo，分布t

广义极值分布gev均匀分布unif

广义帕累托分布gP离散均匀分布unid

超几何分布hyge韦伯分布wbl

对数正态分布logn

MATLAB中求概率密度有一个通用的函数pdf()，通过此函数可以求表I3-7中所示的

23种输入分布的概率密度函数，此函数的调用格式如下：

Y=pdf('name',X,A\,N2,A3)

其中，name为表13-7中分布的字母缩写，X为样本矩阵，、42和N3是分布参数

矩阵，y为概率密度矩阵。

对于某些分布，有些参数矩阵可以不必输入。X、41、N2和43必须是具有同样的大

小的矩阵，当其中之一输入为标量时，程序会将其调整为与其他输入同维的矩阵。

另外，MATLAB对于每一种分布还有一个专用的求其概率密度的函数。它们的基本

用法和含义见表13-8。

Y◄◄◄319

精通MATLAB科学计算（第2回

表13-8概率密度函数

分布名称概率密度函数常用调用格式

Beta分布betapdfr=betapdft%J,fi)

二项分布binopdfy=binopdf(X,N,P)

/分布chi2pdfy=chi2pdRX,F)

混合分布copulapdfy=copulapdf('Gaussian\U,rho)：y=copulapdf('z',(7,rho,Nu)K=copulapdfl；family,t/,alpha)

极值分布cvpdfy=cvpdRX,mu,sigma)

指数分布exppdfy=exppdf(Xmu)

产分布fpdfy^fpdV,H,r2)

7分布gampdfy=gampdRX/l,8)

几何分布gcopdfQgcopdRXP)

广义极值分布gevpdf片gevpdRXK,sigma,mu)

广义帕累托分布gppdfP=gppdf{X,K,sigma,theta)

超几何分布hygcpdfy=hygcpdf(X,M,K,M

对数正态分布lognpdfK^lognpdflA^mu,sigma)

多项式分布mnpdfr=mnpdf(XPROB)

多元正态分布

mvnpdfy=mvnpd*X);H=mvnpdf(X1mu)：K=mvnpdRX,mu,sigma)

多元/分布mvtpdfY=mvtpdfU,C,助

续表

分布名称概率密度函数常用调用格式

负二项分布nbinpdf修nbinpdRXKF)

非中心户分布ncfjxify=ncmd”,NU1,NU2,DELTA)

非中心，分布nctpdfy=nctpdf(A；KDELTA)

非中心/分布ncx2pdfK=ncx2pdf(XV,DELTA)

正态分布normpdf^nonnpdffXmu,sigma)

泊松分布poisspdfr=poisspdf(XLAMBDA)

瑞利分布raylpdfr=raylpd«XS)

学生（分布tpdfy=tpdfIXO

离散均匀分布unidpdfN=unidpdfWV)

连续均匀分布unifpdf^unifpdVJ^)

韦伯分布wblpdf外wblpdRXS/)

上述函数的输出均为对应分布在X处的概率密度。后面的几个输入参数是描述分布的

参数矩阵。需要注意的是，输入的参数的维数必须相等，否则会导致错误结果。

【例13-8］概率密度计算实例。计算坐标［0,2］上的均匀分布在0.4,0.6,1,4,1.6

处的概率密度值。

解：可以分别用通用函数pdf和专用函数unifjxlf求解。在MATLAB中编写程序

exl308.m来实现，具体代码如下所示。

320►►►►

第13章概率统计计算

%exl308.m求[0,2]均匀分布在某点的概率密度值

clearall

x=[0.4,0.6,1.4,1.6];

fl=pdf(»unif'zxz0,2)%用pdf()计算

f2=unifpdf(xz0z2)unifpdf()计算

结果为：

>>exl308

fl=

0.50000.50000.50000.5000

f2=

0.50000.50000.50000.5000

【例13-9]概率密度图绘制实例。作出N(0,l),N(0,4),N(l,4)的概率密度图。

解：通过在MATLAB中编写程序cxl309.m来实现，具体代码如下所示。

%exl309.m绘制正态分布的概率密度图

clearall

x=[-10:0.1:10];

fl=normpdf(xz0,1);%N(0,1)pdf

f2=normpdf(xz0,2);(0,4)pdf

f3=normpdf(x,1,2);%N(1,4)pdf

plot'ro',x,f3,'g+')

legend「X~N(0,1)I3~N(0,4)I”~N(1,4)D

计算结果如图13-1所示。

13.3.2累积与逆累积分布函数

◄◄◄◄321

精通MATLAB科学计算（第2版i-----------------------------------------

1.累积分布函数

若随机变量X的概率密度函数为/（x）,则累积分布函数的定义如下：

X

F（x）=P（X„