多組平均水平的比較

1、Stata軟件基本操作和數(shù)據(jù)分析入門第五講 多組平均水平的比較復(fù)習(xí)和補(bǔ)充兩組比較的統(tǒng)計(jì)檢驗(yàn)配對(duì)設(shè)計(jì)資料(又稱為Dependent Samples)對(duì)于小樣本的情況下，如果配對(duì)的差值資料服從正態(tài)分布，用配對(duì)t檢驗(yàn)(ttest 差值變量=0)大樣本的情況下，可以用配對(duì)t檢驗(yàn)小樣本的情況下，并且配對(duì)差值呈偏態(tài)分布，則用配對(duì)符號(hào)秩檢驗(yàn)(signrank 差值變量=0)成組設(shè)計(jì)(Two Independent Samples)如果方差齊性并且大樣本情況下，可以用成組t檢驗(yàn)(ttest 效應(yīng)指標(biāo)變量,by(分組變量)如果方差齊性并且兩組資料分別呈正態(tài)分布，可以用成組t檢驗(yàn)如果方差不齊，或者小樣本情況下偏態(tài)

2、分布，則用秩和檢驗(yàn)(Ranksum test)groupx0790930910920940770930740910101083073088010209001000810910830106084078087095010111011100111418611061107110719418911041981110189110318911211941951921109198198112011041110多組比較完全隨機(jī)分組設(shè)計(jì)(要求各組資料之間相互獨(dú)立)方差齊性并且獨(dú)立以及每一組資料都服從正態(tài)分布(小樣本時(shí)要求)，則采用完全隨機(jī)設(shè)計(jì)的方差分析方法(即：?jiǎn)我蛩胤讲罘治?，One Way ANOVA)進(jìn)行分析

3、。方差不齊或小樣本情況下資料偏態(tài)，則用Kruskal Wallis 檢驗(yàn)(H檢驗(yàn))例5.1 為研究胃癌與胃粘膜細(xì)胞中DNA含量（A.U）的關(guān)系，某醫(yī)師測(cè)得數(shù)據(jù)如下，試問(wèn)四組人群的胃粘膜細(xì)胞中平均DNA含量是否相同？組別groupDNA含量（A.U）淺表型胃炎19.8112.7312.2912.5312.959.5312.68.912.2714.2610.68腸化生214.6117.5415.11713.3915.3213.7418.2413.8112.6314.5316.17早期胃癌323.2620.820.623.517.8521.9122.1322.0419.5318.4121.4820.

4、24晚期胃癌423.7319.4622.3919.5325.920.4320.7120.0523.4121.3421.3825.70由于這四組對(duì)象的資料是相互獨(dú)立的，因此屬于完全隨機(jī)分組類型的。檢驗(yàn)問(wèn)題是考察四組DNA含量的平均水平相同嗎。如果每一組資料都正態(tài)分布并且方差齊性可以用One way-ANOVA進(jìn)行分析，反之用Kruskal Wallis檢驗(yàn)。STATA數(shù)據(jù)輸入格式gx19.81112.73112.29112.53112.9519.53112.618.9112.27114.26110.68214.61217.54215.1217213.39215.32213.74218.24213

5、.81212.63214.53216.17323.26320.8320.6323.5317.85321.91322.13322.04319.53318.41321.48320.24423.73419.46422.39419.53425.9420.43420.71420.05423.41421.34421.38425.7分組正態(tài)性檢驗(yàn)，=0.05. sktest x if g=1 Skewness/Kurtosis tests for Normality - joint - Variable | Pr(Skewness) Pr(Kurtosis) adj chi2(2) Probchi2-+-

6、x | 0.491 0.485 1.07 0.5861. sktest x if g=2 Skewness/Kurtosis tests for Normality - joint - Variable | Pr(Skewness) Pr(Kurtosis) adj chi2(2) Probchi2-+- x | 0.482 0.541 0.96 0.6201. sktest x if g=3 Skewness/Kurtosis tests for Normality - joint - Variable | Pr(Skewness) Pr(Kurtosis) adj chi2(2) Prob

7、chi2-+- x | 0.527 0.750 0.52 0.7704. sktest x if g=4 Skewness/Kurtosis tests for Normality - joint - Variable | Pr(Skewness) Pr(Kurtosis) adj chi2(2) Probchi2-+- x | 0.260 0.616 1.75 0.4166上述結(jié)果表明每一組資料都服從正態(tài)分布。單因素方差分析的STATA命令：oneway 效應(yīng)指標(biāo)變量 分組變量,t b其中t表示計(jì)算每一組均數(shù)和標(biāo)準(zhǔn)差,b表示采用Bonferroni統(tǒng)計(jì)方法進(jìn)行兩兩比較。本例命令為oneway

