数据的统计分析
在日常生活中我们会在很多事件中收集到一些数据(比如：考试分数、窗口排队人数、月用电量、灯泡寿命、测量误差、产品质量、月降雨量等数据)，这些数据的产生一般都是随机的．这些随机数据乍看起来并没有什么规律，但通过数理统计的研究发现：这些随机数还是符合着某种分布规律的，这种规律被称为统计规律．
本实验旨在通过对概率密度函数曲线的直观认识、对数据分布的形态猜测、对某些概率分布的密度函数的参数估计(以正态为例)以及进行简单的正态假设检验，来揭示生活中的随机数据的一些统计规律．
1.概率密度函数pdf系列．以normpdf( )为例，调用格式：
y=normpdf(x, mu,sigma)，
计算参数为mu和sigma的样本数据x的正态概率密度函数．参数sigma必须为正．其中：mu为均值，sigma为标准差．
2.参数估计fit系列．以normfit( )为例，调用格式：
[muhat, sigmahat, muci, sigmaci] = normfit(x, alpha)，
对样本数据x进行参数估计，并计算置信度为100(1－alpha)%的置信区间．如alpha=0.01时，则给出置信度为99％的置信区间．不写明alpha，即表示alpha取0.05．
3．load( )函数．调用格式：
S = load('数据文件')
将纯数据文件(文本文件)中的数据导入Matlab，S是双精度的数组，其行数、列数与数据文件相一致．
4. hist(x, m)函数：画样本数据x的直方图，m为直方图的条数，缺省值为10．
5. tabulate( )函数：绘制频数表．返回table矩阵，第一列包含x的值，第二列包含该值出现次数，最后一列包含每个值的百分比．
6．ttest(x,m,alpha)函数：假设检验函数．此函数对样本数据x进行显著性水平为alpha的t假设检验，以检验正态分布样本x(标准差未知)的均值是否为m．h=1表示拒绝零假设，h=0表示不能拒绝零假设．
7．normplot(x)或weibplot(x)函数：统计绘图函数，进行正态分布检验．
研究表明：如果数据是来自一个正态分布，则该线为一直线形态；如果它是来自其他分布，则为曲线形态．
完全类似地可探索以下一系列函数的用法与作用：
8．累积分布函数cdf系列，如：normcdf( )．
9．逆累积分布函数inv系列，如：norminv( )．
10．随机数发生函数rnd系列，如：normrnd( )．
11．均值与方差函数stat系列，如：normstat( )．
1.常见的概率分布的密度函数及其图形
1)常见概率分布的密度函数(20个，打√的10个将在后面作介绍)
序号
中文函数名
英文函数名
英文简写
备注
1
Beta分布
Beta
beta
2
二项分布
Binomial
bino
√
3
卡方分布
Chisquare
chi2
√抽样
4
指数分布
Exponential
exp
√
5
F分布
F
f
√抽样
6
Gamma分布
Gamma
gam
7
几何分布
Geometric
geo
√
8
超几何分布
Hypergeometric
hyge
9
对数正态分布
Lognormal
logn
10
负二项式分布
Negative Binomial
nbin
11
非中心F分布
Noncentral F
ncf
12
非中心t分布
Noncentral t
nct
13
非中心卡方分布
Noncentral Chi-square
ncx2
14
正态分布
Normal
norm
√
15
泊松分布
Poisson
poiss
√
16
瑞利分布
Rayleigh
rayl
17
T分布
T
t
√抽样
18
均匀分布
Uniform
unif
√
19
离散均匀分布
Discrete Uniform
unid
√
20
Weibull分布
Weibull
weib
2)常见概率分布的密度函数文字说明与图形演示：
(A)常见连续分布的密度函数
(1)正态分布
若连续型随机变量
的密度函数为：
则称为服从正态分布的随机变量，记作
．特别地，称
时的正态分布
为标准正态分布，其概率分布的密度函数参见图1．一个非标准正态分布的密度函数参见图2中的虚线部分(
)．
正态分布是概率论与数理统计中最重要的一个分布，高斯(Gauss)在研究误差理论时首先用正态分布来刻画误差的分布，所以正态分布又称高斯分布．一个变量如果是由大量微小的、独立的随机因素的叠加效果，那么这个变量一定是正态变量．比如测量误差、产品质量、月降雨量等都可用正态分布描述．
x=-8:0.1:8;
y=normpdf(x, 0, 1);
y1=normpdf(x, 1, 2);
plot(x, y, x, y1, ':' );
图1标准正态分布图2标准正态与非标准正态
(2)均匀分布(连续)
若随机变量的密度函数为
则称服从区间
上的均匀分布(连续)，记作
，其概率分布的密度函数见参见图3
．
均匀分布在实际中经常使用，譬如一个半径为
的汽车轮胎，因为轮胎上的任一点接触地面的可能性是相同的，所以轮胎圆周接触地面的位置是服从
上的均匀分布，这只要看一看报废轮胎四周磨损程度几乎是相同的就可明白均匀分布的含义了．
x=-10:0.01:10;r=1;
y=unifpdf(x, 0, 2*pi*r);
plot(x, y);
图3均匀分布(连续)图4指数分布
(3)指数分布
若连续型随机变量的密度函数为：
其中
，
则称为服从参数为
的指数分布的随机变量，记作
．
在实际应用问题中，等待某特定事物发生所需要的时间往往服从指数分布．如某些元件的寿命；某人打一个电话持续的时间；随机服务系统中的服务时间；动物的寿命等都常假定服从指数分布．
指数分布的重要性还在于它是具有无记忆性的连续型随机变量．即：设随机变量服从参数为的指数分布，则对任意的实数
，有
其概率分布的密度函数参见见图4
．
x=0:0.1:30;
y=exppdf(x, 4);
plot(x, y)
(B)常见离散分布的密度函数
(4)几何分布
在一个贝努里实验中，每次试验成功的概率为
，失败的概率为
，设试验进行到第
次才出现成功，则的分布列为：
容易看到
是几何级数
的一般项，于是人们称它为几何分布，其概率分布的密度函数参见图5
．
x=0:30;
y=geopdf(x, 0.5);
plot(x, y)
图5几何分布图6二项分布
(5)二项分布
如果随机变量的分布列为：
则这个分布称为二项分布，记为
．当
时的二项分布又称为0-1分布，分布律为
0
1
一般的二项分布的密度函数参见图6
.
x=0:50;
y=binopdf(x, 500, 0.05);
plot(x, y);
(6)泊松(Poisson)分布
泊松分布是1837年由法国数学家泊松(Poisson S.D.1781-1840)首次提出的，其概率分布列是：
记为
，其概率分布的密度函数参见图7
．
泊松分布是一种常用的离散分布，它与单位时间(或单位面积、单位产品等)上的计数过程相联系，譬如：单位时间内，电话总机接到用户呼唤次数；1平方米内，玻璃上的气泡数；一铸件上的砂眼数；在单位时间内，某种放射性物质分裂到某区域的质点数等等．
x=0:50;
y=poisspdf(x, 25);
plot(x, y);
注：对比二项分布的概率密度函数图可以发现，当二项分布的
与泊松分布充分接近时，两图拟合程度非常高(图6与图7中的
)，直观地验证了泊松定理(泊松分布是二项分布的极限分布)，请对比图6与图7．
图7泊松分布图8均匀分布(离散)
(7)均匀分布(离散)
如果随机变量的分布列为：
则这个分布称为离散均匀分布，记为
，其概率分布的密度函数参见图8
．
n=20;
x=1:n;
y=unidpdf(x, n);
plot(x, y, 'o-' );
(C)三大抽样分布的密度函数
(8)
分布
设随机变量
相互独立，且同服从正态分布，则称随机变量
服从自由度为
的分布，记作
，亦称随机变量
为变量．其概率分布的密度函数参见图9
、图10
，分布的密度函数解析式参见本章的附录表格．
x=0:0.1:20; x=0:0.1:20;
y=chi2pdf(x, 4); y=chi2pdf(x, 10);
plot(x, y); plot(x, y)
图9 分布图10 
分布
(9)
分布
设随机变量
，
，且与
相互独立，则称随机变量
服从自由度为
的分布，记作
，其概率分布的密度函数参见图11，即
，分布的密度函数解析式参见本章的附录表格．
x=0.01:0.1:8.01;
y=fpdf(x, 4, 10);
plot(x, y)
图11 分布图12 
分布
(10)分布
设随机变量
，且与相互独立，则称随机变量
服从于自由度为的分布，记作
，其概率分布的密度函数参见图12，即
．分布的密度函数解析式参见本章的附录表格．
细心的读者可能已经发现，图12的分布图与图1、图2的正态分布十分相似．可以证明：当
时，分布趋于标准正态分布．
x=-6:0.01:6;
y=tpdf(x, 4);
plot(x, y)
2．对给定数据画频数直方图(Histogram)或频数表(Frequency Table)
假定有若干个给定的数据集，它们满足上述10种分布之一，我们现在的任务就是利用画频数直方图等手段，确定它们到底服从哪一类分布．
例1：某一次书面考试的分数罗列如下，试画频数直方图．
鉴于数据的数量较大(包含有120个数据)，可以先在一个文本文件中输入，保存为data1.txt．
75 69 100 80 70 74 78 59 72 73
63 79 69 81 62 87 80 66 86 75
70 85 85 64 78 65 69 67 78 72
60 50 57 83 77 79 78 74 67 83
71 67 71 74 84 74 83 75 73 74
60 91 65 69 80 63 86 67 73 80
74 68 72 80 95 61 77 85 82 71
80 76 83 69 87 76 72 69 66 86
74 87 59 81 88 75 83 71 77 81
88 67 67 76 71 76 79 79 90 62
80 85 81 75 72 57 94 91 83 78
66 74 79 74 82 79 87 76 81 68
x=load('data1.txt');
x=x(:);
hist(x)
结果参见图13．从图形形态上来看，图13较为接近图2所示的正态分布．
图13例1的频数直方图图14例2的频数直方图
例2：某一次上机考试的分数罗列如下(data2.txt，包含有130个数据)，试画频数直方图．
51 70 95 91 70 83 83 96 66 61
79 79 57 85 95 83 63 71 71 72
91 60 69 100 67 87 72 50 60 63
87 98 71 74 96 55 83 67 92 78
56 62 77 79 84 55 59 61 93 56
82 61 88 97 98 95 73 79 81 87
56 92 53 57 93 89 77 89 56 92
99 86 68 57 91 57 81 65 80 99
79 95 79 86 74 56 70 61 72 81
57 75 98 89 69 61 71 77 72 78
70 73 67 59 62 86 84 93 82 80
90 94 84 89 80 67 97 73 80 94
69 64 51 51 92 62 52 86 67 97
x=load('data2.txt');
x=x(:);
hist(x)
结果参见图14．图14看上去很接近图8所示的均匀分布(离散)．
例3：以下给出上海1998年来的月降雨量的数据(data3.txt，包含有98个数据)：
1184.4 1113.4 1203.9 1170.7 975.4 1462.3 947.8
1416.0 709.2 1147.5 935 1016.3 1031.6 1105.7
849.9 1233.4 1008.6 1063.8 1004.9 1086.2 1022.5
1330.9 1439.4 1236.5 1088.1 1288.7 1115.8 1217.5
1320.7 1078.1 1203.4 1480.0 1269.9 1049.2 1318.4
1192.0 1016.0 1508.2 1159.6 1021.3 986.1 794.7
1318.3 1171.2 1161.7 791.2 1143.8 1602.0 951.4
1003.2 840.4 1061.4 958.0 1025.2 1265.0 1196.5
1120.7 1659.3 942.7 1123.3 910.2 1398.5 1208.6
1305.5 1242.3 1572.3 1416.9 1256.1 1285.9 984.8
1390.3 1062.2 1287.3 1477.0 1011.9 1217.7 1197.1
1143.0 1018.8 1243.7 909.3 1030.3 1124.4 811.4
820.9 1184.1 1107.5 991.4 901.7 1176.5 1113.5
1272.9 1200.3 1508.7 772.3 813.0 1392.3 1006.2
x=load('data3.txt');
x=x(:);
hist(x)
结果参见图15．图15看上去很接近图10所示的分布．
图15例3的频数直方图图16例4的频数直方图
在重复数据较多的情况下，我们也可以利用Matlab自带的函数tabulate( )产生频数表，并以频数表的形式来发掘数据分布的规律．
例4：给出以下数据：(data4.txt，含有46个数据)
2 3 6 4 1 5 1 2 3 1 4 2 3 1 3 3 2 3 1 6 4 6 4
6 5 4 3 6 4 3 3 3 3 4 4 5 6 2 1 2 3 4 5 6 5 4
则：
x=load('data4.txt');
x=x(:);
tabulate(x)
hist(x, 6)
Value Count Percent
1 6 13.04%
2 6 13.04%
3 12 26.09%
4 10 21.74%
5 5 10.87%
6 7 15.22%
结果参见图16．图16看上去好象没有什么规律可循．
例5：现累积有100次刀具故障记录，当故障出现时该批刀具完成的零件数如下：(data5.txt)
459 362 624 542 509 584 433 748 815 505
612 452 434 982 640 742 565 706 593 680
926 653 164 487 734 608 428 1153 593 844
527 552 513 781 474 388 824 538 862 659
775 859 755 49 697 515 628 954 771 609
402 960 885 610 292 837 473 677 358 638
699 634 555 570 84 416 606 1062 484 120
447 654 564 339 280 246 687 539 790 581
621 724 531 512 577 496 468 799 544 645
764 558 378 765 666 763 217 715 310 851
x=load('data5.txt');
x=x(:);
hist(x) %%结果参见图17，很象图2所示的正态分布
figure
histfit(x) %%结果参见图18，加入了较接近的正态分布的密度曲线
图17例5的hist(x)图18例5的histfit(x)
3.参数估计
当我们可以基本确定数据集符合某种分布时，下一步我们就该确定这个分布的参数了．由于正态分布情况发生的比较多，故一般我们首先考虑的分布将是正态分布．考虑最多的也是正态分布情况．
对于未知参数的估计，可分两种情况：点估计与区间估计．
(1)点估计：构造样本与某个统计量有关的一个函数，作为该统计量的一个估计，称为点估计．Matlab统计工具箱中，一般采用最大似然估计法给出参数的点估计．可以证明：
①正态分布
中，
最大似然估计是
，
的最大似然估计是
；
②泊松分布
的最大似然估计是
；
③指数分布
的最大似然估计是
，等等．
例6：已知上述例1的数据服从正态分布，试求出和的值．
解：x=load('data1.txt');
x=x(:);
[mu, sigma] = normfit(x)
mu =
75.3417
sigma =
8.8768
因此，=mu=75.3412，=sigma2=8.87682=78.7982．
(2)区间估计：构造样本与某个统计量有关的两个函数，作为该统计量的下限估计与上限估计，下限与上限一般能够构成一个区间．这个区间作为该统计量的估计，称为区间估计．Matlab统计工具箱中，一般也采用最大似然估计法给出参数的区间估计．
例7：已知上述例1的数据集服从正态分布，试求出和
的置信度为95％的区间估计．
解：x=load('data1.txt');
x=x(:);
[mu, sigma muci, sigmaci] = normfit(x)
mu =
75.3417
sigma =
8.8768
muci =
73.7371
76.9462
sigmaci =
7.8781
10.1678
因此，73.7371
76.9462，7.8781
10.1678．
例8：从自动机床加工的同类零件中抽取16件，测得长度值为(data6.txt)：
12.15 12.12 12.01 12.08 12.09 12.16 12.06 12.13
12.07 12.11 12.08 12.01 12.03 12.01 12.03 12.06
已知零件长度服从正态分布，求零件长度的均值和标准差的置信度为99%的置信区间．
解：x=load('data6.txt');
x=x(:);
[mu, sigma, muci, sigmaci] = normfit(x, 0.01)
mu =
12.0750
sigma =
0.0494
muci =
12.0386
12.1114
sigmaci =
0.0334
0.0892
其中muci(1)、muci(2)分别是平均值在99％置信度下的上下限；而sigmaci(1)、sigmaci(2)分别是标准差在99％置信度下的上下限．
4．正态假设检验
对总体的分布律或分布参数作某种假设，根据抽取的样本观察值，运用数理统计的分析方法，检验这种假设是否正确，从而决定接受假设或拒绝假设，这就是假设检验问题．这里仅以正态假设检验为例，来说明假设检验的基本过程．
正态假设检验的一般过程是：
(1)对比正态分布的概率密度函数图，判断某统计量的分布可能服从正态分布；
(2)利用统计绘图函数normplot( )或weibplot( )进行正态分布检验．
(3)假设检验：利用Matlab统计工具箱给出的常用的假设检验方法的函数ttest(x,m,alpha)，进行显著性水平为alpha的t假设检验，以检验正态分布样本x(标准差未知)的均值是否为m．运行结果中，当h=1时，表示拒绝零假设；当h=0时，表示不能拒绝零假设．
例9：试说明例5所示的刀具的使用寿命服从正态分布，并且说明在方差未知的情况下其均值m取为597是否合理？
解：(1)对比正态分布的概率密度函数图(图17、图18)以及对正态分布的描述(一个变量如果是由大量微小的、独立的随机因素的叠加效果，那么这个变量一定是正态变量．比如测量误差、产品质量等都可用正态分布描述)，可得初步结论：该批刀具的使用寿命可能服从正态分布．
(2)利用统计绘图函数normplot(x)进行分布的正态性检验．由于：
x=load('data5.txt');
x=x(:);
normplot(x)
图19刀具寿命分布正态性检验
结果如图19所示，经观察这100个离散点非常靠近倾斜直线段，图形为线性的，因此可得出结论：该批刀具的使用寿命近似服从正态分布．
(3)利用函数ttest(x,m,alpha)进行显著性水平为alpha的t假设检验．由于：
x=load('data5.txt');
x=x(:);
h=ttest(x,597,0.05)
得：h = 0
检验结果：h=0，表示不拒绝零假设，说明所提出的假设“寿命均值为597”是合理的．
读者可以验证：当执行h=ttest(x,555,0.05)，将得到h = 1，表示拒绝零假设．请读者自行解释此结果的含义．
1．了解本实验中虽已提及但没有详细介绍的其余10种概率分布的密度函数，如Beta分布、Gamma分布、Weibull分布等，写出它们的概率分布的密度函数表达式(本实验的附录中已经列出一部分)，并画出相应的图形．
2．写出本实验所列出的10种概率累积分布函数表达式，并画出相应的概率累积分布函数图形．
3．用tabulate( )函数将例1、例2的分数数据按频数表的方式进行统计，每5分为一个分数段(可参见例4)，观察数据分布有什么规律．
4．用weibplot(x)函数进行例9的正态分布检验，比较与例9的差别．
5．例3给出的上海1998年来的月降雨量的数据(data3.txt)看上去很接近图10所示的
分布，但分布好象没有直接进行参数估计的函数，试寻求对此数据进行参数估计的可能方法．
6．向例3给出的上海1998年来的月降雨量的数据(data3.txt)中“补充”一些数据，使其看上去很接近正态分布，并求此时的均值和标准差的点估计与置信度为97%的区间估计．
7．在第6题基础上，说明在方差未知的情况下，其均值取为1150是否合理？
8．ttest( )函数的完整用法是：[h,sig,ci] = ttest(x,m,alpha,tail)
其中sig为观察值的概率，当sig为小概率时则对零假设提出质疑(这里的零假设为：
．也可以是其它形式，例如：
、
等)；ci为真正均值μ的1-alpha置信区间；不写tail，表示其取值为0．
说明：若h=0，表示在显著性水平alpha下，不能拒绝零假设；若h=1，表示在显著性水平alpha下，可以拒绝零假设．
若tail=0，表示备择(对立)假设为：
(默认，双边检验)；若tail=1，表示备择(对立)假设为：
(单边检验)；若tail=-1，表示备择(对立)假设为：
(单边检验)．
试用该函数求解如下问题：某种电子元件的寿命X(以小时计)服从正态分布，、均未知．现测得16只元件的寿命如下：
159 280 101 212 224 379 179 264 222 362 168 250
149 260 485 170
问当取alpha=0.05时：(1)是否有理由认为元件的平均寿命不大于225(小时)？(2)是否有理由认为元件的平均寿命不大于295(小时)？
9．查看函数ttest2( )的用法，并用于处理Matlab统计工具中的数据文件gas.mat．回答问题：一月份油价price1与二月份油价price2的均值是否相同？
附录：Matlab中的其它部分概率分布函数名及其数学意义列表：
函数名
对应分布
数学意义
batapdf
Beta分布
chi2pdf
卡方分布
fpdf
F分布
gampdf
Gamma分布
，
raylpdf
瑞利分布
，
tpdf
t分布
weibpdf
Weibull分布
