假设检验中的第二类错误

weibin

2018/10/18

一、文献综述

1.假设检验的提出与发展

假设检验是将具体的事件转变为抽象的数学语言,是一种有重要应用价值的统计推断形式。 1900年由卡尔·皮尔逊(K.Pearson)提出来的拟合优度 χ2\chi^2 检验, 主要用于估量数据与模型拟合的程度如何,并由此奠定了假设检验的思想基础。 20世纪20年代费歇尔(R.A.Fisher)提出的显著性检验对其进行了细化。 最终由奈曼(J.Neyman)和艾·皮尔逊(E.S.Pearson)提出了较为完整的假设检验理论。 Neyman-Pearson关于假设检验的理论(NP理论),是建立在概率的频率解释基础之上的、关于假设检验的一套形式比较完美的数学理论。 但NP理论不是另起炉灶,而是对Pearson和Fisher的工作的继承和发展。

2.Pearson的思想,拟合优度检验

为了更好地估量拟合程度如何,K.Pearson提出了拟合优度检验:

先设 XX 只取有限个不同的值 $a_1,\ldots,a_r$ 。 理论分布 FF (完全已知的情况下)集中在 aia_i 点的概率记为 pip_i 。检验假设:

H:P(X=ai)=pi(i=1,,r)H:P(X=a_i)=p_i(i=1,\ldots,r)

其中, pi>0p_i>0 已知, i=1rpi=1\sum_{i=1}^rp_i=1

νi\nu_iX1,,XnX_1,\ldots,X_n 中等于 aia_i 的个数, νi\nu_i 称为 aia_i 的观察频数。 有 i=1rνi=n\sum_{i=1}^r\nu_i=nnpinp_i 称为 aia_i 的理论频数,意即当 XX 的分布确为 FF 时, νi\nu_i “在理论上"应取的值。 K.Pearson引进Pearson χ2\chi^2 统计量,以反映样本与理论分布的偏离:

k=k(X1,,Xn;F)=i=1r(νinpi)2npik=k(X_1,\ldots,X_n;F)=\sum_{i=1}^r\frac{(\nu_i-np_i)^2}{np_i}

K.Pearson在给出 χ2\chi^2 统计量的同时,证明了一个极限定理:若假设 HH 为真,则当样本大小 nn\to\infty 时,统计量 kk 的均匀分布收敛于 χr12\chi_{r-1}^2 ,即自由度为 r1r-1χ2\chi^2 分布。

由此可近似算出概率 P(kk0)P(k \geq k_0) 。 这个数称为拟合优度。因为这个数愈大,则产生像 k0k_0 这么大或更大的偏离机会愈多,因而实际得到 k0k_0 这么大偏离这件事并不稀奇。 一般地,会先给定一个值 α\alpha (如0.01,0.05等),从 χ2\chi^2 分布表上查出满足条件 P(χr12χr12(α))=αP(\chi_{r-1}^2 \geq \chi_{r-1}^2(\alpha))=\alpha 的值 χr12(α)\chi_{r-1}^2(\alpha) ,然后视 k0>χr12(α)k_0>\chi_{r-1}^2(\alpha) 或否,决定否定或接受 HH .这就是Pearson的拟合优度检验。

如果我们需要通过实践去检验一个理论、假说等,只要能设计一种观察或试验,使其结果 XX 当理论或假设成立时有确定的分布 FF ,则关于该假设是否正确的问题可以转化为检验假设(a),其中 X1,,XnX_1,\ldots,X_nXX 的独立观察值。 检验方法可用Pearson χ2\chi^2 检验法,也可以用种种不同的方法定义 X1,,XnX_1,\ldots,X_nFF 之间的偏离,从而引出种种不同的检验法。 这些都可以称为拟合优度检验。由K.Pearson开创的这个方向上的工作,是假设检验的一个重要组成部分。

3.Fisher的思想,显著性检验

1919年,Fisher进入Rothamsted农业试验站,从事统计学和遗传学方面的研究工作。他在其著作《The Design of Experiments》中提出了显著性检验。 Fisher通过提出著名的女士饮茶问题说明了他的观点:

奶茶由牛奶与茶按一定比例混合而成,一位女士声称,她可以分辨出一杯奶茶是先倒茶后加奶(TM)还是先加奶后倒茶(MT)。设计如下试验,来检验她的说法是否可信。 准备8杯饮料,TM和MT各半,把它们随机排成一列,让该女士一次品尝,并告诉她TM和MT各有4杯,然后请她指出哪4杯是TM,设她全说对了。

Fisher的推理为:引进一个假设 HH :该女士并无鉴别能力。 其意义为:当 HH 正确时,不论该女士如何做,她事实上只能从8杯奶茶中随机选择4杯作为TM,共有 C84=70C_8^4=70 种,其中只有1种是全部挑对。 因此,若该女士果真全部挑对,则我们必须承认,下述两种情况必有其一发生:

第二种的情况显然比较稀奇,因而有相当的理由承认第一种情况的可能性。即该女士4杯全对这个结果是一个不利于假设 HH 的显著的证据。据此,我们否定 HH 。这样一种推断过程就称为显著性检验。

根据这个例子,不难把Fisher显著性检验的思想归纳为几点:

4.Neyman-Pearson理论

在继承了K.Pearson和Fisher的思想之后,Neyman和E.S.Pearson提出了NP理论,建立了完整的假设检验理论。

(1)建立假设

设有来自某个参数分布族 {F(x,θ)θΘ0}\{F(x,\theta)|\theta\in\Theta_0\} 的样本 x1,,xnx_1,\ldots,x_n ,其中 Θ\Theta 为参数空间, 设 Θ0Θ\Theta_0\subset\Theta ,且 Θ0\Theta_0\neq\varnothing ,则命题 H0:θΘ0H_0:\theta\in\Theta_0 称为原假设或零假设(null hypothesis),若有另一个 Θ1(Θ1Θ,Θ0Θ1=)\Theta_1(\Theta_1\subset\Theta,\Theta_0\cap\Theta_1=\varnothing) ,则命题 H1:θΘ1H_1:\theta\in\Theta_1 称为 H0H_0 的备择假设(alternative hypothesis)。 于是就有了一对假设:

H0:θΘ0vsH1:θΘ1H_0:\theta\in\Theta_0 vs H_1:\theta\in\Theta_1

如果 Θ0\Theta_0 中只含一点,则可以将原假设写成 H0:θ=θ0H_0:\theta=\theta_0 。此时备择假设通常有三种可能:

H1:θθ0,H1:θ<θ0,H1:θ>θ0H_1':\theta\neq\theta_0,H_1'':\theta<\theta_0,H_1''':\theta>\theta_0

我们称 H0vsH1:θθ0H_0 vs H_1':\theta\neq\theta_0 为双侧假设或双边假设, H0vsH1:θ<θ0H_0 vs H_1'':\theta<\theta_0 以及 H0vsH1:θ>θ0H_0 vs H_1''':\theta>\theta_0 为单侧假设或单边假设。

(2)选择检验统计量,给出拒绝形式

由样本对原假设进行检验总是通过一个统计量完成的,该统计量称为检验统计量。