8、 x group,t b. oneway x g,t b | Summary of x g | Mean Std. Dev. Freq.-+- 1 | 11.686364 1.6884388 11 2 | 15.173333 1.749173 12 3 | 20.979167 1.7668279 12 4 | 22.0025 2.2429087 12-+- Total | 17.583191 4.6080789 47 Analysis of Variance Source SS df MS F Prob F-Between groups 824.942549 3 274.98085 77.87

9、 0.0000 Within groups 151.839445 43 3.53114987- Total 976.781994 46 21.2343912Bartletts test for equal variances: chi2(3) = 1.1354 Probchi2 = 0.769方差齊性的檢驗(yàn)為：卡方=1.1354，自由度=3,P值=0.769，因此可以認(rèn)為方差是齊性的。H0：1=2=3=4 四組總體均數(shù)相同H1：1，2，3，4不全相同=0.05，相應(yīng)的統(tǒng)計(jì)量F=77.87以及相應(yīng)的自由度為3和43，P值 g = 1 - Binom. Interp. - Variable | O

10、bs Percentile Centile 95% Conf. Interval-+- x | 11 50 12.29 9.729564 12.7932- g = 2 - Binom. Interp. - Variable | Obs Percentile Centile 95% Conf. Interval-+- x | 12 50 14.855 13.74745 16.91172- g = 3 - Binom. Interp. - Variable | Obs Percentile Centile 95% Conf. Interval-+- x | 12 50 21.14 19.60552

11、 22.12043- g = 4 - Binom. Interp. - Variable | Obs Percentile Centile 95% Conf. Interval-+- x | 12 50 21.36 20.09042 23.69596得到這4組中位數(shù)分別為：M1=12.29，M2=14.855，M3=21.14和M4=21.36ranksum x if g=1 | g=2,by(g)Two-sample Wilcoxon rank-sum (Mann-Whitney) test g | obs rank sum expected-+- 1 | 11 72 132 2 | 12

12、204 144-+- combined | 23 276 276unadjusted variance 264.00adjustment for ties 0.00 -adjusted variance 264.00Ho: x(g=1) = x(g=2) z = -3.693 Prob |z| = 0.0002P值M1)，并且差別有統(tǒng)計(jì)學(xué)意義。第1組與第3組比較ranksum x if g=1 | g=3,by(g)Two-sample Wilcoxon rank-sum (Mann-Whitney) test g | obs rank sum expected-+- 1 | 11 66 13

13、2 3 | 12 210 144-+- combined | 23 276 276unadjusted variance 264.00adjustment for ties 0.00 -adjusted variance 264.00Ho: x(g=1) = x(g=3) z = -4.062 Prob |z| = 0.0000P值M1)，并且差別有統(tǒng)計(jì)學(xué)意義，其他比較類似進(jìn)行。要注意的問(wèn)題：在方差分析中，要求每一組資料服從正態(tài)分布(小樣本時(shí))，并不是要求各組資料服從一個(gè)正態(tài)分布(因?yàn)檫@就意味各組的總體均數(shù)相同，失去統(tǒng)計(jì)檢驗(yàn)的必要性)，所以不能把各組的資料合在一起作正態(tài)性檢驗(yàn)?？偟闹v，方差分析

14、對(duì)正態(tài)性具有穩(wěn)健性，即：偏態(tài)分布對(duì)方差分析的結(jié)果影響不會(huì)太大，故正態(tài)性檢驗(yàn)的取0.05也就可以了。樣本量較大時(shí)，方差分析對(duì)正態(tài)性要求大大降低(根據(jù)中心極限定理可知：樣本均數(shù)近似服從正態(tài)分布)。并且由于大多數(shù)情況下，樣本資料只是近似服從正態(tài)分布而不是完全服從正態(tài)分布。由于在大樣本情況下，用正態(tài)性檢驗(yàn)就變?yōu)楹苊舾?，?duì)于不是完全服從正態(tài)分布的資料往往會(huì)拒絕正態(tài)性檢驗(yàn)的H0：資料服從正態(tài)分布。因?yàn)檎龖B(tài)性檢驗(yàn)不能檢驗(yàn)資料是否近似服從正態(tài)分布，而是檢驗(yàn)是否服從正態(tài)分布。故在大樣本情況下，考察資料的近似正態(tài)性，應(yīng)用頻數(shù)圖進(jìn)行考察。方差齊性問(wèn)題對(duì)方差分析相對(duì)比較敏感，并且并不是隨著樣本量增大而方差齊性對(duì)方差分

