3.1 假设检验的基本概念
例3.1 已知小麦亩产服从正态分布,传统小麦品种平均亩产800斤,现有新品种产量未知,试种10块,每块一亩,产量为:
775,816,834,836,858,863,873,877,885,901
问:新产品亩产是否超过了800斤?
假设检验就是概率意义上的反证法。要证明命题H1:,可以首先假设H0:。本体中容易计算样本均值超过800了,有没有可能超过800的原因是由于抽样的随机性引起的?是否总体均值根本没有变化?我们看如下的统计量:
容易看出,如果新品种确有增产效应,应偏大,不利于H0,取,查表求临界值,使得,即构造不利于H0,有利于H1的小概率事件,如果在一次试验中该小概率事件发生了,就有理由拒绝H0,认为H1成立。
严格逻辑意义上的反证法思路如下:欲证H1成立,先假设其否命题H0成立,然后找出逻辑意义上的矛盾,从而推翻H0成立,严格证明H1成立。假设检验的思路类似,只不过引出的不是矛盾,而是小概率事件在一次实验中发生。
我们称想要证明的命题H1为备择假设,对立的命题H0称为原假设,面对样本,我们必须表态是接受原假设还是拒绝原假设,这有可能出现两类错误。如果客观上原假设的确成立,面对样本的异常我们拒绝了原假设,这种"以真为假"的错误我们称为第一类错误,发生的概率用表示;如果客观上备择假设成立,我们却接受了原假设,这种"以假为真"的错误我们称为第二类错误,用发生的概率用表示。假设假设检验一般首先控制第一类错误,即:当我们拒绝原假设时有比较充足的理由,犯错误的概率不超过预设的,称为显著性水平。常用的显著性水平有
这种预设显著性水平的假设检验也称为显著性检验,以后我们提到的假设检验都是显著性检验。对于显著性检验,当接受原假设时,可以认为是拒绝的证据不足。
对于例3.1的问题,取,当时拒绝原假设。这里称为检验统计量,所确定的的取值范围称为拒绝域。
x=[775,816,834,836,858,863,873,877,885,901];
T=(mean(x)-800)/(std(x)/sqrt(9)),
ta=tinv(0.95,9),
计算结果T=4.1669>ta=1.8331,故拒绝原假设,认为确有增产。
之所以查表求临界值,是因为当初计算机及数学软件尚未普及,人们利用稀有的计算机资源计算出了一些关键的临界值,供没有计算机的人们膜拜使用。因此上述解题套路是几乎所有教科书上使用的方法,不妨称为"查表法"。
由于计算机及数学软件的普及,统计方法的使用套路也应该更新,如果写作业写论文都用计算机打字,真正数学计算反而要翻书本查表,怎么看也都很滑稽。
其实,Matlab可以计算常用分布在任意一点的分布函数的值,例如对于上述T=4.1669,可以直接计算分布函数在该点的值:
p=tcdf(T,9)
计算结果为0.9988,超过了。或者计算出1-p=0.0012,小于我们预设的显著性水平。面对0.0012这个值,我们拒绝了原假设,就是使用了概率意义上的反证法。我们可以做一个比喻:张三每天上网游戏,期末考试肯定不及格,我们说:"要想张三及格,除非明天太阳从西边出来"。这里原假设是"及格",备择假设"不及格"是我们想证明的东西。其等价的逆否命题是:因为明天太阳不会从西边出来,所以张三一定不及格。这是我们说话的内含逻辑。"太阳从西边出来"是不可能事件,我们使用的是语文上"夸张"的修辞方法以表达对张三的极度鄙视。
现在,面对新品种亩产数据,我们的结论是:要说没有增产效应,除非明天下大雹子。这里没有"夸张",因为1-p=0.0012大约为千分之一,是类似于不可能事件的极小概率事件,和明天下大雹子一样罕见(大约三年才得一见)。我们计算出来的1-p越小,说明备择假设成立的证据越充足。
几十年前,对于自由度为9的分布,我们只能将
1.3830,1.8331 ,2.2622,2.8214
等少数几个值印在书上,现在我们可以计算p=tcdf(T,9)在任意一点分布函数的值。
3.2 正态总体参数的假设检验
3.2.1 单正态总体均值的假设检验
设为来自正态总体简单随机样本,为我们关心的已知的值,原假设为:
H0:
(1)方差已知情形
此时,检验统计量为,H0成立时,依据备择假设的不同提法,分三种情况分别给出拒绝域。
1)双侧检验 备择假设H1: 拒绝域:
这种情形我们关心的是总体均值是否发生了变化,增多减少都是我们同等关注的。例如要研究某种药物的副作用,是否引起血压的变化,变大变小都是副作用,如果实验证明了确有副作用,就该停产或慎用。
2)单侧检验(右侧) 备择假设H1: 拒绝域:
这种情形我们关心的是总体均值是否有增加效应,例如小麦亩产。无增产效应或者减产都是我们不希望看到的,我们希望证明的是增产了。
3)单侧检验(左侧) 备择假设H1: 拒绝域:
这种情形我们希望看到总体均值变小了。每匹布上疵点的个数。新工艺后是否有减少。
(2)方差未知情形
原假设H0:
此时,检验统计量为,H0成立时,依据备择假设的不同提法,分三种情况分别给出拒绝域。
1)双侧检验 备择假设H1: 拒绝域:
2)单侧检验(右侧) 备择假设H1: 拒绝域:
3)单侧检验(左侧) 备择假设H1: 拒绝域:
其实,上一章中区间估计与这里的双侧检验本质上是相同的:区间套中接受原假设,没套中则拒绝原假设。只不过检验统计量的计算更简单些。类似于单侧检验,也可以有单侧区间估计。
3.2.2 单正态总体方差的假设检验
设为来自正态总体简单随机样本,为我们关心的已知的值,原假设为H0:,检验统计量为
当H0成立时,,由此可查临界值表,构造拒绝域。
(1)双侧检验 此时备择假设为H1:,也就是说,我们希望通过样本找到总体方差比较有明显变化的证据,无论变大变小都是我们希望证明的。
此时取临界值与,使得,,拒绝域为:(方差变小了),或者(方差变大了)。
当已经赋值的时候,执行如下Matlab命令可得到临界值。
a=0.05, n=20, c1=chi2inv(a/2,n-1), c2=chi2inv(1-a/2,n-1),
(2)单侧检验(右侧) 此时备择假设为H1:,也就是说,我们关心的是方差是否变大了。此时临界值为满足,可用
c=chi2inv(1-a,n-1)
(3)单侧检验(左侧) 此时备择假设为H1:,也就是说,我们关心的是方差是否变小了。此时临界值为满足,可用
c=chi2inv(a,n-1)
3.2.3 两正态总体均值的假设检验
设为来自正态总体的简单随机样本,为来自正态总体的简单随机样本,且两样本独立。为比较两个总体的期望,提出如下原假设:
H0:
与前面类似,备择假设有双侧、单侧(左侧、右侧)等提法。
(1)方差已知情形
此时检验统计量为,当H0成立时服从标准正态分布,临界值,含义及计算方法同前。
1)双侧检验 H1:,拒绝域:
2)右侧检验 H1:,拒绝域:
3)左侧检验 H1:,拒绝域:
(2)方差未知但相等情形
此时原假设仍为H0:,备择假设同样有三种提法。检验统计量为:
当H0成立时,由此得临界值,。
1)双侧检验 H1:,拒绝域:
2)右侧检验 H1:,拒绝域:
3)左侧检验 H1:,拒绝域:
3.2.4 两正态总体方差的假设检验
设为来自正态总体的简单随机样本,为来自正态总体的简单随机样本,且两样本独立。为比较两个总体的方差,提出如下原假设:
H0:
与前面类似,备择假设有双侧、单侧(左侧、右侧)等提法。此时检验统计量为,当H0成立时,,在Matlab中,如果m,n已经赋值,例如m=8,n=10则
c1=finv(0.025,7,9),c2=finv(0.975,7,9)
分别给出了时的两个临界值,双侧检验的拒绝域为或。
c3=finv(0.05,7,9)
给出了左侧检验临界值,时拒绝原假设,认为备择假设H1:成立。
c4=finv(0.95,7,9)
给出了右侧检验临界值,时拒绝原假设,认为备择假设H1:成立。
3.2.5 大样本非正态总体均值的假设检验
设为来自非正态总体的简单随机样本,设总体均值与总体方差有限,原假设
H0:
此时可以将作为近似的检验统计量,当样本容量很大时(例如100),由中心极限定理知H0成立时近似服从标准正态分布,可以仿照3.2.1小节中的算法检验如下三个备择假设:
H1:; H1:; H1:
设为来自非正态总体的简单随机样本,为来自非正态总体的简单随机样本,且两样本独立。两个总体有有限的均值与方差,均值为与,为比较两个总体的期望,提出如下原假设:
H0:
与前面类似,备择假设有双侧、单侧(左侧、右侧)等提法。此时可以将
近似作为检验统计量,当两个样本容量都很大时(例如100),由中心极限定理知H0成立时近似服从标准正态分布,可以仿照3.2.3小节中的算法检验如下三个备择假设:
H1:; H1:; H1:
3.3 三个常用的非参数检验
大样本情形下,对于非正态总体,可以利用中心极限定理近似用标准正态分布进行假设检验。小样本情形,若总体不是正态分布的,可以使用非参数检验的方法。非参数检验的效率稍差,但适应各种总体类型,应用范围较广。
3.3.1 符号检验
例3.2 已知原来工艺下生产的某种灯泡的中位数为800小时,现改进生产工艺,试产10只灯泡,实验得到每只寿命为:
775,816,834,836,858,863,873,877,885,901
问:新工艺生产的灯泡寿命中位数是否超过了800小时?
H0:
一般情况下,灯泡寿命不是正态分布的,不能用例3.1的方法。符号检验使用的是计数统计量,先设
则有
即记录样本点中大于800的个数。若H0成立,应该大约占样本容量的一半左右,若异常的大,说明备择假设H1:成立。
H0成立时,,可以利用二项分布构造拒绝域:
使得若H0成立时
,
利用二项分布的分布律可以计算出临界值,用如下Matlab函数文件计算。
function t=bt(n,a)
SS=2^n*a;
S=0;
c=1;
k=n+1;
while S<=SS
k=k-1;
S=S+c;
c=c*k/(n-k+1);
end
t=k+1;
以上自定义函数扩展了Matlab的功能,可以替代教科书上的"符号检验临界值表",并且可以使用任意的n及。
在例3.2中,,对于,使用命令t=bt(10,0.05)可以得到临界值9,临界值9,落在拒绝域内,故拒绝原假设,认为新工艺生产的灯泡寿命中位数超过了800小时。
只要去代替,也可以进行双侧符号检验。
例3.3 20个品酒师对A、B两种白酒进行品尝,有17个品酒师认为A品质好,3个品酒师认为B品质好,在的显著性水平下,检验两种白酒品质是否存在差异?
解 ,设原假设为
H0:两种白酒品质无差异
令表示认为A品质好的品酒师的人数,则H0成立时应该在10左右取值,如果值异常大,或者异常小,都说明两种白酒品质有差异。取临界值与,使得,,由于关于对称,故有,因此可用水平为的单侧检验求出临界值。命令
t2=bt(20,0.05/2)
得到,因此,此例中拒绝域为
,或者
落在拒绝域内,可以认为两种白酒品质有显著差异。
有些教科书中没有0.025的临界值,而我们的函数bt.m扩展了功能。
Matlab中有自带的SIGNTEST函数,可以直接用于符号检验。默认的检验是双侧的。
对于配对实验的两总体均值检验问题,也可用符号检验。
3.3.2 Wilcoxon秩和检验
我们要研究的问题是两总体均值的假设检验,设
,
要检验第二个总体是否有增加效应,即检验如下问题:
H0: H1:
Wilcoxon秩和检验的方法是:将两个样本混合为
混合之后样本容量为,每个样本点在样本中从小到大排列的名次称为该样本点的秩,用表示在混合样本中的秩,表示在混合样本中的秩,检验统计量为
例如诸为 1.1,3.3,5.5,7.7,诸为2.2,4.4,6.6,以下列表给出混合样本及秩
混合样本 | 1.1 | 3.3 | 5.5 | 7.7 | 2.2 | 4.4 | 6.6 |
秩 | 1 | 3 | 5 | 7 | 2 | 4 | 6 |
则。
若H0成立,则的值应该适中。注意到每个秩序的平均值为,故H0成立时,,的值在此值附近应该是正常的。若的值异常偏大,说明第二个总体确有增加效应。利用matlab自身的函数
p = ranksum(X,Y)
可以进行双侧的秩和检验。返回的p值小于给定的则拒绝原假设,认为H1:成立。
H0成立时,可以证明关于对称,要检验H1:,只要判定,并且p = ranksum(X,Y)即可。
自定义rsum函数用于求
function W=rsum(x,y)
[s,t]=size(x);
m=max(s,t);
if t<m
x=x\';
end
[s,t]=size(y);
n=max(s,t);
N=m+n;
if t<n
y=y\';
end
xy=[x,y];
[z,I]=sort(xy);
W=0;
for i=1:N
if(I(i))>m
W=W+i;
end
end
为了求出Wilcoxon秩和检验的临界值,我们给出如下定理,证明参见文献[1]。
定理3.1 在H0成立时,的概率分布为
其中表示从中取个数其和恰为的取法的个数。可用如下初始条件及递推公式计算:
当
自己编程tmnd.m计算如下:
function tmn=tmnd(m,n,d)
N=m+n;
nn=n*(n+1)/2;
NN=n*(2*m+n+1)/2;
if m<0 | n<0 | d<nn | d>NN
tmn=0;
elseif m>0 & n==0 & d==0
tmn=1;
elseif m>0 & n==0 & d>0
tmn=0;
elseif m==0 & n>0 & d==nn
tmn=1;
elseif m==0 & n>0 & d~nn
tmn=0;
else
T=zeros(m,n,NN);
for i=1:m
for k=1:i+1;
T(i,1,k)=1;
end
end
for j=1:n
kk=j*(j+1)/2;
KK=(j+1)*(j+2)/2-1;
for k=kk:KK
T(1,j,k)=1;
end
end
for i=2:m
for j=2:n
s=i+j;
for k=1:d
if k<=s
T(i,j,k)=T(i-1,j,k);
else
T(i,j,k)=T(i,j-1,k-s)+T(i-1,j,k);
end
end
end
end
tmn=T(m,n,d);
end
可以证明,H0成立时,的概率分布关于E=n*(m+n+1)/2对称,我们给出单侧检验临界值的求法,以下自定义函数wr.m,其中输入参数m,n,alpha分别是对照组样本容量、实验组样本容量、检验的显著性水平,而输出值c表示右侧临界值,即满足的最小正整数。
function c=wr(m,n,alpha)
% return the min c such that P(W>=c)<=alpha
NN=n*(2*m+n+1)/2;
nn=n*(n+1)/2;
N=m+n;
E=n*(N+1)/2;
a=1;
for k=1:n
a=a*(N+1-k)/k;
end
Alpha=a*alpha;
k=nn;
P=0;
while P<Alpha
P=P+tmnd(m,n,k);
k=k+1;
end
c1=k-1;
c=2*E-c1;
上述函数可用于右侧检验。若左侧检验,c1=2*E-c即为左侧临界值。若双侧检验,先求出c2=wr(m,n,alpha/2),再由c1=2*E-c2即可。
例3.4 某班级共15名同学,某次英语水平考试,分数如下:
男:53,55,59,65,71,77,81
女:56,62,68,76,84,86,90,96
在显著性水平下,能否认为女生英语水平高于男生?要求采用Wilcoxon秩和检验。
解 注意这是一个单侧检验问题,使用matlab命令:
x=[53,55,59,65,71,77,81]
y=[56,62,68,76,84,86,90,96]
rsum(x,y)
c=wr(7,8,0.05)
上述计算中,注意到rsum(x,y)=78,而临界值为c=78,的值落在拒绝域内,故可拒绝原假设,认为女生成绩显著高于男生。
3.3.3 Wilcoxon符号秩检验
设为来自连续总体的简单随机样本,关于点对称,检验假设
H0: H1:
Wilcoxon符号秩检验统计量为:
其中
,即把依照绝对值由小到大排列,的名次。
H0成立时,,故在此值附近取值说明原假设成立。若异常大,则要拒绝原假设,说明H1:成立。
对于双侧检验问题
H0: H1:
Matlab有自带的函数
p=signrank(x,m)
这里x为样本,m代表,若显著性水平为,则时拒绝原假设。
对于单侧检验,H1:,要拒绝原假设需要同时满足两个条件:条件一,;条件二,p=signrank(x,m)<。为计算,自编函数:
function wp=rpsum(x,m);
n=length(x);
x=x-m;
y=abs(x);
[z,I]=sort(y);
wp=0;
for i=1:n
if x(I(i))>0
wp=wp+i;
end
end
保存了上述函数后,即可进行单侧检验。
例3.5 某班级共15名同学,某次英语水平考试,分数如下:
53,55,59,65,71,77,81,56,62,68,76,84,86,90,96
在显著性水平下,能否认为平均成绩高于60分?要求分别用:
(1)符号检验;(2)Wilcoxon符号秩检验。
解 注意这是一个单侧检验问题:
H0: H1:
使用matlab命令:
x=[53,55,59,65,71,77,81,56,62,68,76,84,86,90,96]
(1)符号检验
注意这里n=15,B=11,利用前面自定义的bt.m函数计算:
t=bt(15,0.05)
得到临界值,B=11<,没有落入拒绝域,故接受H0,认为平均成绩没有明显高于60分。
(2)Wilcoxon符号秩检验
E=n*(n+1)/4,
wp=rpsum(x,60),
计算结果发现wp=106> E=60,满足单侧检验条件一,再计算
p=signrank(x,60)
结果得p=0.0071<2,故拒绝原假设,认为平均成绩明显高于60分。
3.4 检验的功效函数
为了简单起见,我们只讨论位置参数的单侧检验:
H0: H1:
其中为总体的中位数。
对于上述检验,当总体为方差已知正态总体时,有检验;当总体为方差未知正态总体时,有检验;当总体为连续对称总体时,有符号检验及Wilcoxon符号秩检验。自然有一个问题,如何评价不同的检验方法的优劣?
对于相同的样本容量,对于相同的显著性水平,一般比较区间时拒绝的概率,此时为犯第二类错误的概率。不同的检验方法犯第一类错误的概率已经被控制了,具有相同的水平,此时比较时的,小者为好;或者等价地说,比较时的,越大越好。
称,为检验的功效函数。功效大的检验就是好的检验。
以下画出正态总体方差已知时检验的功效函数。H0时,不妨设总体服从标准正态分布,已知,均值用m表示。
以下固定样本容量n=20,固定显著性水平a=0.05,此时检验临界值为
u0=norminv(0.95)=1.6449。当m>0时,检验统计量为
容易计算
=1-normcdf(u0-m*sqrt(20))
以下利用Matlab作图功能画出此时的功效函数。
u0=norminv(0.95)
m=0:0.01:1;
w=1-normcdf(u0-m*sqrt(20));
plot(m,w)
结果如图3-1所示。
图图3-1 n=20,=0.05单侧检验功效函数
请读者自己研究,随着样本容量的增加,功效函数的图形会有怎样的变化?
注意,这是水平为的检验的出发点,类似于百米赛跑,此点是起跑点。如果相同起跑点,随着的增加,功效函数越来越大,对于两条功效函数曲线,在备择假设的范围内大者为佳。
上述功效函数容易得到精确的曲线,稍微复杂的情形,拒绝概率的精确值不易计算,可以使用随机模拟的方法得到功效函数。
例如,要研究t检验的功效函数、符号检验的功效函数、Wilcoxon符号秩检验的功效函数,并与检验的功效函数进行对比。首先固定如下四个因素:
(1)总体分布;
(2)样本容量;
(3)显著性水平a=0.05;
(4)取定
前三条都满足时,三种方法的临界值就完全确定了,拒绝域也完全确定了:
t检验:,拒绝域为t0=tinv(0.95,19)=1.7291;
符号检验:大于0样本点个数,拒绝域t=bt(20,0.05)=15;
Wilcoxon符号秩检验:拒绝域
为评价不同的检验,我们可以分别计算功效函数。这可以采用随机模拟的方法,利用万次随机试验中拒绝的频率近似代替拒绝概率。
以下命令文件保存为p123.m
m=0:0.1:1;
p1=zeros(1,11);
p2=zeros(1,11);
p3=zeros(1,11);
t0=tinv(0.95,19);
b0=15;
w0=150;
s20=sqrt(20);
N=10000;
for mm=1:11
for k=1:N
x=randn(1,20)+m(mm);
T=s20*mean(x)/std(x);
if T>=t0
p1(mm)=p1(mm)+1;
end
B=0;
for i=1:20
if x(i)>0
B=B+1;
end
end
if B>=b0
p2(mm)=p2(mm)+1;
end
wp=rpsum(x,0);
if wp>=w0
p3(mm)=p3(mm)+1;
end
end
p1(mm)=p1(mm)/N;
p2(mm)=p2(mm)/N;
p3(mm)=p3(mm)/N;
mm
end
P=[p1;p2;p3]
plot(m,p1,m,p2,m,p3)
计算结果为
P =
0.0482 0.1105 0.2100 0.3522 0.5304 0.6958 0.8292 0.9132 0.9609 0.9874 0.9953
0.0192 0.0460 0.0864 0.1619 0.2551 0.3775 0.5247 0.652
请发表评论