当有了样本之后,按该法则就可以决定是接受 H0H_0 还是拒绝 H0H_0 ,即检验就等价于把样本空间划分为两个互不相交的部分 WWW\overline{W} , 当样本属于 WW 时,拒绝 H0H_0 ;否则接受 H0H_0 。于是我们把 WW 称为拒绝域,而 W\overline{W} 称为接受域。当拒绝域确定了之后,检验的判断准则也随之确定:

(3)两类错误

由于样本是随机的,故当我们应用某种检验做判断时,我们可能作出正确的判断,也可能做出错误的判断。 因此,我们可能犯下两种错误:

在实际中,我们分别称第一类、第二类错误为弃真错误与取伪错误。我们定义犯第一类、第二类错误概率如下:

(4)势函数(power function)

一个检验的势函数为该检验样本观测值 XX 落入拒绝域 WW 内的概率,记为:

g(θ)=P0(XW,θΘ=Θ0Θ1.g(\theta)=P_0(X\in{W}, \theta\in\Theta=\Theta_0\cup\Theta_1.

所以当 θΘ0\theta\in\Theta_0 时, g(θ)=α=α(θ)g(\theta)=\alpha=\alpha(\theta) ,当 θΘ1\theta\in\Theta_1g(θ)=1β(θ)g(\theta)=1-\beta(\theta) 。因此可见,犯两类错误的概率都是参数 θ\theta 的函数, 并可由势函数得到,即:

g(θ)={α(θ),θΘ0,1β(θ),θΘ1,g(\theta) = \left\{ \begin{array}{ll}\alpha(\theta),&\theta\in\Theta_0,\\1-\beta(\theta),&\theta\in\Theta_1,\end{array}\right.

我们可以设一个检验以了解第一类、第二类错误概率的关系。设拒绝域 W={xc},xN(θ,σ2)W=\{\overline{x}\leq{c}\},\overline{x}\sim{N}(\theta,\sigma^2) ,则可以算出该检验的势函数

g(θ)=Pθ(xθσcθσ)=Φ(cθσ)g(\theta)=P_\theta(\frac{\overline{x}-\theta}\sigma\leq\frac{c-\theta}\sigma)=\Phi(\frac{c-\theta}\sigma)

利用这个势函数容易写出其犯两类错误的概率分别为:

α(θ)=Φ(cθσ),θΘ0,\alpha(\theta)=\Phi(\frac{c-\theta}\sigma),\theta\in\Theta_0,β(θ)=1Φ(cθσ),θΘ1.\beta(\theta)=1-\Phi(\frac{c-\theta}\sigma),\theta\in\Theta_1.

由上式可以看出犯两类错误的概率 α,β\alpha,\beta 间的关系:

这说明在样本量一定的条件下不可能找到一个使 α,β\alpha,\beta 都小的检验。 因为无法同时控制检验犯第一类、第二类错误的概率,所以只能采取这种方案,费歇尔提出的显著性检验就仅限制犯第一类错误的概率。

对检验问题 H0:θΘ0vsH1:θΘ1H_0:\theta\in\Theta_0 vs H_1:\theta\in\Theta_1 ,如果一个检验满足对任意的 θΘ0\theta\in\Theta_0 ,都有

g(θ)α,g(\theta)\leq\alpha,

则称该检验是显著性水平为 α\alpha 的显著性检验,简称水平为 α\alpha 的检验。

二、第二类错误计算

1.单个正态总体均值检验

x1,,xnx_1,\ldots,x_n 是来自 N(μ,σ2)N(\mu,\sigma^2) 的样本,总共有如下三种关于 μ\mu 的检验问题:

其中 μ0\mu_0 是已知常数。由于正态总体含两个参数,总体方差 σ2\sigma^2 已知与否对检验有影响。 由NP理论,我们如果想要计算假设检验犯第二类错误的概率,我们需要先计算出检验的拒绝域与接受域。根据接受域以及均值计算出检验犯第二类错误的概率。 下面我们分两种情况讨论在方差是否已知的情况下,在已知显著性条件 α\alpha 下,检验犯第二类错误的概率,以及在限定犯第二类错误概率的情况下保持显著性条件下所需要的样本量。

(1) σ\sigma 已知时的$Z$检验

根据NP理论,在已知 σ\sigma 时,由于 μ\mu 的点估计是 x\overline{x} ,且 xN(μ,σ2/n)\overline{x}\sim{N}(\mu,\sigma^2/n) ,所以我们选用检验统计量:

Z=xμ0σ/nZ=\frac{\overline{x}-\mu_0}{\sigma/\sqrt{n}}

当样本均值 x\overline{x} 不超过设定均值 μ0\mu_0 时,应倾向于接受原假设;当样本均值 x\overline{x} 超过 μ0\mu_0 时,应倾向于拒绝原假设。 可是,在有随机性存在的场合,如果 x\overline{x}μ0\mu_0 大一点就拒绝原假设似乎不当, 只有当 x\overline{x}μ0\mu_0 差距大到一定程度时拒绝原假设才是恰当的,这就存在一个临界值 cc 。设要求显著水平为 α\alpha

1) H1:μμ0H_1:\mu\neq\mu_0

拒绝域与接受域分别为:

W1={(x1,,xn):xμ0c},W_1=\{(x_1,\ldots,x_n):|\overline{x}-\mu_0|\geq{c}\},W1={(x1,,xn):xμ0c}.\overline{W_1}=\{(x_1,\ldots,x_n):|\overline{x}-\mu_0|\leq{c}\}.

cc 满足:

Pμ0(Zcσ/n)=αP_{\mu_0}(|Z|\geq\frac{c}{\sigma/\sqrt{n}})=\alpha

因为在原假设成立时, ZN(0,1)Z\sim{N}(0,1) ,所以

c=σZ1α/2nc=\frac{\sigma Z_{1-\alpha/2}}{\sqrt{n}}

由此根据检验的势函数,我们可以得出犯第二类错误的概率为:

β(μ1)=Φ(c+μ0μ1σ/n)Φ(c+μ0μ1σ/n)\beta(\mu_1)=\Phi(\frac{c+\mu_0-\mu_1}{\sigma/\sqrt{n}}) - \Phi(\frac{-c+\mu_0-\mu_1}{\sigma/\sqrt{n}})

由此,当给定样本均值,显著水平以及犯第二类错误的概率后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

2) H1:μ>μ0H_1:\mu>\mu_0

拒绝域与接受域分别为:

W1={(x1,,xn):xc},W_1=\{(x_1,\ldots,x_n):\overline{x}\geq{c}\},W1={(x1,,xn):xc}.\overline{W_1}=\{(x_1,\ldots,x_n):\overline{x}\leq{c}\}.

cc 满足:

Pμ0(Zcμ0σ/n)=αP_{\mu_0}(Z\geq\frac{c-\mu_0}{\sigma/\sqrt{n}})=\alpha

因为在原假设成立时, ZN(0,1)Z\sim{N}(0,1) ,所以

c=σZ1αn+μ0c=\frac{\sigma{Z_{1-\alpha}}}{\sqrt{n}}+\mu_0

由此根据检验的势函数,我们可以得出犯第二类错误的概率为:

β(μ1)=Φ(cμ1σ/n)\beta(\mu_1)=\Phi(\frac{c-\mu_1}{\sigma/\sqrt{n}})

由此,当给定样本均值,显著水平以及犯第二类错误的概率 β\beta 后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

nσ2(ZβZ1α)2(μ0μ1)2n\geq\frac{\sigma^2(Z_\beta-Z_{1-\alpha})^2}{(\mu_0-\mu_1)^2}

3) H1:μ<μ0H_1:\mu<\mu_0

拒绝域与接受域分别为:

W1={(x1,,xn):xc},W_1=\{(x_1,\ldots,x_n):\overline{x}\leq{c}\},W1={(x1,,xn):xc}.\overline{W_1}=\{(x_1,\ldots,x_n):\overline{x}\geq{c}\}.

cc 满足:

Pμ0(Zcμ0σ/n)=αP_{\mu_0}(Z\leq\frac{c-\mu_0}{\sigma/\sqrt{n}})=\alpha

因为在原假设成立时, ZN(0,1)Z\sim{N}(0,1) ,所以

c=σZαn+μ0c=\frac{\sigma{Z_\alpha}}{\sqrt{n}}+\mu_0

由此根据检验的势函数,我们可以得出犯第二类错误的概率为:

β(μ1)=1Φ(cμ1σ/n)\beta(\mu_1)=1-\Phi(\frac{c-\mu_1}{\sigma/\sqrt{n}})

由此,当给定样本均值,显著水平以及犯第二类错误的概率 β\beta 后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

nσ2(Zβ+Zα)2(μ1μ0)2n\geq\frac{\sigma^2(Z_\beta+Z_\alpha)^2}{(\mu_1-\mu_0)^2}

(2) σ\sigma 未知时的 tt 检验

由于 σ\sigma 未知,我们选择使用 σ\sigma 的无偏估计样本的标准差 ss 来替换 σ\sigma ,这就有了 tt 检验统计量:

t=n(xμ0)s.t=\frac{\sqrt{n}(\overline{x}-\mu_0)}s.

设要求显著水平为 α\alpha ,临界值为 cc

1) H1:μμ0H_1:\mu\neq\mu_0