15、析減少影響的。但是當(dāng)各組樣本量接近相同或相同時(shí)，方差齊性對(duì)方差分析呈現(xiàn)某種穩(wěn)健性。即：只有當(dāng)各組樣本量相同時(shí)，方差齊性對(duì)方差分析結(jié)果的影響大大降低。這時(shí)隨著樣本量增大，影響會(huì)進(jìn)一步降低。相反，如果各組樣本量相差太大時(shí)，方差齊性對(duì)方差分析結(jié)果的影響很大。這時(shí)隨著樣本量增大，影響會(huì)進(jìn)一步加大。隨機(jī)區(qū)組設(shè)計(jì)(處理組之間可能不獨(dú)立)殘差(定義為：，也就是隨機(jī)區(qū)組方差分析中的誤差項(xiàng))的方差齊性且小樣本時(shí)正態(tài)分布，則用隨機(jī)區(qū)組的方差分析(無(wú)重復(fù)的兩因素方差分析,Two-way ANOVA)。不滿足方差齊性或小樣本時(shí)資料偏態(tài)，則對(duì)用秩變換后再用隨機(jī)區(qū)組的方差分析也可以直接用非參數(shù)隨機(jī)區(qū)組的秩和檢驗(yàn)Fredm

16、an test)。例2下表是某湖水中8個(gè)觀察地點(diǎn)不同季節(jié)取樣的氯化物含量測(cè)定值，請(qǐng)問(wèn)在不同季節(jié)該湖水中氯化物的含量有無(wú)差別？表2 某湖水中不同季節(jié)的氯化物含量測(cè)定值（mg/L）location no春夏秋冬121.2818.3317.2714.91222.7819.8116.5514.85320.9018.9316.3616.30419.9021.2317.8615.73521.4919.0915.1117.05622.3817.9216.5714.34721.6719.3917.1916.31822.0619.6516.5814.33顯然同一地點(diǎn)不同季節(jié)的氯化物含量有一定的相關(guān)性，故不能采用

17、完全隨機(jī)設(shè)計(jì)的方差分析方法對(duì)4個(gè)季節(jié)的氯化物含量進(jìn)行統(tǒng)計(jì)分析。可以把同一地點(diǎn)的4個(gè)季節(jié)氯化物含量視為一個(gè)區(qū)組，因此可以用隨機(jī)區(qū)組的方差分析進(jìn)行統(tǒng)計(jì)分析。設(shè)第8個(gè)地點(diǎn)在冬季的氯化物總體均數(shù)為0，同樣在冬季，第i個(gè)地點(diǎn)的氯化物總體均數(shù)與第8個(gè)地點(diǎn)在冬季的氯化物總體均數(shù)相差i，i=1，2，3，4，5，6，7。因此在冬季的這8個(gè)地點(diǎn)在冬季的氯化物總體均數(shù)可以表示為地點(diǎn)編號(hào)12345678冬季氯化物均數(shù)0+10+20+30+40+50+60+70假定在同一地區(qū)，春季的氯化物總體均數(shù)與冬季的氯化物總體均數(shù)相差1，因此春節(jié)和冬季的氯化物總體均數(shù)可以表示為地點(diǎn)編號(hào)12345678冬季氯化物均數(shù)0+10+20+

18、30+40+50+60+70春季氯化物均數(shù)0+1+10+1+20+1+30+1+40+1+50+1+60+1+70如果10說(shuō)明在同一地點(diǎn)，冬季和春季的氯化物總體均數(shù)相同；10說(shuō)明春季的氯化物含量平均高于冬季氯化物含量，反之 F -+- Model | 180.214326 10 18.0214326 17.40 0.0000 | t | 177.344737 3 59.1149122 57.08 0.0000 id | 2.86958916 7 .409941308 0.40 0.8942 | Residual | 21.749618 21 1.0356961 -+- Total | 201.

19、963944 31 6.51496593 處理效應(yīng)H0：1230的檢驗(yàn)對(duì)應(yīng)的統(tǒng)計(jì)量相應(yīng)的P值 F = 0.0000 Residual | 21.749618 21 1.0356961 R-squared = 0.8923-+- Adj R-squared = 0.8410 Total | 201.963944 31 6.51496593 Root MSE = 1.0177- x Coef. Std. Err. t P|t| 95% Conf. Interval-_cons (0) 15.37992 .5966746 25.78 0.000 14.13906 16.62077t1= 1 6.08

