附录1 矩阵基础知识
矩阵的概念:矩阵就是矩形的数表。例如:
Apq
代表由p?q个数字排成的数表,我们称它为p行q列矩阵。矩阵用大写黑体字母表示。其下标表示它所包含的行列数,也可省略不写。用小写字母表示矩阵中的各个数字,如aij表示A矩阵中第i行第j列的那一个数字,称为矩阵的元素。有时也可用(aij)表示矩阵A。
向量是只有一行或一列的矩阵。当p = 1时,矩阵只有一行,称为行向量;当q = 1时,矩阵只有一列,称为列向量。
2. 矩阵的基本运算
(1)相等:两个矩阵A,B,若它们有所有元素对应相等,即对任意i,j,均有aij = bij,则称A与B相等,记为A = B。显然A与B相等的前提条件是它们有相同的行数和列数。
(2)加法:两个矩阵A,B,则A + B = C为一个新的矩阵,其元素为A和B的对应元素相加的和。即:若A = (aij), B = (bij), 则C = (cij) = (aij+bij)。显然加法也要求A,B矩阵有相同的行列数。
(3)乘法:两个矩阵Apq和Bqr, 则A·B = Cpr为一个新矩阵,其第i行第j列的元素cij为A的第i行元素与B的第列元素的乘积和,即:。显然矩阵乘法要求第一个矩阵的列数等于第二个矩阵的行数。
例1
如上面例题中结果的第一行第一列元素–33 = 3×1 + 2×2 + (–8)×5,第二行第一列元素13 = (-4)×1 + 6×2 + 1×5,等等。
注意:一般来说,矩阵乘法不满足交换律即AB≠BA。象上面的例子,BA根本就不能相乘,因为B有三列,而且A只有两行,不满足矩阵乘法的条件。再例如A1n为n阶行向量,Bn1为n阶列向量,则AB为一个数字,而BA为一个n×n阶的矩阵。
(4)转置:把矩阵A以它的主对角线(从左上到右下)为轴旋转180°,它的行变成列,列变成行,称为转置。记为A(。即:
若A = A(,则称A为对称矩阵。
(5)矩阵的行列式:若矩阵A为方阵,则我们可按某种规则从矩阵A计算出一个数作为它的值,这个值称为矩阵的行列式,记为。对于二阶方阵,它的行列式定义为主对角线乘积减去副对角线乘积。主对角线是指从左上到右下的对角线,而副对角线则是指从左下到右上的对角线。
例2 , ,
则 = 5(2–7(3 = –11 则 = 21(9 – (–3)(7 = 210
要计算高阶方阵的行列式,则需引入代数余子式的概念.通过它可把方阵的阶数逐次降低,直到只剩二阶行列式,从而可用上述方法求出最终结果。
子式:对于任意n阶行列式,删除任一元素aij所在的i行j列后所得n - 1阶行列式称为aij的子式。
代数余子式:子式乘以(-1)i+j,称为aij的代数余子式,记为A(ij)。
定理:行列式的值等于它任意一行或任意一列的所有元素与其代数子式的乘积之和。即:
,
称为按i行展开;或
,
称为按j列展开。
反复使用上述公式,直到各子式均变为2阶,然后可用前述方法求出其值。
若= 0,则称A为退化的方阵。
(6)单位阵。它是一个方阵,主对角线(从左上到右下的对角线)上元素均为1,其它元素均为0。记为Inn。它在矩阵乘法中起着类似数字1在数字乘法中的作用,所以称为单位阵。即:设A,I均为n·n方阵,则有AI = IA = A。换句话说,任何矩阵与单位阵(当然阶数必须适当)相乘,均不改变其数值。
(7)逆矩阵。若A为非退化方阵,即,则有与A同阶的方阵A–1存在,使
AA–1 = A–1A = I
A–1称为A的逆矩阵。它的求法为:设A = (aij),则:
其中A(ij)为aij的代数余子式。注意A–1中代数余子式的下标是经过转置的,即第i行第j列位置上是A的第j行i列元素aji的代数余子式。
附录2. 采用微软公司的Excel软件进行常见的统计计算。
Excel是一个功能十分强大的电子表格软件,它是微软公司办公软件Office中的一部分。利用它可以方便地进行许多计算工作,画图工作等,也包括常用的一些统计计算。使用这种通用办公软件的最大优点是普及率高,容易得到;其次是使用简单,不用记许多特殊指令;同时它也能复盖常用的统计方法,可满足一般工作时需要。另一方面,与许多著名的统计软件如SAS等相比,它也有一些明显的缺点,例如自动化程度不高,需要掌握一些基本统计公式;功能也不够强,有些统计计算不能做等。
在本附录中,我们假设读者已对Excel有一定了解,因此不再介绍Excel的基本用法。主要介绍以下几种统计计算:
1. 假设检验。包括正态总体的假设检验,离散分布的假设检验,以及用皮尔逊统计量进行非参数检验。
2. 方差分析。
3. 回归分析,包括简单作图。
§1 假设检验
正态总体单样本假设检验:
1. 统计知识复习:
若要检验方差,则统计假设为:
H0:( = ( 0 HA:( ≠ (0 (双边检验)
或: HA:( > (0 或 ( < (0 (单边检验)
统计量为:
若要检验均值,则统计假设为:
H0:( = ( 0 HA:( ≠ (0 (双边检验)
或:HA:( > (0 或 ( < (0 (单边检验)
统计量的选取则要分为以下两种情况:
a) 总体方差(2已知:u检验
b) 总体方差(2未知:t检验
2. 方差检验的计算方法:
设H0:( = (0,且原始数据在A1:A20位置。
1( 在空单元格(设为B1)中输入公式:
“= Var(A1:A20)*19/(0 ( 2 ↙”
这一步是计算(2统计量,其中Var为Excel的内部函数,功能为求指定数据的方差。“↙”表示回车(Enter)键。
2( 在B2格中输入:
“= chidist (B1, 19) ↙”
这一步是计算统计量所对应的概率,相当于查表。注意函数chidist返回的是单尾概率,即P(X > B1),而不是分布函数,即P(X < B1)。
3( 将B2中数据与(比较来确定是否接受H0:
双边检验:若(/2 < B2 < 1 - (/2,则接受H0;否则接受HA;
单边检验:若HA为:( > (0: 当B2 > (时接受H0;
若HA为:( < (0:当B2 < 1 - ( 时接受H0。
也可把上述公式一次输入:
“= chidist (Var(A1:A20)*19/(0 ( 2, 19) ↙”
上述公式中的19也可换为:Count(A1:A20)-1。Count这一内部函数可自动计算A1至A20中数字的个数。
3. 均值检验方法:
仍采用前述原始数据,零假设为:H0:( = (0。
总体方差(0 已知的情况:
1( 在空单元格(设为C1)中输入:
“= Ztest (A1:A20, (0,(0,) ↙”
内部函数Ztest可以直接算出u统计量所对应的单尾概率值。注意它返回的也是单尾概率,不是分布函数。
2( 仍按前述比较B1与(的同样方法比较C1与(,并决定是否接受H0。
总体方差(0未知:应采用t检验。
1( 在空单元格D1~D20中均填充上(0。
2( 在空格E1中输入:
“= Ttest (A1:A20, D1:D20, tails, 1) ↙”
其中tails为一参数,当进行单尾检验时,把它换成1;进行双尾检验时,换成2。最后一个数字“1”也是一个参数,它的用法我们后面将要介绍,这里应取值1。
3( 把E1格中计算出来的值与(相比,E1 > (时,接受H0;E1 < (时,拒绝H0。
注意:Ttest函数不区分统计量是大于0还是小于0,也不管是上单尾检验还是下单尾检验。因此进行单尾检验时可能出现错误拒绝。如当进行上单尾检验,即HA为( > (0,而观测数据平均值却明显小于(0时;或进行下单尾检验,即HA为( < (0,而观测数据平均值却明显大于(0时;在这两种情况下都会出现错误拒绝现象。使用中务请注意先进行直观检验,不属于以上两种情况时再进行统计检验,以免发生错误。
例1. (即本书例3.2)已知某种玉米平均穗重μ0 = 300g,标准差σ0 = 9.5g,喷药后,随机抽取9个果穗,重量分别为(单位为g):308,305,311,298,315,300,321, 294,320。问这种药对果穗重量是否有影响?
解:如表1,把果穗重原始数据填入A4:A12单元。
检验方差是否变化:在B5单元里输入:
“= Var(A4:A12)*8/9.5 ( 2,8)”
回车后,显示数字0.414234 。由于这一数字在0.025和0.975之间,因此接受H0,认为方差没有变化。
检验均值是否变化:由于方差已知,可采用Z-test。在B8单元里输入:
“= ztest(A4:A12,300,9.5)”
回车后,显示数字0.005763 。由于这一数字小于0.025,大于0.005,因此拒绝H0,喷药前后果穗重差异显著,但未达到极显著。
也可当作方差未知,直接进行T检验:
在C4:C12单元格中,填充数字300。
在D5单元格中输入:
“= ttest(A4:A12,D4:D12,2,1)”
回车后,显示数字0.037208 。由于这一数字小于0.05,大于0.01,因此拒绝H0,喷药造成的差异仍为显著,但未达极显著水平。
两种方法差异的讨论见本书例3.2。
表1. 例1计算结果
例1
μ0
300
σ0
9.5
果穗重
308
Chi-test
300
T-test
305
0.414234
300
0.037208
311
300
298
Z-test
300
315
0.005763
300
300
300
321
300
294
300
320
300
二、正态总体双样本假设检验:
1. 统计知识复习:
若要检验方差,统计假设为:H0:(1 = (2;HA:(1 ( (2。一般均为双边检验。统计量为:
其中m和n分别为第一和第二样本的样本容量。
若要检验均值,零假设为:H0:(1 = (2;
备择假设为:HA:(1 ( (2 (双边检验)
或:HA:(1 > (2 或 (1 < (2 (单边检验)
同时,还可能出现以下几种情况:
总体方差已知:u检验
(2)总体方差未知,但相等(即通过了F检验):t检验。
(3)总体方差未知,且不等(即未通过F检验):近似t检验。
近似服从t(df)
其中df = ,
(4)配对检验:用于两总体间明显正相关时。
令di = X1i – X2i, 对H0:(d = 0作单样本检验。
2. 方差检验方法:F检验,H0:(1 = (2;HA:(1 ( (2
假设两组数据分别位于A1:A10, B1:B10。
1( 在空格C1中输入:
“= Ftest (A1:A10, B1:B10) ↙”
注意Ftest返回的是F统计量的双尾概率,因此下一步可直接与α比较。
2°比较:C1<α则拒绝H0;C1>α则接受H0。
3. 均值检验方法:需区分几种情况:
两总体方差,已知:U检验。
1°在空格D1中输入:
“=(Average(A1:A10)-Average(B1:B10))/
sqrt(/count(A1:A10)+/count(B1:B10)) ↙”
这一步计算统计量的值,用了以下几个函数:Average:计算平均数;sqrt:计算平方根;count:计算指定区域中数字的个数。,应直接输入数值,或存贮该数值的位置。
2°在D2中输入:
“= Normsdist(D1) ↙”
这一步计算统计量对应的分布函数概率值。它返回的是分布函数取值(即P(X < x)),而不是尾区概率(一般为P(X > x))。注意在Excel中函数Normsdist是计算标准正态分布的取值,而Normdist是计算一般正态分布的取值。这里由于D1计算过程中已进行了标准化,因此应使用Normsdist。
3°将D2的数值与α比较:
双边检验:α/2 < D2 < 1-α/2时接受H0,否则拒绝H0。
单边检验:上单尾:HA:(1 > (2:当D2 < 1 - ( 时接受H0;
下单尾:HA: (1 < (2:当D2 > (时接受H0。
注意:由于Normsdist函数返回的是分布函数,而不是尾区概率,因此这里单边检验的接受域与使用Chidist和Ztest函数时正好相反。使用时请特别注意所用函数返回的倒底是分布函数还是尾区概率,否则单边检验时很容易出错误。
(2)两总体方差未知:由于Ttest函数中已考虑了方差未知时的各种可能,因此使用中很方便,只需改变一个参数的取值就可以了。
1( 在空格E1中输入:
“= Ttest (A1:A10, B1:B10, tails, type) ↙”
这一函数中后两个参数的取值与意义为:
tails = 1:单尾检验;tails = 2:双尾检验。
type = 1:配对检验;type = 2:方差相等;type = 3:方差不等。
使用时直接把参数换为相应的数值即可。由于函数返回的数值为尾区概率,因此可直接与(相比。
2( 把E1的数值与(比较,E1 > ( 时,接受H0,否则拒绝H0。
注意:单尾检验中不管两个均值谁大Ttest给出的概率都是相同的。因此在上单尾检验(HA:(1 > (2)中第一个样本均值偏小或下单尾检验(HA:(1 < (2)中第一个样本均值偏大都有错误拒绝H0的可能,使用时需要特别注意。
例2. (即本书例3.3) 两发酵法生产青毒素的工厂,其产品收率的方差分别为,现甲工厂测得25个数据,,乙工厂测得30个数据,,问它们的收率是否相同?
解:由于两总体方差已知,可采用正态分布进行检验。在空格E3中输入:
“= normsdist((3.71 – 3.46) / sqrt(0.46 / 25 + 0.37 / 30))”
回车后,显示数字0.923073。由于这一数字在0.025和0.975之间,因此接受H0,认为这两个工厂的收率相同。
例3. 新旧两个小麦品系进行对比试验,旧品系共收获25个小区,新品系收获20个小区,产量(公斤)如下表。问新品系是否值得推广?
旧品系
34.6
38.1
40.5
36.2
39.5
34.1
39.5
38.0
37.9
38.4
39.5
32.9
37.2
新品系
37.1
38.9
39.1
36.2
39.8
40.8
41.2
38.7
40.3
41.5
40.3
37.7
40.9
旧品系
30.8
38.1
38.3
39.3
34.9
31.8
34.5
35.9
38.2
39.7
33.9
36.0
新品系
38.7
37.2
41.9
38.6
39.2
38.2
40.6
解:首先检验方差是否相等:在空格中输入:
“= ftest(E3:E27,F3:F22)”
回车后,显示数字0.024704。由于这一数字小于0.05,因此拒绝H0,认为方差不相等。应采用近似检验。
检验均值是否相等:根据题意,应为单侧检验。在另一空格输入:
“= ttest(E3:E27,F3:F22,1,3)”
回车后,显示数字0.000095。由于这一数字小于0.01,因此拒绝H0,认为新品系极显著地优于旧品系,值得推广。
例4.(即本书例3.6) 10名病人服药前后血红蛋白含量如下:
病人号
1
2
3
4
5
6
7
8
9
10
服药前(x)
11.3
15.0
15.0
13.5
12.8
10.0
11.0
12.0
13.0
12.3
服药后(y)
14.0
13.8
14.0
13.5
13.5
12.0
14.7
11.4
13.8
12.0
问该药是否引起血红蛋白含量变化?
解:根据题意,应采用配对检验。在空格输入:
“= ttest(I3:I12,J3:J12,2,1)”
回车后,显示数字0.223742。由于这一数字大于0.05,因此接受H0,认为服药前后血红蛋白含量没有显著变化。
三、非参数检验:皮尔逊(Pearson)统计量。
统计知识复习:
皮尔逊定理:当P1,P2,…… Pr为总体的真实概率分布时,统计量
随n增加而渐近于自由度为r – 1的(2分布。
若令Oi = ni, Ti = npi,则上式变为:
用途:吻合度检验,列联表独立性检验。
限制条件:各Ti ( 5。
2. 列联表独立性检验:
对列联表进行独立性检验首先应计算理论值。对列联表独立性检验来说,理论值计算公式为:
下面结合例题,介绍计算过程。
例5.(即本书例3. 22) 下表是对某种药的试验结果:
给药方式
有效
无效
口服
58
40
注射
64
31
问给药方式对药效果是否有影响?
解: 表2. 例5的计算结果
有效
无效
口服
58
40
Chi-test
注射
64
31
0.238468
理论值
61.94819
36.05181
60.05181
34.94819
如上表,原始数据在区域M3:N4。计算步骤为:
1( 首先计算理论值:在空格M6输入:
“=SUM($M3:$N3)*SUM(M$3:M$4)/SUM($M$3:$N$4)”
回车后,显示数字61.94819。把M6复制到M7和N6、N7,得到各理论值。请注意上式中美元符号的位置,只有位置正确才能保证复制结果正确。
2( 进行统计检验:在P4单元格输入:
“= CHITEST(M3:N4,M6:N7)”
回车后,显示数字0.238468。把P4的值与α相比:当P4 >α时接受H0,即列联表的行与列相互独立;否则拒绝H0,即行与列不独立。由于这一数字大于0.05,因此接受H0,认为给药方式与药效无关。
此函数的第一个参数为观测值所在区域,第二个参数为理论值所在区域。这两个矩形区域行列数必须相同。返回值为皮尔逊统计量对应的(2分布的尾区概率,其自由度为(r–1)(c–1),其中r, c分别为数据区的行数和列数。如果数据区只有一行或一列,则自由度为数据个数减1。这正是列联表独立性检验所需的自由度。
例6.(即本书例3. 24) 为检测不同灌溉方式对水稻叶片衰老的影响,收集如下资料:
灌溉方式
绿叶数
黄叶数
枯叶数
深水
浅水
湿润
146
183
152
7
9
14
7
13
16
问叶片衰老是否与灌溉方式有关?
解: 表3. 例6的计算结果
灌溉方式
绿叶数
黄叶数
枯叶数
深水
146
7
7
Chi-test
浅水
183
9
13
0.229248
湿润
152
14
16
理论值
140.6947
8.775137
10.53016
180.2651
11.24314
13.49177
160.0402
9.981718
11.97806
如表3,原始数据在区域Q3:S5。首先计算理论值:在空格Q7输入:
“=SUM($Q3:$S3)*SUM(Q$3:Q$5)/SUM($Q$3:$S$5)”
回车后,显示数字140.6947。把Q7复制到区域Q7:S9,得到各理论值。请注意上式中美元符号的位置,只有位置正确才能保证复制结果正确。
在U4单元格输入:
“=CHITEST(Q3:S5,Q7:S9)”
回车后,显示数字0.229248。由于这一数字大于0.05,因此接受H0,认为叶片衰老与灌溉方式无关。
3. 吻合度检验:
对吻合度检验来说,理论值的计算显然与理论分布的类型有关,(2检验的自由度也可能发生变化。例如对正态分布的吻合度检验,如果总体参数(,σ2已知,则统计量自由度为数据个数减1;但若总体参数未知,用样本均值与方差S2代替,则统计量自由度也要再减2。此时直接用Chitest得到的尾区概率就不对了,需要再作一下变换(见例7)。现以正态分布为例介绍吻合度检验计算步骤。
例7.(即本书例3.19) 调查了某地200名男孩身高,得,分组数据见下表。男孩身高是否符合正态分布?
组号
区间
Oi
1
(-∞, 126)
8
2
[126, 130)
13
3
[130, 134)
17
4
[134, 138)
37
5
[138, 142)
55
6
[142, 146)
33
7
[146, 150)
18
8
[150, 154)
10
9
[154, +∞)
9
解:计算结果如表4。计算过程为:
1°在C3至C11中填入身高区间的上界。最后一个应为无穷大,填入足够大的数即可。
2°在D3格中输入:
“=NORMDIST(C3,139.5,7.42,1)”
这一步是计算正态分布值。第一个参数为区间上限;第二个参数为均值;第三个参数为标准差;第四个参数为0时计算密度函数,为1时计算分布函数。
把D3复制到D4:D11。
3°计算各区间的概率。在E3中输入“=D3”,在E4中输入“=D4-D3”,并复制E4到E5:E11。
4°计算理论值:在G3输入
“=E3*200”
并复制G3到G4至G11。
5°计算统计量:在H3输入:
“=(F3-G3)∧2/G3”
把H3复制到H4至H11,并在H12输入:
“=Sum(H3:H11)”
另一种计算统计量的方法为:在I3输入:
“=Chitest(F3:F11, H3:H11)”
在I6输入:
“=Chiinv(I3, 8)”
可见I6的数值与H12是相同的。
6°计算统计量对应的尾区概率:在I9输入:
“=chidist(I6,6) ↙”
7°将I9与α相比,当I9 > α时,接受H0,所观察数据符合正态分布;当I9 ≤ α时,拒绝H0,数据不符合正态分布。在本题中,I9的数值为0.085446 > α,因此应接受H0,可认为男孩身高符合正态分布。计算结果如下表。
表4. 例7的计算结果
组号
区间
边界
正态分布
概率
观察值
理论值
(Oi-Ti)2/Ti
Chi-test
1
<126
126
0.034425
0.034425
8
6.884924
0.180597
0.196303
2
[126,130)
130
0.100216
0.065791
13
13.15823
0.001903
3
[130,134)
134
0.229274
0.129058
17
25.81163
3.008134
统计量
4
[134,138)
138
0.419897
0.190623
37
38.12467
0.033178
11.09629
5
[138,142)
142
0.631914
0.212017
55
42.40336
3.742049
6
[142,146)
146
0.809488
0.177574
33
35.51478
0.17807
P
7
[146,150)
150
0.92148
0.111992
18
22.39832
0.863689
0.085446
8
[150,154)
154
0.97466
0.05318
10
10.63609
0.038041
9
>154
100000
1
0.02534
9
5.068004
3.050627
和
11.09629
分位数
12.59158
本来Chitest函数返回的就是尾区概率,但它使用的自由度为数据个数减1,而现在应使用数据个数减3为自由度,因此要使用函数Chiinv先把尾区概率变回统计量的值,然后再用Chidist求出正确自由度下的尾区概率。
注意使用不同概率模型时,自由度的变化是不同的。一般来说,模型中使用几个统计量代替未知参数,自由度就要在原来的基础上再减少几个。例如上面的例题用了样本期望和方差代替未知参数,因此自由度比正常的Pearson统计量少2;本书中例3.20,统计模型中没有未知参数,因此自由度没有变化;例3.21有一个参数需用统计量代替,因此自由度需再减一。
四、常用离散分布的统计计算:
离散分布统计计算中关键一点是正确建立尾区。尾区是从观察值开始,向对H0成立不利的方向求和。例如水质检验要求大肠杆菌不大于2个/毫升,取2毫升检验,发现5个细菌,问是否判断超标。此时H0为:μ≤4,对H0成立不利的方向应是细菌数增加,因此尾区概率应为:。其中pi为2毫升水样中出现i个细菌的概率。
尾区建立以后用Excel提供的函数求概率是很容易的。然后根据是单尾或双尾检验与α或α/2比较,若尾区概率大于α或α/2,则接受H0;否则拒绝。我们先介绍一下有关函数所需参数的意义,然后结合例题说明使用方法。
1. 二项分布有关函数:
Binomdist (n, N, p, C)
用于计算二项分布的概率或累积概率。其中n:成功次数;N:总实验次数;p:成功概率;C:参数,取值为1时计算从0到n的累积概率,取值为0时计算成功n次的概率。
Critbinom (N, p, α)
用于求二项分布累积概率大于指定临界概率时的最小成功次数。其中参数意义为:N:总实验次数;p:成功概率,α:临界概率。
2. 超几何分布有关函数:
Hypgeomdist (k, n, M, N)
用于计算超几何分布概率。其参数意义为:k:样本中的成功数;n:样本数;M:总体中的成功数;N:总体中个体数。
3. 负二项分布有关函数:
Negbinomdist (x, r, p)
用于计算负二项分布概率。其参数意义为:x:失败次数;r:成功次数;p:成功概率。其最后一次实验必定是成功的。
4. 泊松分布函数:
Poisson (x, λ,c)
用于计算泊松分布概率或累积概率。参数意义为:x:成功次数;λ:平均数;c:参数,取值为1时计算成功次数小于等于x的累积概率;取值为0时计算成功x次的概率。
例8. 产品废品率小于等于0.03为合格。抽检20个样品发现2个废品,该批产品是否合格?若发现3个废品呢?
解:1°在空格B5格中输入:
“=1-Binomdist (1, 20, 0.03, 1)”
回车后,显示数字0.119838。由于尾区是从2累加到20,而Binomdist函数是从0累加到指定值,因此这里应指定第一个参数为1。
2°将B5与α相比:由于B5>α=0.05,故接受H0,发现2个废品可认为合格。
3°在空格B6格中输入:
“=1-Binomdist(2, 20, 0.03, 1)”
回车后,显示数字0.021008。
4°将B6与α相比:由于B6<α=0.05,故拒绝H0,发现3个废品应认为不合格。
例9. 水质检验要求每毫升水中大肠杆菌不得超过3个。现取1毫升检验,发现6个细菌,水质是否合格?若2毫升发现12个细菌呢?
解:1°在空格B12中输入:
“=1-Poisson(5, 3, 1)”
回车后,显示数字0.083918。与前一题类似,H0:λ≤3;故尾区应向多的方向累加。对1毫升发现6个细菌,尾区为:。因此第一个参数应取为5。
2°将B12与α=0.05相比,由于B12>α,故接受H0,1毫升发现6个细菌应认为合格。
3°在空格B13中输入:
“=1-Poisson(11,6,1)”
回车后,显示数字0.020092。由于现改为检测2毫升,故λ应取为6;尾区为:,因此第一个参数应取为11。
4°将B13与α相比,由于B13<α,故拒绝H0,2毫升发现11个细菌应认为不合格。
§2 方差分析
方差分析是重要的统计方法之一,它主要用于比较多组数据的平均数是否相同。Excel有一个用于进行方差分析的宏,但必须进行安装才能使用,同时也不太完善,例如不能区分因素类型等。因此本节中既介绍利用Excel的统计函数,手工进行方差分析的方法,也介绍利用宏自动计算,然后根据需要再对结果加以调整的方法。
一、统计知识复习:
方差分析中的因素可分为固定因素和随机因素,不同因素类型对方差分析的影响主要表现在应选用不同统计量及对结果解释不同。因此进行方差分析应注意区分因素类型。
1. 单因素方差分析:总平方和及自由度可作如下分解:
总平方和:SST = SSA + SSe
自由度: an–1 = (a–1) + a(n–1)
统计量: F = MSA / MSe ~ F(a–1, a(n–1))
当H0不成立,即各水平的平均数有差异时,F统计量有偏大的趋势,因此可进行上单尾检验。若因素为固定因素,结论只适用于参加检验的几个水平;若为随机因素,则可推广到一切水平。
2. 双因素交叉分组方差分析:平方和及自由度分解为:
平方和:SST = SSA + SSB + SSAB + SSe
自由度:abn–1 = (a–1) + (b–1) + (a–1)(b–1) + ab(n–1)
统计量的选择依赖于因素类型:
a) 固定效应模型:
FA = MSA / MSe ~ F(a–1, ab(n–1))
FB = MSB / MSe ~ F(b–1, ab(n–1))
FAB = MSAB / MSe ~ F((a–1)(b–1), ab(n–1))
b) 随机效应模型:
FA = MSA / MSAB ~ F(a–1, (a–1)(b–1))
FB = MSB / MSAB ~ F(b–1, (a–1)(b–1))
FAB = MSAB / MSe ~ F((a–1)(b–1), ab(n–1))
c) 混合模型:(A固定,B随机)
FA = MSA / MSAB ~ F(a–1, (a–1)(b–1))
FB = MSB / MSe ~ F(b–1, ab(n–1))
FAB = MSAB / MSe ~ F((a–1)(b–1), ab(n–1))
均为上单尾检验。固定因素的结果不能推广,随机因素则可推广到一切水平。
3. 双因素系统分组方差分析:
系统分组与交叉分组的不同点在于对应于一级因素的不同水平,系统分组的二级因素各水平可取不同值。此时SSB与SSAB无法分离。其平方和与自由度的分解为:
平方和: SST = SSA + SSB + SSe
自由度:abn–1 = (a–1) + a(b–1) + ab(n–1)
统计量:FB = MSB / MSe ~ F(a(b–1), ab(n–1))
检验因素A的统计量则取决于因素B的类型:
B固定:FA = MSA / MSe ~ F(a–1, ab(n–1))
B随机:FA = MSA / MSB ~ F(a–1, a(b–1))
结果解释仍为固定因素不可推广,随机因素可推广。
4. 多重比较:由于Excel中没有Duncan法,Newman-Q法等所需要的系数表,因此无法使用这些多重比较方法。这里我们只介绍可用的最小显著差数法:
统计量:
一般为双尾检验。其中,为两个处理的平均数,n为重复数。
二、方差分析的手工计算方法:
例10.(即本书例4.1) 用4种不同的配合饲料饲养30日龄小鸡,10天后计算平均日增重,得以下数据:4种饲料效果是否相同?
饲料
日增重值Xij
1
55 49 62 45 51
2
61 58 52 68 70
3
71 65 56 73 59
4
85 90 76 78 69
解:把数据输入Excel,如表5。
表5. 例10的计算结果
饲料1
饲料2
饲料3
饲料4
重复1
55
61
71
85
重复2
49
58
65
90
重复3
62
52
56
76
重复4
45
68
73
78
重复5
51
70
59
69
平均
52.4
61.8
64.8
79.6
方差
41.8
54.2
54.2
66.3
S2
146.02895
SST
2774.55
MSA
636.18333
F
11.753965
S2xi.
127.23667
SSA
1908.55
MSe
54.125
P
0.0002556
SSe
866
计算步骤为:
1°计算各饲料日增重平均值:在B8中输入:
“=Average (B3:B7)”
回车后,显示数字52.4。把B8复制到C8:E8,得到各平均值。
2°计算SST:在B11中输入:
“=Var(B3:E7) ”
回车后,显示数字146.02895。这是全部原始数据的样本方差。在D11中输入:
“=19*B11”
回车后,显示数字2774.55。这就是总平方和SST。公式中19 = an - 1,在本题中,a = 4,n = 5。
3°计算SSA:在B12中输入:
“=Var(B8:E8)”
回车后,显示数字127.23667。这是各平均值的样本方差。在D12中输入:
“=15*B12”
回车后,显示数字1908.55。这就是平方和SSA。公式中15 = n(a-1)。
4°计算SSe: 在D13中输入:
“=D11-D12 ”
回车后,显示数字866。这就是平方和SSe。
5°计算MSA,MSe:在F11中输入:
“=D12/3”
回车后,显示数字636.18333。这就是MSA,其中3=a-1;在F12中输入:
“=D13/16”
回车后,显示数字54.125。这就是MSe,其中16=a(n-1)。
6°计算统计量及其对应概率:在H11中输入:
“=F11/F12”
回车后,显示数字11.753965。这就是F统计量。在H12中输入:
“=Fdist(H11, 3, 6)”
回车后,显示数字0.0002556。这就是F统计量对应的概率值。其中3为统计量分子自由度a-1, 16为分母自由度a(n-1)。由于H12<α=0.01,应拒绝H0,各饲料有极显著差异。
本题属于固定模型,因此可进一步进行多重比较。结果见表6。具体步骤为:
表6. 例10的多重比较
饲料4
饲料3
饲料2
饲料1
79.6
64.8
61.8
52.4
饲料1
52.4
2.48314E-05
0.0169434
0.060429
饲料2
61.8
0.001490289
0.5282252
饲料3
64.8
0.005807619
7°复制平均数,并进行排序:
把B8:E8复制到K3:N3,用“选择性粘贴”,“数值”。
把B2:E2复制到K2:N2。
用鼠标选择K2:N3区域,然后对它进行排序:用“数据”菜单下的“排序”命令,点击“选项”,在出现的菜单中选择“按行排序”,点击“确定”关闭“选项”菜单;并指定关键字为“按行3”,“递减”。点击“确定”后,就完成了对平均数从大到小的排序。排序过程中,平均数和它对应的处理是连在一起排序的,这样有助于判断是那些处理之间有显著差异。
再把L2:N3复制到I4:J6,用“选择性粘贴”,“转置”命令。用鼠标选择I4:J6区域,再按列J递增排序。
8°计算各平均数间T统计量所对应的概率:在K4中输入:
“=Tdist ((K$3-$J4)/SQRT(2*$F$12/5),16,2)”
上式中有些行或列号前有“$”号,是为了在以后的复制中使相应的行号或列号不变化;SQRT为求平方根函数;$F$12为MSe的存贮地址;它前边的数字2为公式中的常数;后边的数字5为本题中的常数n;再后边的数字16是本题中MSe的自由度;最后的数字2为函数Tdist的参数,表示计算双尾概率。
把K4复制到K5,K6,L4,M4;再把K5复制到L5。
9°把上面计算出的各概率与0.05相比,小于0.05的为差异显著;再与0.01相比,小于0.01的为差异极显著。可用不同颜色分别表示。在本题中,饲料4与饲料1,2,3差异均达极显著;而饲料3与饲料1差异显著;其他差异不显著。从以上结果看,饲料4的增重最大,应是最好的。
例11.(即本书例4.3) 为选择最适发酵条件,用三种原料、三种温度进行了实验,得结果如下表。请进行统计分析。
原料
种类
(A)
温 度(B)
30℃
35℃
40℃
1
41
49
23
25
11
13
25
24
6
22
26
18
2
47
59
50
40
43
38
33
36
8
22
14
18
3
35
53
50
43
38
47
44
55
33
26
29
30
解:把数据输入Excel表,并排列如下:
表7. 例11方差分析结果
因素A
温度30
温度35
温度40
因素B
原料1
原料2
原料3
原料1
原料2
原料3
原料1
原料2
原料3
重复1
41
47
35
24
43
38
6
8
33
重复2
49
59
53
25
38
47
22
22
26
重复3
23
50
50
13
33
44
26
14
19
重复4
25
40
43
11
36
55
18
18
30
平均
34.5
49
45.25
18.25
37.5
46
18
15.5
27
温度平均
42.9167
33.9167
20.1667
原料平均
23.5833
34
39.4167
SST
7170
SSST
5513.5
SSA
3150.5
MSA
1575.25
FA
5.67E-07
SSB
1554.17
MSB
777.083
FB
0.000132
SSAB
808.833
MSAB
202.208
FAB
0.025322
SSE
1656.5
MSE
61.3519
按以下步骤进行计算:
1°计算各处理平均数:在B8输入:
“=Average(B4:B7)”
回车后,显示数字34.5。把B8复制到C8:J8,得到各处理的平均数。
2°计算因素A,即温度的各水平平均值:在B9输入:
“=Average (B8:D8)”
回车后,显示数字42.9167。再用鼠标标记B9:D9,点“跨列居中”键。再把B9复制到E9,H9。这样就得到了各温度的平均值。
3°计算因素B,即原料各水平的均值:在B10输入:
“=Average (B8,E8,H8)”
回车后,显示数字23.5833。把B10复制到C10,D10。这样就得到了各原料的平均值。
4°计算总平方和SST:在B12输入:
“=35*Var(B3:J6)”
回车后,显示数字7170。其中35=abn-1,由于本题中a=b=3, n=4,故总自由度为35。
5°计算次总平方和SSST:在B13输入:
“=32*Var(B8:J8)”
回车后,显示数字5513.5。其中32=n(ab-1)。
6°计算SSA:在B14输入:
“=24*Var(B9,E9,H9)”
回车后,显示数字3150.5。其中24=bn(a-1)。
7°计算SSB:在B15输入:
“=24*Var(B10:D10)”
回车后,显示数字1554.17。其中24=an(b-1)。
8°计算SSAB:在B16输入:
“=B13-B14-B15”
回车后,显示数字808.8333。
9°计算误差平方和SSe:在B17输入:
“=B12-B13”
回车后,显示数字1656.5。
10°计算各个均方:在D14输入:
“=B14/2”
回车后,显示数字1575.25,为MSA。其中2=a-1。在D15输入:
“=B15/2”
回车后,显示数字777.0833,为MSB。其中2=b-1。在D16输入:
“=B16/4”
回车后,显示数字202.2083,为MSAB。其中4=(a-1)(b-1)。在D17输入:
“=B17/27”
回车后,显示数字61.35185。为MSe,其中27=ab(n-1)。
11°计算各统计量对应的尾区概率:在F14输入:
“=Fdist(D14/D17,2,27)”
回车后,显示数字5.67×10-7,为统计量FA对应的概率值。其中D14/D17为FA统计量的值,2,27分别为其分子,分母自由度。在F15输入:
“=Fdist (D15/D17, 2, 27)”
回车后,显示数字0.000132,为统计量FB对应的概率值。其中D15/D17为FB统计量的值,2,27分别为其分子分母自由度。在F16输入:
“=Fdist(D16/D17, 4, 27)”
回车后,显示数字0.025322,为统计量FAB对应的概率值。其中D16/D17为FAB统计量的值,4,27分别为其分子分母自由度。
12°将F14,F15,F16中的数值分别与α比较,若大于α,则接受H0,认为该因素影响不显著;否则影响显著。对于本题来说,A,B两因素影响均达极显著水平,它们的交互作用也达到了显著水平。
以上是认为A,B均为固定因素的检验方法。若认为有一个或两个因素为随机因素,则应相应改变统计量及自由度:若认为两因素均为随机因素,则应在检验主效应时改用MSAB为分母,即将F14,F15中输入的公式分别改为:
“=Fdist(D14/D16, 2, 4) ”
“=Fdist(D15/D16, 2, 4) ”
其他不变,但现在结果可推广到A,B因素的一切水平。若只有一个因素为随机,设A固定,B随机,则F15公式同固定模型,F14同随机模型,即分别为:
“=Fdist(D14/D16, 2, 4) ”,
“=Fdist(D15/D17, 2, 27) ”
比较方法仍不变,但A因素结果不能推广,B因素则可以。多重比较在各处理的平均数之间进行,方法同单因素方差分析,本例题仅给出结果(见表8),不再重复计算步骤。
表8中第一列是处理条件,即具体温度和原料种类;第二列是该处理平均数;第三列是平均数排序的序号。从下表可知,X1~X4和X6~X9两组内各平均数之间差异除X1与X4及X6与X9之外均不显著;而这两组间差异大多达到显著或极显著。两组中的X1,X2,X3以及X7,X8,X9更是没有多少差异。因此在实践中可根据实际问题要求选平均数大的还是小的,从这两组中选取一组;再根据其他条件如成本,原料来源,操作方便等从中选取需要的处理。
表8. 例11多重比较结果
温度40
原料2
温度40
原料1
温度35
原料1
温度40
原料3
温度30
原料1
温度35
原料2
温度30
原料3
温度35
原料3
温度30
原料2
15.5
18
18.25
27
34.5
37.5
45.25
46
49
x9
x8
x7
x6
x5
x4
x3
x2
x1
x1
1.86E-06
6.15E-06
6.94E-06
4.76E-04
1.43E-02
4.75E-02
5.04E-01
5.92E-01
x2
7.83E-06
2.62E-05
2.96E-05
1.95E-03
4.75E-02
1.36E-01
8.93E-01
x3
1.12E-05
3.78E-05
4.26E-05
2.75E-03
6.28E-02
1.73E-01
x4
4.76E-04
1.55E-03
1.74E-03
6.87E-02
5.92E-01
x5
1.95E-03
6.05E-03
6.75E-03
1.87E-01
x6
4.75E-02
1.16E-01
1.26E-01
x7
6.24E-01
9.64E-01
x8
6.55E-01
以上是交叉分组方差分析的做法。系统分组方差分析与交叉分组的最大不同点是SSB与SSAB不可分离,因此计算变得较为简单。下面以例12说明具体计算步骤。
例12.(即本书例4.8) 为比较4种酶在不同温度下的催化效率,特设计如下实验:由于文献记载各酶最适温度分别为30℃,25℃,37℃,40℃,现设定温度水平如下:最适温-5℃,最适温,最适温+5℃。其他条件均保持一致。保温2小时后,测定底物消耗量(毫克)。全部实验重复三次,得结果如下:
温度
酶种类
A1
A2
A3
A4
偏低
适宜
偏高
14.4, 15.2, 13.5
15.9, 15.1, 14.4
13.8, 12.9, 14.6
13.5, 14.4, 15.2
15.1, 16.4, 15.8
15.7, 14.8, 16.0
14.5, 16.3, 15.4
16.4, 18.1, 16.7
15.8, 14.7, 14.1
11.2, 9.8, 10.5
12.5, 10.9, 11.6
10.3, 11.4, 9.9
请进行统计分析
解:把原始数据输入Excel如表9中A2:M6区域,然后计算如下:
表9. 例12计算结果
因素A
酶A1
酶A2
酶A3
酶A4
因素B
温度偏低
温度适宜
温度偏高
温度偏低
温度适宜
温度偏高
温度偏低
温度适宜
温度偏高
温度偏低
温度适宜
温度偏高
重复1
14.4
15.9
13.8
13.5
15.1
15.7
14.5
16.4
15.8
11.2
12.5
10.3
重复2
15.2
15.1
12.9
14.4
16.4
14.8
16.3
18.1
14.7
9.8
10.9
11.4
重复3
13.5
14.4
14.6
15.2
15.8
16
15.4
16.7
14.1
10.5
11.6
9.9
平均
14.3667
15.1333
13.7667
14.3667
15.7667
15.5
15.4
17.0667
14.8667
10.5
11.6667
10.5333
方差
0.72333
0.56333
0.72333
0.72333
0.42333
0.39
0.81
0.82333
0.74333
0.49
0.64333
0.60333
A平均
14.4222
15.2111
15.7778
10.9
S2
4.614921
SST
161.5222
Sxij
4.43037
SSST
146.2022
Sxi.
4.797119
SSA
129.5222
MSA
43.17407
PA
7.504E-12
SSB
16.68
MSB
2.085
PB
0.0116
SSe
15.32
MSe
0.638333
1°计算各处理平均数:在B8输入:
“=Average (B4:B6)”
回车后,显示数字14.3667。把B8复制到C8:M8,得到各处理平均数。
2°计算A因素,即各酶种的平均值:在B10输入:
“=Average (B8:D8)”
回车后,显示数字14.4222。然后标记B10:D10区域,点“跨列居中”键,再把B10复制到E10,H10,K10。
由于不同酶种所需温度不同,再求温度平均数已无意义。
3°计算总平方和SST:在B12输入:
“=Var(B4:M6)”
回车后,显示数字4.614921,为全部原始数据的方差。在B13输入:
“=(4*3*3-1)*B12”
回车后,显示数字161.5222,为总平方和SST。公式中4,3,3分别为本例题中a, b, n的取值,下同。
4°计算次总平方和SSST:在B14输入:
“=Var(B8:M8)”
回车后,显示数字4.43037,为各处理平均数的方差。在B15中输入:
“=3*(4*3-1)*B14”
回车后,显示数字146.2022,为次总平方和SSST。公式中第一个3为n,另外两数分别为a,b。
5°计算SSA:在B16输入:
“=Var(B9,E9,H9,K9)”
回车后,显示数字4.797119,为各酶种平均数的方差。在B17输入:
“=3*3*(4-1)*B16”
回车后,显示数字129.5222,为SSA。
6°计算SSB(这里实际相当交叉分组的SSB+SSAB):在B18输入:
“=B15-B17”
7°计算SSe: 在B19输入:
“=B13-B15”
8°计算各因素均方,在E17,E18,E19中分别输入:
“=B17/3”
“=B18/8”
“=B19/24”
显示数字分别为:43.17407,2.085,0.63833。公式中3,8,24分别为各平方和的自由度,其表达式分别为a-1, a(b-1), ab(n-1)。
9°计算统计量FA和FB所对应的尾区概率:在I17,I18输入:
“=Fdist (E17/E19, 3, 24)”
“=Fdist (E18/E19, 8, 24)”
回车后,显示数字分别为:7.504×10-12,0.0116。以上是B为固定因素时的计算公式。若B为随机因素,则I17中的公式应改为:
“=Fdist (E17/E18, 3, 8)”
回车后,显示数字为:0.00397。
10°将I17,I18中的数值与α相比,大于α时接受H0,该因素影响不显著;小于α时拒绝H0,该因素影响显著。在本题中,A因素即酶的种类影响极显著,B因素即温度(包括交互效应)影响显著,但未达极显著。
若需要也可对各处理平均数进行多重比较,方法与前相同,不再重复。
三、采用Excel中的宏进行方差分析。
采用宏进行方差分析的优点是计算都可自动完成,但它只能进行交叉分组固定模型的分析,如果是其他模型则可利用其中间结果再重新计算。
要利用这种方法,首先要加载宏:点击“工具”菜单下的“加载宏”命令,出现一对话框,在其中选取“分析工具库”和“分析工具库-VBA函数”,再点击“确定”钮。然后,在“工具”菜单下就会出现“数据分析”命令,点击后出现对话框,其中有方差分析,相关系数,协方差分析,指数平滑等多种分析工具可用。现在我们就介绍一下用它进行单因素和双因素方差分析的方法。
1. 单因素方差分析
在单因素方差分析中,因素类型对分析过程没有影响,因此不用重新计算。只需把数据输入Excel,就可利用宏进行计算。
例13.仍采用与例10相同的数据:用4种不同的配合饲料饲养30日龄小鸡,10天后计算平均日增重,得以下数据:4种饲料效果是否相同?
饲料
日增重值Xij
1
55 49 62 45 51
2
61 58 52 68 70
3
71 65 56 73 59
4
85 90 76 78 69
解:如下表,把原始数据输入:
表10. 例13方差分析结果
饲料1
饲料2
饲料3
饲料4
重复1
55
61
71
85
重复2
49
58
65
90
重复3
62
52
56
76
重复4
45
68
73
78
重复5
51
70
59
69
方差分析:单因素方差分析
SUMMARY
组
计数
求和
平均
方差
饲料1
5
262
52.4
41.8
饲料2
5
309
61.8
54.2
饲料3
5
324
64.8
54.2
饲料4
5
398
79.6
66.3
方差分析
差异源
SS
df
MS
F
P-value
F crit
组间
1908.55
3
636.1833
11.75396
0.000256
3.238867
组内
866
16
54.125
总计
2774.55
19
计算过程为:
1°如上表,把原始数据输入B2:E7区域。
2°点击“工具”,“数据分析”命令,在对话框中选“单因素方差分析”,然后点“确定”。
3°在出现的对话框中,指定输入区域为B2:E7,分组方式为“列”选定“标志位于第一行”,再指定输入区域为A9,然后点“确定”。
4°结果解释:结果包括两个表,第一个summary表中有各水平的重复数,总和,平均数,方差等数据;第二个方差分析表,给出组内和组间平方和,自由度,均方,F统计量,尾区概率,F分位数值等数据。尾区概率(表中标记为“P-Value”)小于α时拒绝H0,大于时接受H0。本题尾区概率为0.000256 < 0.01,应拒绝,差异极显著。这一结果与例10手工计算相同。
2. 双因素方差分析:
宏提供有重复和无重复双因素方差分析,但都是针对固定模型,交叉分组。若为其他模型则需利用中间结果重新计算。仍采用例11说明使用方法。
例11.(即本书例4.3) 为选择最适发酵条件,用三种原料、三种温度进行了实验,得结果如下。请进行统计分析。
原料
种类
(A)
温 度(B)
30℃
35℃
40℃
1
41
49
23
25
11
13
25
24
6
22
26
18
2
47
59
50
40
43
38
33
36
8
22
14
18
3
35
53
50
43
38
47
44
55
33
26
29
30
解:要用宏进行双因素方差分析,原始数据必须排列成以下形状(区域设为A2:D14):
原料1
原料2
原料3
温度30
41
47
35
49
59
53
23
50
50
25
40
43
温度35
24
43
38
25
38
47
13
33
44
11
36
55
温度40
6
8
33
22
22
26
26
14
19
18
18
30
计算步骤为:
1°把原始数据按上表形式输入Excel:每个因素A的水平各占一列,因素B的各水平在这一列中依次排列,相同处理的各重复要排在一起。
2°点击“工具”,“数据分析”命令,在对话框中选定输入区域为“A2:D14”,每一样本的行数为“4”(即重复数为4),输出区域为“A16”,然后点击“确定”。出现的结果如下表。
3°结果解释:在出现的SUMMARY表中有因素B各水平分别列出的重复数,和,平均值,方差等数值;最后的总计表中也有因素A各水平的相应数值;在方差分析表中,列出样本(即因素B)、列(即因素A)、交互(交互作用)、内部(即误差)的SS(平方和)、df(自由度)、MS(均方)、F(统计量)、P-Value(尾区概率)、F crit(F统计量的分位数)等数据。把各尾区概率与α相比,大于α时接受H0,否则拒绝。
表11. 例11采用宏计算的输出结果
方差分析:可重复双因素分析
SUMMARY
原料1
原料2
原料3
总计
温度30
计数
4
4
4
12
求和
138
196
181
515
平均
34.5
49
45.25
42.91667
方差
158.3333
62
64.25
118.8106
温度35
计数
4
4
4
12
求和
73
150
184
407
平均
18.25
37.5
46
33.91667
方差
52.91667
17.66667
50
179.9015
温度40
计数
4
4
4
12
求和
72
62
108
242
平均
18
15.5
27
20.16667
方差
74.66667
35.66667
36.66667
66.69697
总计
计数
12
12
12
求和
283
408
473
平均
23.58333
34
39.41667
方差
142.9924
242.1818
125.3561
方差分析
差异源
SS
df
MS
F
P-value
F crit
样本
3150.5
2
1575.25
25.67567
5.67E-07
3.354131
列
1554.167
2
777.0833
12.66601
0.000132
3.354131
交互
808.8333
4
202.2083
3.29588
0.025322
2.727766
内部
1656.5
27
61.35185
总计
7170
35
把三个P-value与手工计算结果相比,可见它们是相同的。
4°若模型不是交叉分组固定模型,则应进行如下计算:
a)有随机因素:利用各因素的MS(均方)数据,采用与例11同样的公式重新计算尾区概率,并与α相比。
b)系统分组:若A为一级因素,则把B(即样本)与交互的SS(平方和)相加,df(自由度)也相加,令它们相除为MSB,然后用与例9同样的公式重新计算尾区概率。若B为一级因素,则合并列与交互。得到尾区概率后再与α相比。
例14. 假设例11中原料(因素B)为随机因素,温度(因素A)为固定因素,则可进行以下补充计算:由于为混合模型,只需重算PA:在空格中输入:
“= fdist(1575.25 / 202.2083, 2, 4)”
回车后,显示数字0.041732。由于这一数字小于0.05,可认为温度间差异显著,但未达极显著。
上述公式中的数字也可用它们的位置代替。
例15. 假设例11中不同原料需用不同的温度水平,即应选用系统分组模型,且原料为一级因素,温度为二级因素。此时应进行以下补充计算:
1o 计算温度的平方和与自由度:在空格F3,F4分别输入:
“= 3150.5 + 808.8333”
“= 2 + 4”
上述数字分别为宏输出的样本平方和、交互平方和、样本自由度、交互自由度。也可用它们的位置代替。
2o 计算温度的均方:在空格F5输入:
“= F3 / F4”
3o 计算温度的F统计量对应的概率:在空格F6中输入:
“= fdist(F5 / 61.35185, F4, 27)”
回车后,显示数字4.1122(10(6。由于它小于0.01,应认为温度间的差异达到极显著。上式中的数字61.35185和27分别为宏输出内部均方和自由度。
从上面的过程可知,利用宏确实可以大大化简计算过程,不过要注意分析模型的类型,必要时进行所需的补充计算。
§3 回归分析
一、统计知识复习。
一元线性回归的统计模型为:yi = ( + (xi + (i。
目的:求出参数(,(的估计值a, b。
方法:最小二乘法。即令残差达到最小。
结果:
其中
2. 对回归方程进行统计检验的方法:
1°对回归系数b, a作t检验:
tb = b / Sb
, ta = a / Sa
一般只对b作检验,自由度均为n - 2。
2°对相关系数r进行统计检验:
检验方法一般为查表。
3°方差分析:
上单尾检验。
这三种检验实际是等价的,只要采用一种即可。
3. 预测值的置信区间:
线性回归的用途之一是预测,即对新的x,计算。新的x取值最好接近,至少不得超出各xi的变化范围。
1°条件均值的置信区间:
点估计:a+bx0
方差:
95%置信区间:
2°下一次观察值y0的置信区间:
点估计:a+bx0
方差:
95%置信区间:
4. 多元线性回归:
原理完全相同,仍采用最小二乘法;自由度有变化。设数据组数为n,自变量个数为m,则回归平方和SSR的自由度为m,残差平方和SSe的自由度为n – m – 1。
二、有关内部函数介绍:
Excel提供的有关回归分析的内部函数有:Linest, Intercept, Slope, Steyx Trend, correl,以及用于指数回归的Logest。注意当这些函数中需要输入因变量和自变量时, 都是因变量在前, 自变量在后。这与一般先x后y的习惯不同。
1. Linest(y数组,x数组,c, s)
其中C取值为true,1或省略,则函数计算截距a;若为False或0,则函数强制a=0。S取值若为true或1,则返回全部统计值;若为False,0或省略,则只返回a, b。
这一函数可用于多元回归,其输入数据以数组形式提供,因变量y只占1列,m个自变量占m列,数据组数(即行数)为n。它的输出数据也是一个数组,为m+1列,5行。输出数据的排列方式为:
bm bm-1 … b2 b1 a
Sem Sem-1 … Se2 Se1 Sea
r2 Sey
F df
SSR SSe
各符号的意义为:
b1, b2, … bm:分别为自变量x1, x2, … xm的回归系数;
a:截距;
Sei:bi的子样标准差;
Sea:a的子样标准差;
Sey:σ的估计值,计算公式为,即;
r2:相关系数的平方。注意平常查表使用的为r;
F:检验y与全体x相关性的统计量。其自由度为:(m, df);
df:SSe的自由度,为n – m – 1;
SSR:回归平方和,自由度为m;
SSe:剩余平方和,又称残差平方和。
由于此函数输出为数组,必须按数组函数方法输入,步骤为:
1°选定输出数据所占的区域;
2°输入公式,例如“=Linest(A2:A10, B2:D10,1,1)”
3°左手按住Ctrl+Shift键,右手再按Enter键。
由于输出数组为一整体,其中任一数字均不能被单独修改。若要修改公式,可将光标移入这一输出区域的任一单元格,则编辑区都会出现公式。修改后再按Ctrl+Shift+Enter即可。若要删除,则可在光标移入区域后,按Ctrl + “/” 选定整个区域,再删除或移动。
对方程进行统计检验的方法:
1°对bi作检验:H0:βi = 0, HA:βi ≠ 0。在一空格输入:
“= Tdist(ABS(bi/Sei),df, 2) ↙”
把返回数字与α比较,大于α则接受H0,否则拒绝。
对a的检验与上述步骤相同。
2°对y与全体x的相关性作检验:在空格中输入:
“= Fdist(F, m, df) ↙”
把返回数据与α比较,大于α则接受H0,即认为y与全体x无关;否则认为相关。
2. 若只需要斜率,截距,相关系数,σ的估计值这几个数中的一个,则可分别采用以下函数:
Slope(y数组,x数组):返回斜率;
Intercept(y数组,x数组):返回截距;
Correl(y数组,x数组):返回相关系数r (不是r2);
Steyx(y数组,x数组):返回σ的估计值,其公式为:。
这几个函数的共同特点是它们只返回单独一个数字,因此可用于更复杂的计算公式中。
3. 若需要预测新的x所对应的y预测值,可用函数:
Trend(y数组,x数组,新x数组,C)
其中C为参数,其取值决定函数如何计算截距。当c取值为1(或逻辑值“true”)时,计算a;当c取值为0(或逻辑值“False”)时,令a=0。
使用这一函数还应注意以下几点:
1°它返回一个数组,因此应先指定返回区域,输入公式后,同时按“Ctrl+Shift+Enter”三个键。
2°若省略新x数组,函数采用原来的x数组计算y的预测值;若连原来的x数组也省略,函数自动认为1,2,3,…… n为自变量。
3°新x数组应与原x数组有相同的列数,即自变量个数相同;但可有不同的行数,即可有不同的数据组数。
4. 指数回归:函数为:
Logest(y数组,x数组,c,s)
这一函数的输入变量与前面介绍的线性回归函数Linest完全相同,输出信息与使用方法也相同,只是Logest函数的回归公式为:
它实际是把上式两边取对数后按线性回归来做,返回的统计量都是线性化(即取对数)后的数据的统计量。
5. 预测值置信区间的建立:
条件均值的95%置信区间:
其理论公式为:
其中,可用Trend(y数组,x数组,x0, c)计算;
分位数t0.975(n–2) = Tinv(0.05, n–2);
在Linest返回的数值中,有,因此可用以下公式计算置信区间的上下限:
Trend(y数组, x数组, x0, c)±Tinv(0.05, n-2)
下一次观察值y0的95%置信区间:
理论公式为:
与上类似,可得计算公式为:
Trend(y数组,x数组,x0,c)±Tinv(0.05, n-2)*
三、计算步骤:
例16. (即本书例5.6) 江苏武进县测定1959-1964年间3月下旬至4月中旬平均温度累积值x和一代三化螟蛾盛发期y的关系如下表(盛发期以5月10日为起算日):试作回归分析。
年代
1956
1957
1958
1959
1960
1961
1962
1963
1964
累积温x
35.5
34.1
31.7
40.3
36.8
40.2
31.7
39.2
44.2
盛发期y
12
16
9
2
7
3
13
9
–1
解: 表12. 线性回归的输出结果
累积温x
盛发期y
b
-1.09962
48.54932
a
31.7
9
Sb
0.271567
10.12779
Sa
31.7
13
r2
0.700801
3.265989
MSe1/2
34.1
16
F
16.3958
7
df
35.5
12
SSR
174.8888
74.66678
SSE
36.8
7
39.2
9
PF
0.004876
40.2
3
Pt
0.004876
40.3
2
44.2
-1
具体计算步骤为:
1°把原始数据输入A3:B11区域,如表12。
2°用鼠标选定E2:F6区域。输入公式:
“=Linest(B3:B11, A3:A11, 1, 1)”
3°同时按下“Ctrl+Shift+Enter”键,返回数据的排列为:
E2:b F2: a
E3:Sb F3: Sa
E4:r2 F4: Sey=
E5:F F5: df
E6:SSR F6: SSe
为清楚起见,我们在D2:D6,G2:G6中标上了各数据的统计意义。
5°对回归方程进行F检验:在E8输入:
“=Fdist (E5, 1, F5)”
回车后,显示数字0.004876。把E8的返回值与α比较,若α> 0.05,则认为回归失败,即接受β=0;若α< 0.05,则认为回归成功,即β≠ 0。本例题E8<<0.05,回归是成功的。
为进行比较,我们在E9中给出了对b进行t检验的结果。在E9输入:
“=TDIST(ABS(E2/E3),F5,2)”
回车后,显示的数字与E8是完全相同,即F检验和t检验的尾区概率是完全相同的,说明只进行一种检验即可。
上式中由于Tdist 函数要求输入的统计量值x为正数,故增加了内部函数ABS,其功能是取绝对值。F5为自由度,2表示为双尾检验。
以下步骤是为了画出包括回归线、观测值,以及条件均值和下次观察值95%置信区间的图。数据见表13。由于此表共有140余行,这里显示的仅是前19行。
表13. 回归分析绘图数据(部分)
累积温x
预测值
观测值
均值下限
均值上限
预测值下限
预测值上限
31
14.46104
9.785646
19.13643
5.433224
23.48885
31.1
14.35107
9.729153
18.97299
5.350837
23.35131
31.2
14.24111
9.672384
18.80984
5.268075
23.21415
31.3
14.13115
9.615328
18.64697
5.184935
23.07736
31.4
14.02119
9.557975
18.4844
5.101413
22.94096
31.5
13.91123
9.500315
18.32213
5.017506
22.80494
31.6
13.80126
9.442337
18.16019
4.933211
22.66932
31.7
13.6913
9
9.384028
17.99857
4.848523
22.53408
31.8
13.58134
13
9.325378
17.8373
4.763441
22.39924
31.9
13.47138
9.266373
17.67638
4.67796
22.26479
32
13.36141
9.207
17.51583
4.592077
22.13075
32.1
13.25145
9.147246
17.35566
4.505789
21.99711
32.2
13.14149
9.087097
17.19588
4.419092
21.86389
32.3
13.03153
9.026538
17.03652
4.331984
21.73107
32.4
12.92157
8.965553
16.87758
4.244461
21.59867
32.5
12.8116
8.904127
16.71908
4.156519
21.46669
32.6
12.70164
8.842243
16.56104
4.068156
21.33513
32.7
12.59168
8.779884
16.40347
3.979368
21.20399
6°在I3:I143中填充31~45的数值,间隔为0.1。
7°计算预测值(也是回归线):标记J3:J143,输入公式:
“=Trend (B3:B11, A3:A11, I3:I143)”
同时按下“Ctrl+Shift+Enter”三键。在J列出现预测值。其中J3的值为14.46140。
8°在K列中适当地方入原始观察y值,使它与I列的X值对应。由于原始数据中x等于31.7的有两组,我们把它们分别放在对应31.7和31.8的地方。
9°计算条件均值置信区间下限:在L3中输入公式:
“=$J3-TINV(0.05,$F$5)*SQRT($F$4^2/($F$5+2)+
$E$3^2*($I3-AVERAGE($A$3:$A$11))^2)”
回车后,显示数字9.785646。把L3复制到L4:L143。
式中加了“$”号的地址是为了在复制过程中使它不改变,也可用相应单元中的数值代替,例如$F$5可换为7。
10°计算条件均值置信区间上限:由于这几个置信区间的公式大同小异,故可把L3的公式复制到M3:O3,然后加以修改。
把L3复制到M3:O3。
把光标移到M3,在编辑栏中把$J3后边的“-”号改为“+”号,并按“Enter”键。显示数字为19.13643。
11°计算下次观察值下限:把光标移到N3后,在编辑栏中原公式的“$F$4^2”后边加上:“*($F$5+3)”,完整公式为:
“=$J3-TINV(0.05,$F$5)*SQRT($F$4^2*($F$5+3)/($F$5+2)+
$E$3^2*($I3-AVERAGE($A$3:$A$11))^2)”
按“Enter”键。显示数字为5.433224。
12°计算下次观察值上限:把N3复制到O3,光标移到O3,并把编辑栏中的公式中开始处的“$J3-”改为“$J3+”,其他不变,按“Enter”键。显示数字为23.48885。
把L3:O3复制到L4:O143。
13°在L2:O2中加上各列数据名称,如“累积温”,“预测值”,“观测值”,……等,如表13。
以上就完成了全部计算。下面来画图。
14°把光标移入空格,例如Q3。
15°按工具栏中的“图表向导”钮,或按“插入”,“图表”。然后按照图表向导的指引,选择所需图类型、数据区域等,并对图形中不满意的部分进行修正。由于这一部分不属于统计内容,而是Excel的基本使用技巧,在此不再详细介绍,有兴趣同志可参阅有关书籍。最后图形见图1。
图1. 例16的回归线及置信区间图
§4 Excel中常用统计函数简介。
Excel函数名主要由字母组成, 输入时大小写均可。函数指南中共有9大类函数,其中列入“统计”类的共71个,而我们常用的约有50余个,可分为以下三大类:
一、对数据进行统计处理的函数:
这类函数主要用于从原始数据计算一些常用统计量,如均值,方差,相关系数等。主要有:
1. 求平均数:Average (x1, x2, … xn)
函数的输入可以是数值,也可以是存贮数据的地址或区域。其他函数输入类似,不再重复。
求相关系数:Correl(数组1,数组2)
数组可输入地址,也可直接输入数据。若直接输入数据, 每个数组要用一个大括号“{}”围住。以下其他函数输入数组的方法相同。
3. 计算参数组中的数字个数:Count (x1, x2, …xn)
4. 计算参数组中非空单元格数:CountA (x1, x2, …xn)
5. 计算协方差:Covar(数组1,数组2)
6. 计算:Devsq(x1, x2, … xn)
7. 计算最大值:max(x1, x2, … xn)
8. 计算中位数:median(x1, x2, … xn)
9. 计算最小值:min(x1, x2, … xn)
10. 计算出现频率最高的值:mode(x1, x2, … xn)
11. 计算组合数;Permut (n, m)
12. 计算一个数在数列中的排序:Rank(数,数列,次序)
其中次序为一个参数,它取值为0或被省略时,按递减排序;取值为1时递增排序。
13. 标准化(公式为):standardize (x,μ,σ)
14. 计算样本标准差:stdev (x1, x2, … xn)
15. 计算样本方差:Var(x1, x2, … xn)
16. 计算:Avedev (x1, x2, … xn)
17. 计算频率分布:Frequency(数据,间隔点)
它返回的数据以垂直数组形式给出,其个数比间隔点数多1,计算时会忽略空白与文字。
18. 计算正数的几何平均数:Geomean(x1, x2, … xn)
19. 计算调和平均数H(公式为):Harmean(x1, x2, … xn)
20. 找出第K个最大值:Large(数组,K)
21. 计算服从指定离散分布的随机变量落入某一区间的概率:
Prob (数组x, 数组p,下限,上限)
其中x为离散分布的取值,p为对应的概率。区间上限可省略,省略后只计算下限一个点的概率。
22. 计算样本峭度(峰度):Kurt(x1, x2, … xn)
其计算公式为:
23. 计算样本偏度:Skew(x1, x2, … xn)
公式为:
注:计算峭度与偏度函数的输入最多为30个数字。
24. 取第K个最小值:small(数组,K)
25. 计算总体标准差(公式为):Stdevp(x1, x2, … xn)
26. 计算总体方差(公式为):Varp(x1, x2, … xn)
二、计算常用统计分布的函数:
注意这些函数中正态分布、指数分布、离散分布一般计算分布函数(即P(X<x)),χ2分布、t分布、F分布一般计算尾区概率(即P(X>x)),而以“dist”结尾的和以“inv”结尾的互为反函数。函数的输入变量中x一般表示统计量的取值;p表示概率;df表示自由度。
1. 计算χ2分布的单尾概率:Chidist (x, df)
2. 计算χ2分布的单尾分位数:Chiinv (p, df)
3. 计算F分布的单尾概率:Fdist (x, df1, df2)
df1为分子自由度,df2为分母自由度。
4. 计算F分布的单尾分位数:Finv (p, df1, df2)
5. 计算超几何分布概率:Hypgeomdist (k, n, M, N)
其中k:样本中成功数; n:样本含量;
M:总体中成功数; N:总体中个体数。
计算公式为:
6. 计算负二项分布概率:Negbinomdist (x, r, p)
其中x:失败次数;r:成功次数;p:成功概率。
计算公式为:
7. 计算指数分布:Expondist (x, λ, c)
其中c取值为true或1,则计算累积分布;为False或0,则计算密度函数;
累积分布公式为:F(x, λ)=;密度函数公式为:f(x,λ)=
8. 计算二项分布概率:Binomdist (x, n, p, c)
其中x:成功次数;n:总次数;p:成功概率,
c为true或1:累积概率;c为False或0:成功x次的概率。
计算公式为:
9. 计算正态分布概率:Normdist (x, μ,σ,c)
其中μ:数学期望;σ:标准差;
c为true或1:计算分布函数;c为False或0:计算密度函数。
10. 计算正态分布分位数:Norminv (p, μ,σ)
11. 计算标准正态分布分布函数:Normsdist (x)
12. 计算标准正态分布分位数:Normsinv (p)
13. 计算泊松分布概率:Poisson (x, mean, C)
其中mean为均值;
C为true或1:计算[0, x]中的累积概率;C为False或0:计算x点的概率。
14. 计算t分布尾区概率:Tdist (x, df, tails)
其中x必须大于0;df为自由度;tails为1:单尾;为2:双尾。
15. 计算t分布双尾尾区的分位数:Tinv (p, df)
其中p:双尾尾区概率。
16. 计算累积二项分布的逆:Critbinom (N, p, α)
其中N:总实验次数;p:成功概率;α:临界值。
函数返回值为累积概率大于等于临界值的最小成功次数x。
三、直接进行某些统计检验的函数:
这些函数大部分已在前边的例题中出现过,为方便查阅,我们再次把它们列在下面:
1. 列联表独立性检验:Chitest(观测值,期望值)
其中观测值、期望值应存在两个行列数相同的矩形区域内,函数返回Pearson统计量所对应的分布的上单尾概率,自由度为(r-1) (C-1), r, C为矩形区域的行,列数。如为1行或1列,则自由度为r-1或C-1。
2. 计算标准正态分布均值置信区间宽度的一半: Confidence (α,σ,n)
其中α:显著性水平;σ:标准差;n:样本含量。
3. F检验:Ftest(数组1,数组2)
用于检验两组数据方差是否相等,返回双边尾区概率。
4.求线性回归截距:Intercept(y数组,x数组)
5.线性回归:Linest(y数组,x数组,C,S)
C为true,1或省略:计算a;为False或O:令a=0;
S为true或1:返回统计值;为False或O:只返回a,b。
y数组为单列,x数组则可有多列,相当于多元回归。两数组应有相同的行数。
返回统计值的排列为:(设x数组有m列)。
bm bm-1 … b2 b1 a
Sem Sem-1 … Se2 Se1 Sea
r2 Sey
F df
SSR SSe
其中各Se为标准误差,于上一列的估计值对应;
6. 指数回归:Logest(y数组,x数组,C,S)
回归模型:
把上式两边取对数后按线性回归做。参数及返回数值排列同线性回归,返回统计量为线性方程统计量。
7. 计算线性回归斜率:slope(y数组,x数组)
8. 计算线性回归中:Steyx(y数组,x数组)
9. 计算线性回归预测值:Trend(y数组,x数组,新x数组,C)
先用y,x做线性回归,再利用回归方程预测新x所对应的y值。
C为true或1:计算截距a;为False或0:令a = 0。
若新x省略则用旧x数组,若旧x也省略则用1,2,… n为自变量。新、旧x应有相同的列数。
10. t检验:Ttest(数组1,数组2,tails, type)
tails为1:单尾;为2:双尾。
Type为1:配对检验;为2:方差相等;为3:方差不等。
用于检验两数组均值是否相等。返回值为尾区概率。
11. Z检验(本书中一般称为U检验):Ztest (数组,μ,σ)
用于已知标准差σ的情况下检验数组期望是否为μ。若σ未知一般应使用t检验。
公式为: 1-Normsdist
其中为数组均值,n为数组中数据个数,σ为总体标准差,若省略则用子样标准差S代替。返回值为单边尾区概率。
12. 计算相关系数:Correl(数组1,数组2)
注意:1°所有与回归有关的函数中,均把因变量y放在前边,自变量x放在后边,不要搞错次序。
2°返回值为一个数组时要按数组方式输入,即先选定输出区域,输入公式后,同时按下“Ctrl+Shift+Enter”键。