μ=μ0\mu=\mu_0tt(n1)t\sim t(n-1) ,从而检验 H1:μμ0H_1:\mu\neq\mu_0 的拒绝域为:

W1={(x1,,xn):xμ0c},W_1=\{(x_1,\ldots,x_n):|\overline{x}-\mu_0|\geq{c}\},W1={(x1,,xn):xμ0c}.\overline{W_1}=\{(x_1,\ldots,x_n):|\overline{x}-\mu_0|\leq{c}\}.

cc 满足:

Pμ0(tcs/n)=αP_{\mu_0}(|t|\geq\frac{c}{s/\sqrt{n}})=\alpha

因为在原假设成立时, tt(n1)t\sim t(n-1) ,所以

c=st1α/2(n1)nc=\frac{s t_{1-\alpha/2}(n-1)}{\sqrt{n}}

由此根据检验的势函数,我们可以得出犯第二类错误的概率为:

β(μ1)=tn1(c+μ0μ1s/n)tn1(c+μ0μ1s/n).\beta(\mu_1)=t_{n-1}(\frac{c+\mu_0-\mu_1}{s/\sqrt{n}})-t_{n-1}(\frac{-c+\mu_0-\mu_1}{s/\sqrt{n}}).

由此,当给定样本均值,显著水平以及犯第二类错误的概率后,我们可以由上式求得犯第二类错概率的最小样本量 nn

2) H1:μ>μ0H_1:\mu>\mu_0

μ=μ0\mu=\mu_0tt(n1)t\sim t(n-1) ,从而检验 H1:μ>μ0H_1:\mu>\mu_0 的拒绝域为:

W={tt1α(n1)}W=\{t\geq t_{1-\alpha}(n-1)\}$$\overline{W}=\{t 所以 c=st1α(n1)n+μ0c=\frac{s t_{1-\alpha}(n-1)}{\sqrt{n}}+\mu_0 再根据检验的势函数即可求得犯第二类错误的概率为:

β(μ1)=tn1(cμ1s/n)\beta(\mu_1)=t_{n-1}(\frac{c-\mu_1}{s/\sqrt{n}})

由此,当给定样本均值,显著水平以及犯第二类错误的概率后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

3) H1:μ<μ0H_1:\mu<\mu_0

μ=μ0\mu=\mu_0tt(n1)t\sim t(n-1) ,从而检验 H1:μ<μ0H_1:\mu<\mu_0 的拒绝域为:

W={ttα(n1)}W=\{t\leq t_\alpha(n-1)\}W={t>tα(n1)}\overline{W}=\{t>t_\alpha(n-1)\}

所以 c=stα(n1)n+μ0c=\frac{s t_\alpha(n-1)}{\sqrt{n}}+\mu_0 再根据检验的势函数即可求得犯第二类错误的概率为:

β(μ1)=1tn1(cμ1s/n)\beta(\mu_1)=1-t_{n-1}(\frac{c-\mu_1}{s/\sqrt{n}})

由此,当给定样本均值,显著水平以及犯第二类错误的概率后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

2.两个正态总体均值差的检验

x1,,xmx_1,\ldots,x_m 是来自正态总体 N(μ1,σ12)N(\mu_1,\sigma_1^2) 的样本, y1,,yny_1,\ldots,y_n 是来自另一个正态总体 N(μ2,σ22)N(\mu_2,\sigma_2^2) 的样本,两个样本相互独立。 设原假设为 H0:μ1=μ2H_0:\mu_1=\mu_2 下面对常见的两种情况进行讨论。

(1) σ1,σ2\sigma_1,\sigma_2 已知时的两样本 ZZ 检验

此时 μ1μ2\mu_1-\mu_2 的点估计 xy\overline{x}-\overline{y} 的分布完全已知,

xyN(μ1μ2,σ12m+σ22n).\overline{x}-\overline{y}\sim{N}(\mu_1-\mu_2,\frac{\sigma_1^2}{m}+\frac{\sigma_2^2}{n}).

由此可采用 ZZ 检验的方法,构造检验统计量为:

Z=(xy)/σ12m+σ22n.Z=(\overline{x}-\overline{y})/\sqrt{\frac{\sigma_1^2}m+\frac{\sigma_2^2}n}.

在原假设成立时,有 ZN(0,1)Z\sim{N}(0,1) 。检验的拒绝域则取决于备择假设的具体内容。 下面根据备择假设的不同分情况讨论。

1) H1:μ1μ2H_1:\mu_1\neq\mu_2

对于该检验,检验的拒绝域与接受域分别为:

W={xyc}W=\{|\overline{x}-\overline{y}|\geq c\}W={xyc}\overline{W}=\{|\overline{x}-\overline{y}|\leq c\}

cc 满足:

P0(Zcσ12/m+σ22/n)=αP_0(|Z|\geq\frac{c}{\sqrt{\sigma_1^2/m+\sigma_2^2/n}})=\alpha

所以 c=Z1α/2σ12m+σ22nc=Z_{1-\alpha/2}\sqrt{\frac{\sigma_1^2}m+\frac{\sigma_2^2}n} 再根据势函数即可求得犯第二类错误的概率为:

