史海芳, 姬永刚
(1. 吉林大学 数学研究所, 长春 130012; 2. 中国民航大学 理学院, 天津 300300)
考虑单向方差分析模型:
yij=μi+ij,j=1,2,…,ni,i=1,2,…,k.
(1)
其中: 响应变量yij表示第i种试验第j个个体的值;μi为第i种试验的平均值; 随机误差ij为相互独立的服从正态分布N(0,σ2)的随机变量.
在很多实际问题中, 通过一些先验信息可知总体均值满足一定的序约束. 如简单半序μ1≤…≤μk, 简单树半序μ1≤μi(i=2,3,…,k), 伞型半序μ1≤…≤μp≥μp+1≥…≥μk等. 文献[1-2]给了在这些先验信息下一些感兴趣参数的估计和检验. 此外, 文献[3-5]分别考虑了构造总体均值置信区间的问题. 进一步, 文献[6]介绍了用U-I(union-intersection)检验法考虑似然比检验统计量的精确分布, 通过计算p值, 解决了简单半序约束下的检验问题. 但上述的频率方法在理论上一般较复杂, 实际操作不太方便. 因此, 近年人们试图从Bayes角度考虑序约束下的估计和检验问题, 如文献[7-8]分别用Bayes方法考虑了如下等值检验问题:
H0:μ1=…=μkv.s.H1:μ1≤…≤μk,μ1<μk;
(2)
的估计, 同时考虑了如下等值检验问题:
(3)
最后用Gibbs抽样和Metropolis-Hastings方法给出了数值模拟. 与非Bayes方法相比, Bayes方法在处理中、 小样本问题时具有很大优势. 特别地, 如果拒绝原假设, 根据Bayes方法还可以知道哪些总体均值和标准差之比是不同的.
[δhρh,θh]=ρhI(δh=0)+(1-ρh)g(δhθh)I(δh>0),h=1,2,…,k-1,
(4)
(5)
综上, Bayes分层模型为:
由Bayes分层模型知, 所有未知参数包括δj,ρj,θj,ξ1,σi(j=1,2,…,k-1;i=1,2,…,k).
由Bayes公式知
因此, 给定y,{σi},ρj,θj,ξ1,δ-j下,δj的满条件后验分布是一个混合分布:
其中
(7)
由Bayes公式知
由Bayes公式知
ρi的满条件后验分布可以记作
当α0=β0=1时, 有
由Bayes公式知
θi的满条件后验分布可以记作
由Bayes公式, 有
显然,σi的满条件后验不是一个常见分布, 所以本文用Metropolis-Hastings方法给出σi的数值模拟.
Gibbs抽样步骤如下:
2) ① 先通过下式计算fj和ej:
(8)
解方程组(8)得
② 计算r.
③ 从(0,1)上的均匀分布抽取一随机数u:
下面通过数值模拟检验本文的方法. 考虑总体个数k=4时的情况, 已知δj=μj+1/σj+1-μj/σj,j=1,2,3. 共模拟4组总体. 为方便, 令每组总体中均值和标准差比的间隔δj都分别相等, 分别为0.5,0.75,1,1.25. 因此产生如下4组数据:
表1 不同样本数下检验问题(3)的势
表2~表4分别列出了一次抽样中参数的估计值和δj=0的后验概率, 其中Bayes估计是本文方法给出的估计, Frequentist估计是利用文献[10]中频率方法得到的参数估计值. 由表2可见, 对于参数δ1,δ2和δ3, 当n=30,100时, 本文方法都给出了比频率方法更精确的估计. 尽管当n=30时, 频率方法给出的参数ξ估计值比本文方法更接近真值, 但当n=100 时, 本文方法给出的估计值更接近真值. 由表3可见, 本文方法给出的参数ξ,δ1,δ3估计值优于频率方法给出的估计值, 但对于参数δ2, 频率方法给出了比本文方法更精确的估计. 由表4可见, 只有参数δ1当n=30时和参数δ3当n=100时, 频率方法优于本文方法, 其他情况都是本文方法优于频率方法. 总之, 两种方法给出的估计值类似, 都与真值很接近. 对于检验问题(3), 从表中数据可见应用本文给出的Bayes方法得到的结果较理想. 其中在表3中, 尽管当间隔δj较小且样本量不大(n=30)时出现了与实际不符的检验结果(P(δ2=0Y)=0.611 7), 但随着样本的增多或间隔的增大(表4), 从后验概率都小于0.5易见4个正态总体均值与标准差比都不同. 这主要是由于间隔较小时很难正确识别出所有不同的均值和标准差的比.
表2 当μ=(1,1,1,1), σ=(1,1,1,1)时一次抽样中参数的估计值和δj=0的后验概率
表3 当μ=(0,1,2,3), σ=(1,1,1,1)时一次抽样中参数的估计值和δj=0的后验概率
表4 当μ=(0,2,4,6), σ=(1,1,1,1)时一次抽样中参数的估计值和δj=0的后验概率
进一步, 当原假设成立时, 本文对犯第一类错误的概率进行了模拟, 结果列于表5. 考虑4组总体, 数据按如下总体产生:
其中σ=(σ1,σ2,σ3,σ4)分别取(1,1,1,1),(3,3,3,3),(6,6,6,6),(10,10,10,10).
表5 不同样本数和方差下犯第一类错误的概率
由表5可见, 本文方法能很好地控制犯第一类错误的概率, 甚至当标准差很大(σ=10)时, 犯第一类错误的概率仍然很小, 因此本文方法是一个相对保守的方法. 此外, 在表5中, 当σ=(1,1,1,1)时似乎随着样本量的增加犯第一类错误的概率也增加, 为了检验这一结论本文做了很多模拟, 结果显示当n=49时, 犯第一类错误的概率为0.005 0, 当n=60,70,80,90,110,120时, 犯第一类错误的概率均为0, 因此, 并不存在犯第一类错误的概率随样本量增加而增大的趋势. 当间隔δ不是太小时, 本文的Bayes方法在处理检验问题(3)时, 在控制犯第一类错误很小的同时, 也能达到一个合理的检验势.
[1] Robertson T, Wright F T, Dykstra R L. Order Restricted Statistical Inference [M]. New Jersey: John Wiley & Sons Inc, 1988.
[2] Silvapulle M J, Sen P K. Constrained Statistical Inference: Inequality, Order and Shape Restrictions [M]. New Jersey: John Wiley & Sons Inc, 2005.
[3] Hayter A J. A One-Sided Studentized Range Test for Testing against a Simple Order Alternative [J]. Journal of the American Statistical Association, 1990, 85(411): 778-785.
[4] LIU Lin. Simultaneous Statistical Inference for Monotone Dose-Response Means [D]. Newfoundlan: Memorial Universtiy of Newfoundland, 2001.
[5] LIU Lin, Lee C C, PENG Jian-an. Max-Min Multiple Comparison Procedure for Isotonic Dose-Response Curves [J]. Journal of Statistical Planning and Inference, 2002, 107(1/2): 133-141.
[6] 史宁中. 统计检验的理论与方法 [M]. 北京: 科学出版社, 2008.
[7] SHANG Jun-feng, Cavanaugh J E, Wright F T. A Bayesian Multiple Comparison Procedure for Order-Restricted Mixed Models [J]. International Statistical Review, 2008, 76(2): 268-284.
[8] Oh M S, Shin D W. A Unified Bayesian Inference on Treatment Means with Order Constraints [J]. Computational Statistics & Data Analysis, 2011, 55(1): 924-934.
[9] Chou Y M, Owen D B. A Likelihood Ratio Test for the Equality of Proportions of Two Normal Populations [J]. Communications in Statistics: Theory and Methods, 1991, 20(8): 2357-2374.
[10] LI Shu-you, SHI Ning-zhong, ZHANG Bao-xue. Testing Ratios of Means to Standard Deviations from Normal Populations under Order Restrictions with Generalizedp-Values [J]. Chinese Journal of Applied Probability and Statistics, 2009, 25(1): 77-85. (李树有, 史宁中, 张宝学. 正态总体均值与标准差比在序约束下的广义p-值检验 [J]. 应用概率统计, 2009, 25(1): 77-85.)