20、0619 .5088458 11.95 0.000 5.022417 7.1388222= 2 3.816041 .5088458 7.50 0.000 2.757838 4.8742443= 3 1.20765 .5088458 2.37 0.027 .1494472 2.265853 4 (dropped)id1= 1 -.2092595 .7196166 -0.29 0.774 -1.705784 1.2872652= 2 .3387067 .7196166 0.47 0.643 -1.157818 1.8352313= 3 -.034339 .7196166 -0.05 0.962 -

21、1.530864 1.4621864= 4 .5231357 .7196166 0.73 0.475 -.973389 2.019665= 5 .0314307 .7196166 0.04 0.966 -1.465094 1.5279556= 6 -.353369 .7196166 -0.49 0.628 -1.849894 1.1431567= 7 .4840407 .7196166 0.67 0.509 -1.012484 1.980565 8 (dropped)其中，對(duì)應(yīng)的假設(shè)檢驗(yàn)H0：1=0的統(tǒng)計(jì)量t=11.95，P值0.001，95%可信區(qū)間為(5.022，7.139)，因此可以認(rèn)為

22、春季的氯化物平均高于冬季，差別有統(tǒng)計(jì)學(xué)意義。，對(duì)應(yīng)的假設(shè)檢驗(yàn)H0：2=0的統(tǒng)計(jì)量t=7.50，P值2、1=2和1 F = 0.0000相應(yīng)的統(tǒng)計(jì)量F=26.28，P值2的估計(jì)值，所以可以認(rèn)為春季氯化物平均高于夏季的氯化物含量。同理檢驗(yàn)H0：13 vs H1：13，只需輸入命令test bt1=_bt3檢驗(yàn)H0：23 vs H1：23，只需輸入命令test bt2=_bt3此處不在詳細(xì)敘述了。由于隨機(jī)區(qū)組方差分析要求殘差()服從正態(tài)分布，再輸入regress以后，只要輸入predict 殘差變量名,residual，就可以得到殘差計(jì)算值。本例用e表示殘差變量名，因此輸入predict e,res

23、idual就可以得到殘差計(jì)算值e，然后對(duì)殘差進(jìn)行正態(tài)性檢驗(yàn)(sktest 殘差變量名)本例輸入命令為: sktest e 結(jié)果如下（H0：殘差服從正態(tài)分布 vs H1：殘差偏態(tài)分布, 0.05） Skewness/Kurtosis tests for Normality - joint - Variable | Pr(Skewness) Pr(Kurtosis) adj chi2(2) Probchi2-+- e | 0.699 0.586 0.46 0.7948P值=0.93349，因此可以認(rèn)為資料近似服從正態(tài)分布。(大樣本時(shí)，可以不考慮正態(tài)性問(wèn)題)如果資料呈偏態(tài)分布，可以對(duì)資料進(jìn)行秩變換(

24、Rank Transform)后，然后把變換后的秩視為原始數(shù)據(jù)進(jìn)行隨機(jī)區(qū)組的方差分析。秩變換的STATA命令為 egen 秩變量名rank(觀察變量名),by(區(qū)組變量)為了說(shuō)明上述操作分析的過(guò)程，故借用本例資料進(jìn)行秩變換操作說(shuō)明如下(本例資料正態(tài)分布，無(wú)需用秩變換，只是說(shuō)明操作而言).設(shè)用r表示秩變量名，則本例操作為egen r=rank(x) ,by(id) 產(chǎn)生秩ranova 命令anova r t id 結(jié)果如下 Number of obs = 32 R-squared = 0.9125 Root MSE = .408248 Adj R-squared = 0.8708 Source

25、| Partial SS df MS F Prob F-+- Model | 36.50 10 3.65 21.90 0.0000 | t | 36.50 3 12.1666667 73.00 0.0000 id | 0.00 7 0.00 0.00 1.0000 | Residual | 3.50 21 .166666667 -+- Total | 40.00 31 1.29032258 命令regress 結(jié)果如下regress Source | SS df MS Number of obs = 32-+- F( 10, 21) = 21.90 Model | 36.50 10 3.65 Prob F = 0.0000 Residual | 3.50 21 .166666667 R-squared = 0.9125-+- Adj R-s