β(μ1μ20)=Φ(cμ1+μ2σ12/m+σ22/n)Φ(cμ1+μ2σ12/m+σ22/n)\beta(\mu_1-\mu_2\neq 0)=\Phi(\frac{c-\mu_1+\mu_2}{\sqrt{\sigma_1^2/m+\sigma_2^2/n}})-\Phi(\frac{-c-\mu_1+\mu_2}{\sqrt{\sigma_1^2/m+\sigma_2^2/n}})

2) H1:μ1>μ2H_1:\mu_1>\mu_2

对于该检验,拒绝域与接受域分别为:

W={xyc},W=\{\overline{x}-\overline{y}\geq{c}\},W={xyc}.\overline{W}=\{\overline{x}-\overline{y}\leq{c}\}.

cc 满足:

P0(Zcσ12/m+σ22/n)=αP_0(Z\geq\frac{c}{\sqrt{\sigma_1^2/m+\sigma_2^2/n}})=\alpha

所以 c=Z1ασ12/m+σ22/nc=Z_{1-\alpha}\sqrt{\sigma_1^2/m+\sigma_2^2/n} 由此根据检验的势函数,我们可以得出犯第二类错误的概率为:

β(μ1>μ2)=Φ(cμ1+μ2σ12/m+σ22/n)\beta(\mu_1>\mu_2)=\Phi(\frac{c-\mu_1+\mu_2}{\sqrt{\sigma_1^2/m+\sigma_2^2/n}})

3) H1:μ1<μ2H_1:\mu_1<\mu_2

对于该检验,拒绝域与接受域分别为:

W={xyc},W=\{\overline{x}-\overline{y}\leq{c}\},W={xyc}.\overline{W}=\{\overline{x}-\overline{y}\geq{c}\}.

cc 满足:

P0(Zcσ12/m+σ22/n)=αP_0(Z\leq\frac{c}{\sqrt{\sigma_1^2/m+\sigma_2^2/n}})=\alpha

所以 c=Zασ12/m+σ22/nc=Z_\alpha\sqrt{\sigma_1^2/m+\sigma_2^2/n} 由此根据检验的势函数,我们可以得出犯第二类错误的概率为:

β(μ1<μ2)=1Φ(cμ1+μ2σ12/m+σ22/n)\beta(\mu_1<\mu_2)=1-\Phi(\frac{c-\mu_1+\mu_2}{\sqrt{\sigma_1^2/m+\sigma_2^2/n}})

(2) σ1=σ2=σ\sigma_1=\sigma_2=\sigma 未知时的两样本 tt 检验

σ12=σ22=σ2\sigma_1^2=\sigma_2^2=\sigma^2 但未知时,首先 xyN(μ1μ2,(1m+1n)σ2)\overline{x}-\overline{y}\sim N(\mu_1-\mu_2,(\frac{1}m+\frac{1}n)\sigma^2) ,其次,由于

1σ2i=1m(xix)2χ2(m1),1σ2i=1n(yiy)2χ2(n1)).\frac{1}{\sigma^2}\sum_{i=1}^m(x_i-\overline{x})^2\sim \chi^2(m-1),\frac{1}{\sigma^2}\sum_{i=1}^n(y_i-\overline{y})^2\sim \chi^2(n-1)).

1σ2((xix)2+(yiy)2)χ2(m+n2)\frac{1}{\sigma^2}(\sum{(x_i-\overline{x})^2}+\sum{(y_i-\overline{y})^2})\sim \chi^2(m+n-2) ,记

sw2=1m+n2[i=1m(xix)2+i=1n(yiy)2].s_w^2=\frac{1}{m+n-2}[\sum_{i=1}^m(x_i-\overline{x})^2+\sum_{i=1}^n(y_i-\overline{y})^2].

于是有:

t=(xy)(μ1μ2)sw1m+1nt(m+n2)t=\frac{(\overline{x}-\overline{y})-(\mu_1-\mu_2)}{s_w\sqrt{\frac{1}m+\frac{1}n}}\sim t(m+n-2)

构造 tt 统计量:

t=xysw1m+1n.t=\frac{\overline{x}-\overline{y}}{s_w\sqrt{\frac{1}m+\frac{1}n}}.

设显著性要求为 α\alpha ,临界值为 cc 。下面根据备择假设的不同,分类讨论。

1) H1:μ1μ2H_1:\mu_1\neq\mu_2

对于该检验,拒绝域与接受域分别为:

W={μ1μ2c}W=\{|\mu_1-\mu_2|\geq{c}\}W={μ1μ2c}\overline{W}=\{|\mu_1-\mu_2|\leq{c}\}

cc 满足

P0(tcsw1m+1n)=αP_0(|t|\geq\frac{c}{s_w\sqrt{\frac{1}m+\frac{1}n}})=\alpha

所以 c=sw1m+1nt1α/2(m+n2)c=s_w\sqrt{\frac{1}m+\frac{1}n}t_{1-\alpha/2(m+n-2)} 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(μ1μ2)=tm+n2(cμ1+μ2sw1m+1n)tm+n2(cμ1+μ2sw1m+1n)\beta(\mu_1\neq\mu_2)=t_{m+n-2}(\frac{c-\mu_1+\mu_2}{s_w\sqrt{\frac{1}m+\frac{1}n}})-t_{m+n-2}(\frac{-c-\mu_1+\mu_2}{s_w\sqrt{\frac{1}m+\frac{1}n}})

2) H1:μ1>μ2H_1:\mu_1>\mu_2

对于该检验,拒绝域与接受域分别为:

W={tt1α(m+n2)},W=\{t\geq t_{1-\alpha}(m+n-2)\},W={tt1α(m+n2)}.\overline{W}=\{t\leq t_{1-\alpha}(m+n-2)\}.

所以 c=sw1m+1nt1α(m+n2)c=s_w\sqrt{\frac{1}m+\frac{1}n}t_{1-\alpha}(m+n-2) 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(μ1>μ2)=tm+n2(cμ1+μ2sw1m+1n)\beta(\mu_1>\mu_2)=t_{m+n-2}(\frac{c-\mu_1+\mu_2}{s_w\sqrt{\frac{1}m+\frac{1}n}})

3) H1:μ1<μ2H_1:\mu_1<\mu_2

对于该检验,拒绝域与接受域分别为:

W={tt1α(m+n2)},W=\{t\leq t_{1-\alpha}(m+n-2)\},W={tt1α(m+n2)}.\overline{W}=\{t\geq t_{1-\alpha}(m+n-2)\}.

所以 c=sw1m+1ntα(m+n2)c=s_w\sqrt{\frac{1}m+\frac{1}n}t_\alpha(m+n-2) 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(μ1<μ2)=1tm+n2(cμ1+μ2sw1m+1n)\beta(\mu_1<\mu_2)=1-t_{m+n-2}(\frac{c-\mu_1+\mu_2}{s_w\sqrt{\frac{1}m+\frac{1}n}})

3.正态总体方差的检验

(1)单个正态总体方差的 χ2\chi^2 检验

x1,,xnx_1,\ldots,x_n 是来自 N(μ,σ2)N(\mu,\sigma^2) 的样本,在这里我们采用 χ2\chi^2 检验统计量:

χ2=(n1)s2/σ02\chi^2=(n-1)s^2/\sigma_0^2

当原假设 H0:σ2=σ02H_0:\sigma^2=\sigma_0^2 时,有 χ2χ2(n1)\chi^2\sim\chi^2(n-1) ,于是若取显性水平为 α\alpha ,临界值为 c,c1,c2c,c_1,c_2 。下面根据备择假设的不同,分类讨论。

1) H1:σ2σ02H_1:\sigma^2\neq\sigma_0^2

对于该检验,拒绝域与接受域分别为:

W={χ2χα/22(n1)orχ2χ1α/22(n1)}W=\{\chi^2\leq\chi_{\alpha/2}^2(n-1)or\chi^2\geq\chi_{1-\alpha/2}^2(n-1)\}W={χα/22(n1)χ2χ1α/22(n1)}\overline{W}=\{\chi_{\alpha/2}^2(n-1)\leq\chi^2\leq\chi_{1-\alpha/2}^2(n-1)\}

则有

P0(χ2nc12σ02)+P0(χ2nc22σ02)=αP_0(\chi^2\leq\frac{nc_1^2}{\sigma_0^2})+P_0(\chi^2\geq\frac{nc_2^2}{\sigma_0^2})=\alpha

所以 c12=σ02χα/22(n1)n,c22=σ02χ1α/22(n1)nc_1^2=\frac{\sigma_0^2\chi_{\alpha/2}^2(n-1)}n,c_2^2=\frac{\sigma_0^2\chi_{1-\alpha/2}^2(n-1)}n 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(σ2σ02)=χn12(nc22σ12)χn12(nc12σ12)\beta(\sigma^2\neq\sigma_0^2)=\chi_{n-1}^2(\frac{nc_2^2}{\sigma_1^2})-\chi_{n-1}^2(\frac{nc_1^2}{\sigma_1^2})

由此,当给定样本标准差,显著水平以及犯第二类错误的概率后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

2) H1:σ2σ02H_1:\sigma^2\geq\sigma_0^2

对于该检验,拒绝域与接受域分别为:

W={χ2χ1α2(n1)}W=\{\chi^2\geq\chi_{1-\alpha}^2(n-1)\}W={χ2χ1α2(n1)}\overline{W}=\{\chi^2\leq\chi_{1-\alpha}^2(n-1)\}

则有

P0((n1)s2σ02nc2σ02)=αP_0(\frac{(n-1)s^2}{\sigma_0^2}\geq\frac{nc^2}{\sigma_0^2})=\alpha

所以 c2=σ02χ1α2(n1)nc^2=\frac{\sigma_0^2\chi_{1-\alpha}^2(n-1)}n 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(σ2σ02)=χn12(nc2σ12)\beta(\sigma^2\geq\sigma_0^2)=\chi_{n-1}^2(\frac{nc^2}{\sigma_1^2})

由此,当给定样本标准差,显著水平以及犯第二类错误的概率后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

3) H1:σ2σ02H_1:\sigma^2\leq\sigma_0^2

对于该检验,拒绝域与接受域分别为:

W={χ2χα2(n1)}W=\{\chi^2\leq\chi_\alpha^2(n-1)\}W={χ2χα2(n1)}\overline{W}=\{\chi^2\geq\chi_\alpha^2(n-1)\}

则有

P0((n1)s2σ02nc2σ02)=αP_0(\frac{(n-1)s^2}{\sigma_0^2}\leq\frac{nc^2}{\sigma_0^2})=\alpha

所以 c2=σ02χα2(n1)nc^2=\frac{\sigma_0^2\chi_\alpha^2(n-1)}n 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(σ2σ02)=1χn12(nc2σ12)\beta(\sigma^2\leq\sigma_0^2)=1-\chi_{n-1}^2(\frac{nc^2}{\sigma_1^2})

由此,当给定样本标准差,显著水平以及犯第二类错误的概率后,我们可以根据上式求得满足犯第二类错误概率的最小样本量 nn 为多少。

(2)两个正态总体方差比的 FF 检验

x1,,xmx_1,\ldots,x_m 是来自 N(μ1,σ12)N(\mu_1,\sigma_1^2) 的样本, y1,,yny_1,\ldots,y_n 是来自 N(μ2,σ22)N(\mu_2,\sigma_2^2) 的样本。 通常 μ1,μ2\mu_1,\mu_2 均未知,记 sx2,sy2s_x^2,s_y^2 分别是由 x1,,xmx_1,\ldots,x_m 算得的 σ12\sigma_1^2 的无偏估计和由 y1,,yny_1,\ldots,y_n 算得的 σ22\sigma_2^2 的无偏估计(两个都是样本方差),则可建立如下检验统计量:

F=sx2sy2.F=\frac{s_x^2}{s_y^2}.

当原假设 H0:σ12=σ22H_0:\sigma_1^2=\sigma_2^2 成立时,有 FF(m1,n1)F\sim F(m-1,n-1) ,于是若取显性水平为 α\alpha ,临界值为 c1,c2,c1,c2c_1,c_2,c_1',c_2' 。下面根据备择假设的不同,分类讨论。

1) H1:σ12σ22H_1:\sigma_1^2\neq\sigma_2^2

对于该检验,拒绝域与接受域分别为:

W={FFα/2(m1,n1)orFF1α/2(m1,n1)}W=\{F\leq F_{\alpha/2}(m-1,n-1)orF\geq F_{1-\alpha/2}(m-1,n-1)\}W={Fα/2(m1,n1)FF1α/2(m1,n1)}\overline{W}=\{F_{\alpha/2}(m-1,n-1)\leq F\leq F_{1-\alpha/2}(m-1,n-1)\}

则有

P0(F(n1)mc12(m1)nc22)+P0(F(n1)mc12(m1)nc22)=αP_0(F\leq\frac{(n-1)mc_1^2}{(m-1)nc_2^2})+P_0(F\geq\frac{(n-1)m{c_1'}^2}{(m-1)n{c_2'}^2})=\alpha

所以 c12c22=(m1)nFα/2(m1,n1)(n1)m,c12c22=(m1)nF1α/2(m1,n1)(n1)m\frac{c_1^2}{c_2^2}=\frac{(m-1)n F_{\alpha/2}(m-1,n-1)}{(n-1)m},\frac{{c_1'}^2}{{c_2'}^2}=\frac{(m-1)n F_{1-\alpha/2}(m-1,n-1)}{(n-1)m} 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(σ12σ22)=Fm1,n1((n1)mσ22c12(m1)nσ12c22)Fm1,n1((n1)mσ22c12(m1)nσ12c22)\beta(\sigma_1^2\neq\sigma_2^2)=F_{m-1,n-1}(\frac{(n-1)m\sigma_2^2{c_1'}^2}{(m-1)n\sigma_1^2{c_2'}^2})-F_{m-1,n-1}(\frac{(n-1)m\sigma_2^2c_1^2}{(m-1)n\sigma_1^2c_2^2})

2) H1:σ12σ22H_1:\sigma_1^2\geq\sigma_2^2

对于该检验,拒绝域与接受域分别为:

W={FF1α(m1,n1)}W=\{F\geq F_{1-\alpha}(m-1,n-1)\}W={FF1α(m1,n1)}\overline{W}=\{F\leq F_{1-\alpha}(m-1,n-1)\}

则有

P0(Fmm1c12nn1c22)=αP_0(F\geq\frac{\frac{m}{m-1}c_1^2}{\frac{n}{n-1}c_2^2})=\alpha

所以 c12c22=(m1)nF1α(m1,n1)(n1)m\frac{c_1^2}{c_2^2}=\frac{(m-1)n F_{1-\alpha}(m-1,n-1)}{(n-1)m} 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(σ12σ22)=Fm1,n1((n1)mσ22c12(m1)nσ12c22)\beta(\sigma_1^2\geq\sigma_2^2)=F_{m-1,n-1}(\frac{(n-1)m\sigma_2^2c_1^2}{(m-1)n\sigma_1^2c_2^2})

3) H1:σ12σ22H_1:\sigma_1^2\leq\sigma_2^2

对于该检验,拒绝域与接受域分别为:

W={FFα(m1,n1)}W=\{F\leq F_\alpha(m-1,n-1)\}W={FFα(m1,n1)}\overline{W}=\{F\geq F_\alpha(m-1,n-1)\}

则有

P0(Fmm1c12nn1c22)=αP_0(F\leq\frac{\frac{m}{m-1}c_1^2}{\frac{n}{n-1}c_2^2})=\alpha

所以 c12c22=(m1)nFα(m1,n1)(n1)m\frac{c_1^2}{c_2^2}=\frac{(m-1)n F_\alpha(m-1,n-1)}{(n-1)m} 由此根据检验的势函数,我们可以得出犯第二类错误的概率:

β(σ12σ22)=1Fm1,n1((n1)mσ22c12(m1)nσ12c22)\beta(\sigma_1^2\leq\sigma_2^2)=1-F_{m-1,n-1}(\frac{(n-1)m\sigma_2^2c_1^2}{(m-1)n\sigma_1^2c_2^2})

R语言实现及例题验证

R语言实现

根据上述公式,我们发现有些的第二类错误并不容易算出来。 好在统计软件的快速发展,可以使我们快速的得到各个分布的值以及犯第二类错误的概率达到某一值所需要的最小n。 本文将通过R语言来实现以上所有的公式,同时为了便于计算,我利用公式编写了一系列函数用来计算犯第二类错误的概率以及指定$\beta$的最小的样本量。其原代码附在文后。

之后我将利用一些建设检验的相关例题,结合编写的R语言的函数来计算相关数值。

例题验证

1.单个正态总体均值检验( σ\sigma 已知)

1)设某产品的指标服从正态分布,它的标准差 σ\sigma 已知为150,今抽了一个容量为26的样本,计算得平均值为1637。问在5%的显著水平下,认为这批产品的指标的期望值 μ\mu 为1600,其犯第二类错误的概率为多少,控制犯第二类错误在0.1之内所需要的样本量为?

由题,原假设 H0:μ=1600H_0:\mu=1600 ,备择假设为 H1:μ1600H_1:\mu\neq 1600 。使用R语言计算得到结果:

result = secerror(rnorm(26,1637,150),mu=1600,mu1 = 1637,
                  alpha = 0.05,beta = 0.1,sigma = 150,alternative = "two-side")
## 犯第二类错误的概率为:P= 0.7580785 
##  如要控制犯第二类错误的概率在 0.1 之内 
##  最小样本量为: 179

2)某纺织厂在正常的运转条件下,平均每台布机每小时经纱断头数为O.973根,各台布机断头数的标准差为O.16根,该厂进行工艺改进,减少经纱上浆率,在200台布机上进行试验,结果平均每台每小时经纱断头数为O.994根,标准差为0.16根。新工艺上浆率推广( α=0.05\alpha = 0.05 )犯第二类错误的概率以及控制犯第二类错误在0.1之内所需要的样本量?

由题,原假设 H0:μ>0.973H_0:\mu>0.973 ,备择假设为 H1:μ0.973H_1:\mu\leq 0.973 。使用R语言计算得到结果:

result = secerror(x=rnorm(200,0.973,0.162),mu=0.973,mu1 = 0.994,
                  alpha = 0.05,beta = 0.1,sigma = 0.16,alternative = "less")
## 犯第二类错误的概率为:P= 0.9997682 
##  如要控制犯第二类错误的概率在 0.1 之内 
##  最小样本量为: 498

2.单个正态总体均值检验( σ\sigma 未知)

3)已知某种元件的寿命服从正态分布,要求该元件的平均寿命不低于1000h,现从这批元件中随机抽取25 知,测得平均寿命 x=980h\overline{x}=980h ,标准差 σ=65h\sigma=65h ,试在水平 α=0.05\alpha=0.05 下, 推断符合条件后犯第二类错误的概率为多少?控制第二类错误的概率为0.1的最小样本量为多少?

result = secerror(x=rnorm(25,980,65),mu=1000,mu1 = 980,
                  alpha = 0.05,beta = 0.1,sigma1 = 65,alternative = "less")
## 犯第二类错误的概率为:P= 0.5677237 
##  如要控制犯第二类错误的概率在 0.1 之内 
##  最小样本量为: 95

3.两个正态总体均值差的检验

4)下面给出了两个文学家马克·吐温(Mark Twain)的8篇小品文以及斯诺·特格拉斯(Snodgrass)的10篇小品文中由3格字母组成的词比例.

马克·吐温: 0.225,0.262,0.217,0.240,0.230,0.229,0.235,0.217
斯诺·特格拉斯:0.209,0.205,0.196,0.210,0.202,0.207,0.224,0.223,0.220,0.201

设两组数据分别来自正态分布,且两总体方差相等,两样本相互独立, 推断两个作家所写的小品文中包含由3格字母组成的词的比例有显著性的差异( α=0.05\alpha=0.05 )犯第二类错误的概率为多少?

result = secerror(x=c(0.225,0.262,0.217,0.240,0.230,0.229,0.235,0.217),
                  y=c(0.209,0.205,0.196,0.210,0.202,0.207,0.224,0.223,0.220,0.201),
                  alpha = 0.05,alternative = "two-side")
## 犯第二类错误的概率为:P= 0.04889964

4.单个样本方差检验

5)某厂生产的一中电池,其寿命长期以来服从方差 σ2=5000\sigma^2=5000 小时的正态分布。 现有一批这种电池,从生产的情况来看,寿命的波动性有所改变,现随机地抽取26只电池,测得寿命的样本方差 σ2=9200\sigma^2=9200 小时。 根据这一数据推断这批电池寿命的波动性较以往有显著性的变化(取 α=0.02\alpha=0.02 ),犯第二类错误的概率为多少?控制第二类错误的概率为0.1的最小样本量为多少?

result = var_secerror(x = rnorm(26,0,sqrt(9200)),alpha = 0.02,beta = 0.1,
                      sigma = sqrt(5000),sigma1 = sqrt(9200),alternative = "two-side")
## 犯第二类错误的概率为:P= 0.4854096 
##  如要控制犯第二类错误的概率在 0.1 之内 
##  最小样本量为: 68

5.双样本方差检验

6)在十块地上同时试种甲、乙两种品种作物,设每种作物的产量服从正态分布,并计算得 x=30.97,y=21.79,sx=26.7,sy=12.1x=30.97,y=21.79,s_x=26.7,s_y=12.1 。 这两种品种的产量在 α=O.01\alpha=O.01 下,推论犯第二类错误的概率为多少?

var_secerror(x=rnorm(10,30.97,26.7),y=rnorm(10,21.79,12.1),alpha = 0.01,
             sigma = 1,sigma1 = 26.7,sigma2 = 12.1,alternative = "two-side")
## 犯第二类错误的概率为:P= 0.666364
## NULL

总结

本文通过着力于分析假设检验的计算原理,利用势函数,求出在各种情况下的检验推论犯第二类错误的概率。 同时针对单个样本,我们结合其给定第二类概率,可以计算出最小样本量 nn 。 我们可以发现,尽管样本情况繁多,但是,只要结合势函数,以及第二类错误的基本定义,我们依然可以方便的算出检验犯第二类错误的概率。 本次研究为了可以得到快速且准确的解,我大量的使用了R语言的编程功能,利用自己编写的函数,快速地计算出各种情况不同样本的检验情况。 该函数覆盖了所有双样本以内正态总体的假设检验的情况,可在日后需要时方便快速使用。

本次研究注重于公式推导以及计算机程序语言实现。借助于日益强大的计算机,我相信,日后统计学会越来越多的依靠统计软件。 同时计算机的运用,也会使得统计计算变成一件简单的事。

函数程序源代码

Pb = function(params,sigmaknown){
  if(sigmaknown){
    c = params[5]*qnorm(1-params[2]/2,0,1)/sqrt(params[1])
    Pb = abs(pnorm((c+params[3]-params[4])*sqrt(params[1])/params[5])-
               pnorm((-c+params[3]-params[4])*sqrt(params[1])/params[5])-params[6])
    return(Pb)
  }else{
    c = params[5]*qt(1-params[2]/2,0,1)/sqrt(params[1],n-1)
    Pb = abs(pt((c+params[3]-params[4])*sqrt(params[1])/params[5],n-1)-
               pt((-c+params[3]-params[4])*sqrt(params[1])/params[5])-params[6],n-1)
    return(Pb)
  }
}
Pb_1 = function(params){
  c = params[5]*qt(params[2],params[1]-1)/sqrt(params[1])+params[3]
  Pb_1 = abs(1-pt((c-params[4])/params[5]*sqrt(params[1]),params[1]-1)-params[6])
  return(Pb_1)
}
Pb_2 = function(params){
  c = params[5]*qt(1-params[2],params[1]-1)/sqrt(params[1])+params[3]
  Pb_2 = abs(pt((c-params[4])/params[5]*sqrt(params[1]),params[1]-1)-params[6])
  return(Pb_2)
}
n_min = function(n,alpha,mu0,mu1,sigma,beta,sigmaknown){
  rex = nlminb(c(n,alpha,mu0,mu1,sigma,beta),Pb,sigmaknown=sigmaknown,
             lower = c(n,alpha,mu0,mu1,sigma,0),
             upper = c(Inf,alpha,mu0,mu1,sigma,beta))
  return(rex$par)
}
n_min_1 = function(n,alpha,mu0,mu1,sigma,beta){
  rex = nlminb(c(n,alpha,mu0,mu1,sigma,beta),Pb_1,
               lower = c(n,alpha,mu0,mu1,sigma,0),
               upper = c(Inf,alpha,mu0,mu1,sigma,beta))
  return(rex$par)
}
n_min_2 = function(n,alpha,mu0,mu1,sigma,beta){
  rex = nlminb(c(n,alpha,mu0,mu1,sigma,beta),Pb_2,
               lower = c(n,alpha,mu0,mu1,sigma,0),
               upper = c(Inf,alpha,mu0,mu1,sigma,beta))
  return(rex$par)
}
secerror = function(x,y=NULL,alpha=0.05,beta=0.05,
                    mu=0,mu1=NULL,mu2=NULL,sigma=NULL,sigma1=NULL,sigma2=NULL,
                    alternative=c("two-side","less","greater"),varequal=TRUE){
  n = length(x)
  ifelse(is.null(mu1),mu1<-mean(x),mu1<-mu1)
  if(is.null(y)){if(is.null(sigma)){              #don't know sigma  
    sigma = sigma1
    if(alternative=="two-side"){                  #H_1:mu<>mu0
      c=sigma*qt(1-alpha/2,n-1)/sqrt(n)
      P = pt((c+mu-mu1)*sqrt(n)/sigma,n-1)-pt((-c+mu-mu1)*sqrt(n)/sigma,n-1)
      n = n_min(n,alpha,mu,mu1,sigma,beta,FALSE)[1]
      n = ceiling(n)
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
    if(alternative=="less"){                      #H_1:mu<mu0  
      c = sigma*qt(alpha,n-1)/sqrt(n)+mu
      P = 1-pt((c-mu1)*sqrt(n)/sigma,n-1)
      n = n_min_1(n,alpha,mu,mu1,sigma,beta)[1]
      n = ceiling(n)
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
    if(alternative=="greater"){                   #H_1:mu>mu0
      c = sigma*qt(1-alpha,n-1)/sqrt(n)+mu
      P = pt((c-mu1)*sqrt(n)/sigma,n-1)
      n = n_min_2(n,alpha,mu,mu1,sigma,beta)[1]
      n = ceiling(n)
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
    return(result)
    }
  }else{                                         #sigma have known
    if(alternative=="two-side"){                 #H_1:mu<>mu0
      c=sigma*qnorm(1-alpha/2)/sqrt(n)
      P = pnorm((c+mu-mu1)*sqrt(n)/sigma)-pnorm((-c+mu-mu1)*sqrt(n)/sigma)
      n = n_min(n,alpha,mu,mu1,sigma,beta,TRUE)[1]
      n = ceiling(n)
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
    if(alternative=="less"){                      #H_1:mu<mu0
      c = sigma*qnorm(alpha)/sqrt(n)+mu
      P = 1-pnorm((c-mu1)*sqrt(n)/sigma)
      n = ceiling(sigma^2*(qnorm(beta)+qnorm(alpha))^2/(mu1-mu)^2)
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
    if(alternative=="greater"){                  #H-1:mu>mu0
      c = sigma*qnorm(1-alpha)/sqrt(n)+mu
      P = pnorm((c-mu1)*sqrt(n)/sigma)
      n = ceiling(sigma^2*(qnorm(beta)-qnorm(1-alpha))^2/(mu-mu1)^2)
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
  }
  }else{                                         #double sample sigma
    m=n;n=length(y)
    ifelse(is.null(mu2),mu2<-mean(y),mu2<-mu2)
    if(is.null(sigma1)|is.null(sigma2)){         #double sample sigma don't know
      sw = sqrt(((m-1)*var(x)+(n-1)*var(y))/(m+n-2))
      if(alternative=="two-side"){               #H_1:mu1<>mu2
      c=sw*sqrt(1/m+1/n)*qt(1-alpha/2,m+n-2)
      P = pt((c+mu2-mu1)/sw/sqrt(1/m+1/n),m+n-2)-
        pt((-c+mu2-mu1)/sw/sqrt(1/m+1/n),m+n-2)
      result = cat("犯第二类错误的概率为:P=",P,"\n")
      return(result)
    }
    if(alternative=="less"){                     #H_1:mu1<mu2
      c = sw*sqrt(1/m+1/n)*qt(alpha,m+n-2)
      P = 1-pt((c-mu1+mu2)/sqrt(1/m+1/n)/sw,m+n-2)
      result = cat("犯第二类错误的概率为:P=",P,"\n")
      return(result)
    }
    if(alternative=="greater"){                 #H_1:mu>mu2
    c = sw*sqrt(1/m+1/n)*qt(1-alpha,m+n-2)
    P = pt((c-mu1+mu2)/sqrt(1/m+1/n)/sw,m+n-2)
    result = cat("犯第二类错误的概率为:P=",P,"\n")
    return(result)
    }
    }else{                                      #double sample sigma have known
      if(alternative=="two-side"){              #H_1:mu1<>mu2
        c=sigma*qnorm(1-alpha/2)/sqrt(sigma1^2/m+sigma^2/n)
        P = pnorm((c-mu1+mu2)/sqrt(sigma1^2/m+sigma^2/n))-
          pnorm((-c-mu1+mu2)/sqrt(sigma1^2/m+sigma^2/n))
        result = cat("犯第二类错误的概率为:P=",P,"\n")
        return(result)
      }
      if(alternative=="less"){                    #H_1:mu1<mu2
        c = qnorm(alpha)/sqrt(sigma1^2/m+sigma^2/n)
        P = 1-pnorm((c-mu1+mu2)/sqrt(sigma1^2/m+sigma^2/n))
        result = cat("犯第二类错误的概率为:P=",P,"\n")
        return(result)
      }
      if(alternative=="greater"){                 #H_1:mu>mu2
        c = qnorm(1-alpha)/sqrt(sigma1^2/m+sigma^2/n)
        P = pnorm((c-mu1+mu2)/sqrt(sigma1^2/m+sigma^2/n))
        result = cat("犯第二类错误的概率为:P=",P,"\n")
        return(result)
      }
    }
    }
}

var_secerror = function(x,y=NULL,alpha=0.05,beta=0.1,
                        sigma=NULL,sigma1=NULL,sigma2=NULL,
                        alternative=c("two-side","less","greater")){
  if(is.null(sigma)){
    return("Please enter sigma")
  }
  ifelse(is.null(sigma1),sigma1<-sd(x),sigma1<-sigma1)
  if(is.null(y)){                               #single sample
    n = length(x)
    if(alternative=="two-side"){
      c_1 = sigma^2*qchisq(alpha/2,n-1)/n;
      c_2 = sigma^2*qchisq(1-alpha/2,n-1)/n
      P = pchisq(n*c_2/sigma1^2,n-1)-pchisq(n*c_1/sigma1^2,n-1)
      Pb = function(params){
        c_1 = params[2]^2*qchisq(params[4]/2,params[1]-1)/params[1];
        c_2 = params[2]^2*qchisq(1-params[4]/2,params[1]-1)/params[1]
        Pb = abs(pchisq(params[1]*c_2/params[3]^2,params[1]-1)-
                   pchisq(params[1]*c_1/params[3]^2,params[1]-1)-params[5])
        return(Pb)
      }
      n_min_sd = function(n,sigma,sigma1,alpha,beta){
        res = nlminb(c(n,sigma,sigma1,alpha,beta),Pb,
                     lower = c(n,sigma,sigma1,alpha,beta),
                     upper = c(Inf,sigma,sigma1,alpha,beta))
        return(res$par[1])
      }
      n = ceiling(n_min_sd(n,sigma,sigma1,alpha,beta))
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
    if(alternative=="less"){
      c = sigma^2*qchisq(alpha,n-1)/n
      P = 1-pchisq(n*c/sigma1^2,n-1)
      Pb = function(params){
        c = params[2]^2*qchisq(1-params[4],params[1]-1)/params[1]
        Pb = abs(pchisq(params[1]*c/params[3]^2,params[1]-1)-params[5])
        return(Pb)
      }
      n_min_sd = function(n,sigma,sigma1,alpha,beta){
        res = nlminb(c(n,sigma,sigma1,alpha,beta),Pb,
                     lower = c(n,sigma,sigma1,alpha,beta),
                     upper = c(Inf,sigma,sigma1,alpha,beta))
        return(res$par[1])
      }
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
    if(alternative=="greater"){
      c = sigma^2*qchisq(1-alpha,n-1)/n
      P = pchisq(n*c/sigma1^2,n-1)
      Pb = function(params){
        c = params[2]^2*qchisq(params[4],params[1]-1)/params[1]
        Pb = abs(1-pchisq(params[1]*c/params[3]^2,params[1]-1)-params[5])
        return(Pb)
      }
      n_min_sd = function(n,sigma,sigma1,alpha,beta){
        res = nlminb(c(n,sigma,sigma1,alpha,beta),Pb,
                     lower = c(n,sigma,sigma1,alpha,beta),
                     upper = c(Inf,sigma,sigma1,alpha,beta))
        return(res$par[1])
      }
      result = cat("犯第二类错误的概率为:P=",P,"\n",
                   "如要控制犯第二类错误的概率在",beta,"之内","\n",
                   "最小样本量为:",n,"\n")
      return(result)
    }
  }else{                                        #double sample
    m=length(x);n=length(y)
    ifelse(is.null(sigma2),sigma2<-sd(y),sigma2<-sigma2)
    if(alternative=="two-side"){
      c_11_c_21 = (m-1)*n*qf(alpha/2,m-1,n-1)/(n-1)/m;
      c_12_c_22 = (m-1)*n*qf(1-alpha/2,m-1,n-1)/(n-1)/m
      P = pf((n-1)*m*sigma2^2/((m-1)*n*sigma1^2)*c_12_c_22,m-1,n-1)-
        pf((n-1)*m*sigma2^2/((m-1)*n*sigma1^2)*c_11_c_21,m-1,n-1)
      result = cat("犯第二类错误的概率为:P=",P,"\n")
      return(result)
    }
    if(alternative=="less"){
      c_1_c_2 = (m-1)*n*qf(alpha,m-1,n-1)/(n-1)/m
      P = 1-pf((n-1)*m*sigma2^2/((m-1)*n*sigma1^2)*c_1_c_2,m-1,n-1)
      result = cat("犯第二类错误的概率为:P=",P,"\n")
      return(result)
    }
    if(alternative=="greater"){
      c_1_c_2 = (m-1)*n*qf(1-alpha,m-1,n-1)/(n-1)/m
      P = pf((n-1)*m*sigma2^2/((m-1)*n*sigma1^2)*c_1_c_2,m-1,n-1)
      result = cat("犯第二类错误的概率为:P=",P,"\n")
      return(result)}}}

参考文献

[1]《数理统计学简史》.陈希孺.湖南教育出版社.2002年

[2]《概率论与数理统计教程》. 茆诗松.程依明.濮晓龙.高等教育出版社.2011年

[3]《数理统计学教程》.陈希孺.中国科学技术大学出版社.2009年

[4]《论假设检验中的两类错误》.蔡越江.数理统计与管理.18卷3期1999年3月