分類數(shù)據(jù)的統(tǒng)計分析及SAS編程

1、一一 . 概概 述述 n分 析 策 略 分類數(shù)據(jù)分析策略可以分成假設(shè)檢驗和建立模型。1.假設(shè)檢驗法是建立一個關(guān)于聯(lián)系（association）的假設(shè)。通常研究用隨機化的方法進行。 n例如： 把病人隨機分成兩組，檢驗組別與療效之間（列聯(lián)表的行與列之間）是否有關(guān)。這種聯(lián)系的無效假設(shè)為變量間沒有聯(lián)系，而備擇假設(shè)一般有3種情況： 有一般聯(lián)系（general association）。 平行均分有差別（row mean scores differ）。 非零相關(guān)（nonzero correlation）。 在以后討論中我們將對各種不同的聯(lián)系進行說明。2.建立模型的方法可求得各參數(shù)值，說明各因素的作用。通常

2、用最大似然估計或加權(quán)最小二乘法估計。 二二 . 22 表表 n臨床試驗實例 例2.1 在研究老年慢性支氣管炎的中草藥療效時，將病人隨機分配到江剪刀草組和膽麻片組。結(jié)果療效如下表所示： 表表2-1 兩組有效率的比較兩組有效率的比較 這是這是2.1中提到的第一種情況。根據(jù)無效假設(shè)，數(shù)據(jù)為超幾何分布。中提到的第一種情況。根據(jù)無效假設(shè)，數(shù)據(jù)為超幾何分布。 2.2.1 Mantel-Haenszel 卡方Q如果用以下符號表示表2-1中的數(shù)據(jù) 表表2-2 兩組有效率的比較兩組有效率的比較 則每一格出現(xiàn)數(shù)據(jù)為 的概率為： n如上例江剪刀草組無效（263例）的期望值為： 方差 ： 本例 ： 樣本足夠大時 近似

3、正太分布，由此得： 結(jié)論為兩組療效不同 。2.2.2 Pearson 卡方QpPearson 卡方即一般常用的 。這一公式即一般習(xí)慣用的 ： n計算Pearson卡方先要計算出4個期望值 ，一個 已算出為239.87。其他則可由合計值減出。Mantel-Haenszel 與Pearson 的關(guān)系如下： =本例 ： = 2.2.3 SAS 程序 程序2-1 data bronchit; input treat $ outcome $ count; cards; jiang y 70 jiang n 263 dan y 102 dan n 180 ; proc freq order = data;

4、 weight count; tables treat * outcome/chisq; run; 檢驗結(jié)果 ：STATISTICS FOR TABLE OF TREAT BY OUTCOME 表中各列依次為統(tǒng)計量、自由度（DF）、卡方值及P值。 2.7 配對資料四格表 2.7.1 例2.3 205份標本分別接種于甲、乙兩種培養(yǎng)基，共有4種結(jié)果，即甲、乙都生長；甲生長乙不生長，乙生長甲不生長及甲、乙都不生長。試比較兩種培養(yǎng)基的效果，結(jié)果如表2-5 ： 表2-5 兩種培養(yǎng)基的結(jié)果 2.7.2 計算方法 2.7.3 程序 在SAS第6.10版以后，在proc freq 中tables語句中加入了a

5、gree選項作McNemar檢驗。本例程序如下： 程序2-2 data culture; input a $ b $ count; cards; y y 36 y n 34 n y 0 n n 135 ； proc freq order = data; weight count; tables a * b/agree; run; 2.7.4 輸出 除頻數(shù)部分外，輸出的結(jié)果如下：STATISTICS FOR TABLE OF A BY B McNemars Test 其中，Kappa系數(shù)是反映一致性的指標。三三. 多層多層 22 表表3.1 概 述 在醫(yī)學(xué)研究中經(jīng)常遇到分層研究，如果每個層都有一

6、個22表，則有多個22表（sets of 22 table）。例如在多中心臨床試驗中，每個醫(yī)院隨機地把病人分為試驗組和對照組，療效為有效和無效，則每個醫(yī)院的數(shù)據(jù)形成了一個22表數(shù)據(jù)。 3.2 Mantel-Haenszel 檢驗 3.2.1 實例3.1 在三所醫(yī)院中考察某感冒藥A（江剪刀草組）對治療流鼻涕的效果，與對照藥B（膽麻片組）對比記錄其療效如表3-1：表表3-1 某感冒藥某感冒藥A與對照藥與對照藥B治療流鼻涕的效果治療流鼻涕的效果 3.2.2 符號 對于多層四格表，除四格表中數(shù)字以 代表外，以h表示層次，則第h層的四格表可用以下符號表示（表3-2）:表表3-2 多層四格表的符號表示多層

7、四格表的符號表示 3.2.3 QMH 統(tǒng)計量 在無效假設(shè)為兩組療效相同時，第h層 的期望值 及方差 的計算公式為： = (3.1) = (3.2) 校正中心因素的兩組之間療效的差異可用Mantel-Haenszel 1959年提出的統(tǒng)計量表示。 = 式中q為層數(shù)，本例q3。 當(dāng)q層的i行合計數(shù) 的合計數(shù) 大于30時近似服從自由度等于1的卡方分布。 nMantel-Haenszel方法，消除了層次因素的干擾而提高了檢出關(guān)聯(lián)的把握度。 當(dāng)各層次的陽性率的方向一致時，如每層的治愈率都是試驗組高于對照組，則QMH檢驗效果很好。如果各層的方向不一致則可能檢驗不出其關(guān)聯(lián)。 3.2.4 程序 程序3-1 d

8、ata park2; input center treat $ response $ count ; cards; 1 test y 12 1 test n 18 1 placebo y 15 1 placebo n 15 2 test y 31 2 test n 9 3 placebo y 34 2 placebo n 6 3 test y 16 3 test n 14 3 placebo y 15 3 placebo n 15 ； proc freq order = data; weight count; tables center * treat * response/chisq CMH

9、; run; tables語句中要把需要校正的因素，即center放在最前面，CMH選項可輸出Mantel-Haenszel統(tǒng)計量。 3.2.5 輸出 輸出內(nèi)容除包括各層的四格表頻數(shù)、百分數(shù)、X2檢驗、Fisher精確檢驗結(jié)果及相關(guān)等內(nèi)容外，還輸出了小結(jié)統(tǒng)計量（summary statistics）、OR及RR。 這里沒有列出四格表，但由結(jié)果可見，中心1江剪刀草組與膽麻片組療效比為40.00：50.00，以膽麻片組療效較好（表3-3）；中心2為77.50：85.00，也是以膽麻片組療效較差。就是各組療效比的方向在各中心間并不一致，這樣的資料QMH檢驗效果較差。表表3-3 各中心療效方向比較各中

10、心療效方向比較 CMH選項輸出了下列結(jié)果：SUMMARY STATISTICS FOR TREAT BY RESPONSECONTROLLING FOR CENTERn首先輸出Cochran-Mantel-Haenszel統(tǒng)計量，其中檢驗了非零相關(guān)，行平均得分差及一般關(guān)聯(lián)三個備選假設(shè)。這三者在四格表資料中結(jié)果是相同的，因而不必考慮應(yīng)當(dāng)選擇哪一個結(jié)果。這里是檢驗治療方法與療效的關(guān)系，但校正了中心的作用。有表中數(shù)據(jù)可見，QMH=0.586，P=0.44，不能得出兩種療法療效有所不同的結(jié)論。 n輸出的總的OR及RR的可信區(qū)間都包括了1在內(nèi)，因而并無統(tǒng)計意義，結(jié)論與上述相同。n輸出的Breslow-D

11、ay檢驗是檢驗治療組與療效的關(guān)聯(lián)在各中心是否一致。本例QMH=0.821,DF=2,P=0.663,即雖然我們看到中心1，中心2與中心3治療組與療效的關(guān)系有所不同，但還不足以說明其各中心的OR不同。如果Breslow-Day QMH有統(tǒng)計意義則不能用總的OR表示總的療效比較。四四. 2r表和多層表和多層2r表表 4.1 2r表 2行多列的列聯(lián)表為2r表；多行2列的表為s2表。如為分層分析，則可有多層2r表或s2表。4.1.1 例4.1在一種感冒藥片的研究中對咳嗽的療效如表4-1所示：表表4-1 兩組病人咳嗽療效的比較兩組病人咳嗽療效的比較4.1.2 行平均得分 這種資料用一般的卡方檢驗只能得出

12、兩組構(gòu)成是否相同的結(jié)論，并不能得出哪組療效較好的結(jié)論。例如，當(dāng)試驗組無效和痊愈較多，對照組則好轉(zhuǎn)和顯效者較多，這時X2則會較大，P則較大，說明其構(gòu)成不同，但不能說明何者療效好。 一種辦法是人為地給各療效一個分數(shù)，如無效為1、好轉(zhuǎn)為2、顯效為3、痊愈為4，則可計算其均數(shù)，稱為行平均得分（row mean scores）。公式為： (4.1) 式中為各療效的得分， 為第一行的各療效頻數(shù)， 為第一行合計。由此求得第一行的平均得分 。本例試驗組行平均得分為： 當(dāng)然也可求出第二行平均得分，但這里作檢驗時只需第一行平均得分。4.1.3 行平均得分差檢驗 以第一行平均得分與期望得分之差的平方除以方差 ,可得

13、檢驗統(tǒng)計量 : = (4.2) 式中 為平均期望得分， 為方差，可由列合計得出即： (4.3) (4.4) 本例 在以上計算中略去了小數(shù)點后的一些數(shù)字。如果計算更精確，計算結(jié)果應(yīng)為： 近似服從自由度為1的卡方分布。 4.1.4 SAS程序 上例可編程如下： 程序中在tables語句后增加了MCH選項，即求Cochran-Mantel-Haenszel統(tǒng)計量。4.1.5 輸出 以上程序輸出了列聯(lián)表中頻數(shù)及各百分比，與前相同，這里從略。 以下輸出為及相關(guān)系數(shù)。這里的Mantel-Haensel 卡方即 值。 ，與前面精確計算結(jié)果相同 4.2 結(jié)果尺度為次數(shù)的數(shù)據(jù)結(jié)果尺度為次數(shù)的數(shù)據(jù) 當(dāng)2r表中列數(shù)

14、據(jù)為次數(shù)時由于沒有理由認為次數(shù)是正太分布和方差齊性，因而不能用方差分析，而需用CMH分析。這是一個反應(yīng)變量為另一種數(shù)據(jù)類型的例子。 4.2.1 例 4.2 調(diào)查小學(xué)生一學(xué)期中感冒情況，得病學(xué)生感冒次數(shù)如表4-2所示，試研究性別間感冒次數(shù)有無不同。表表4-2 學(xué)生一學(xué)期中感冒次數(shù)學(xué)生一學(xué)期中感冒次數(shù) 4.2.2 程序 程序4-2 data cold; input sex cold count ; cards; male 0 45 male 1 64 male 2 71 female 0 84 female 1 124 female 2 82 ； proc freq order = data; w

15、eight count; tables sex * cold/all; run; 程序中all選項包括了chisq及CMH項。4.2.3 輸出 輸出格式與上例相同。結(jié)果Mantel-Haenszel 卡方及row mean scores differ卡方都是4.249，DF1，P0.039。性別之間感冒次數(shù)差別由統(tǒng)計意義。4.3 多層多層2r表表在分層分析時，每個層有一個2r表則形成多層2r表。4.3.1 例 4.3 在以上研究感冒咳嗽的療效這中，實際是一個多中心研究。每個醫(yī)院（中心層）有一個2r表.結(jié)果如表4-3：表表4-3 三中心感冒咳嗽療效三中心感冒咳嗽療效 4.3.2 符號 以h代表層

16、的符號，共有g(shù)層，則每層2r表的符號為（表4-1）：表表4-4 2r表的符號表示表的符號表示 4.3.3 擴展的Mantel-Haenszel 平均得分統(tǒng)計量各層第一行總得分的合計 由下式求出： (4.7) 第h層第一行平均得分為： (4.8) 在 （兩組行平均得分相等）成立的前提下 的期望值 為： (4.9) 式中 為各層第一行得分期望值 (4.10)h層第一行平均得分的方差為： (4.11) 總方差 ： (4.12) 如果各層各行的合計數(shù) 足夠大，則 近似正態(tài)分布，統(tǒng)計量為 ： (4.13) 服從自由度為1的卡方分布。 4.3.4 程序本例可編程序如下: 程序4-3 data cough2

17、r; input center $ treat $ response $count ; cards; one test none 7 one test improve 11 one test notable 8 one test cure 4 one control none 2 one control improve 13 one control notable 13 one control cure 2 two test none 2 two test improve 7 two test notable 22 two test cure 9 two control none 0 two

18、control improve 6 two control notable 20 two control cure 14 three test none 5 three test improve 9 three test notable 13 three cure 3 three control none 6 three control improve 9 three control notable 15 three control cure 0 ； proce freq order = data; weight count; tables center * treat * response/

19、cmh; run; 程序中輸入中心，治療方法、結(jié)果及計算，用order=data保證輸出時與輸入順序相同，計數(shù)用weight count;語句。Tables語句中有三個變量其中第一個center為分層變量，treat * response為各層2r表。分析用CMH。 4.3.5 輸出 各中心皆輸出一個2r列聯(lián)表的頻數(shù)及百分比（略）。然后輸出CMH，統(tǒng)計量如下: 值即表中的Row Mean Scores Differ。本例 ， 。在校正中心效應(yīng)后兩組的平均得分差沒有統(tǒng)計學(xué)意義。輸出的CMH統(tǒng)計量后面的注明為以表中得分為基礎(chǔ)。也即是4.1.6中所述的整數(shù)分數(shù)。 如果反應(yīng)變量的級別不是等間隔的則可以

20、選用修改的ridit法計算，可寫成以下tables語句： tables center * treat * response/CMH scores = modridit; 與前述多層2r表相同，各層間效應(yīng)的方向一致時 效果較好。 五五. 非參數(shù)方法的非參數(shù)方法的CMH解法解法 非參數(shù)統(tǒng)計方法，如Wilcoxon-Mann-Whitney，Kruskal-Wallis，Spearman Correlation，F(xiàn)reidman和Durbin檢驗，能用Mantel-Haenszel過程，并且可以解決傳統(tǒng)方法不能解決的問題（例如分中心資料的非參數(shù)統(tǒng)計）等。6.1 兩樣本秩和檢驗兩樣本秩和檢驗（Wilc

21、oxon-Mann-Whitney法）法）6.1.1 對于兩個樣本差異統(tǒng)計意義的非參數(shù)檢驗 1945年首先由Wilcoxon提出，主要是針對兩樣本含量相同的情況。1947年，Mann和Whitney又在考慮到不等樣本的情況下補充了這一方法。它的計算公式為: （6.1） m、n分別為兩樣本的含量，規(guī)定mn T1兩樣本混合由小到大排列編秩，較小樣本的秩和 T2兩樣本混合有大到小排列編秩，較小樣本的秩和 編秩時，要求把兩個樣本混合由小到大排列統(tǒng)一編秩，相同的數(shù)據(jù)一律為平均秩次。當(dāng)m15，n28時，計算結(jié)果可以和附表（等級總和數(shù)臨界值）進行比較。如m、n超過附表范圍，便需要用公式（6.3）進行u檢驗。

22、 對于差異具有統(tǒng)計意義的兩組成正態(tài)分布的樣本來說，W-M-W檢驗相對于兩樣本的t檢驗的漸進效率是0.955，而對于總體非正態(tài)（例如非對稱分布），W-M-W檢驗比兩樣本t檢驗的效率高得多，事實上這時的漸進相對效率能高達無窮小，所以W-M-W的方法對于兩樣本的檢驗是十分適用的。6.1.2 例 6.1 表表6-1 不同飼料的兩組雌鼠在不同飼料的兩組雌鼠在8周內(nèi)所增體重周內(nèi)所增體重 先將兩組數(shù)據(jù)混合從小到大排列，并注明組別與等級（如數(shù)字相同則取先將兩組數(shù)據(jù)混合從小到大排列，并注明組別與等級（如數(shù)字相同則取平均等級）。平均等級）。表表6-2 兩樣本等級和檢驗用表兩樣本等級和檢驗用表 較小樣本為低蛋白組較

23、小樣本為低蛋白組m7 以與兩者中的較小者與附表以與兩者中的較小者與附表1數(shù)據(jù)進行比較，當(dāng)數(shù)據(jù)進行比較，當(dāng)m=7，n=12時查附時查附表：等級總和數(shù)臨界值，得表：等級總和數(shù)臨界值，得P 0.05臨界值為臨界值為46，則，則P0.05，差別，差別無統(tǒng)計意義。無統(tǒng)計意義。6.2 完全隨機化設(shè)計資料的檢驗（完全隨機化設(shè)計資料的檢驗（H檢驗）檢驗） Kruskal-Wallis檢驗是一個推廣兩樣本W(wǎng)-M-W檢驗到3個或更多檢驗的方法，又稱為H檢驗方法的前提是假定抽樣總體是連續(xù)的和相同的，檢驗其各樣本之間的分布位置是否相同。對適合方差分析的資料，K-W檢驗也同樣適用。它的計算公式為: K樣本數(shù) 各樣本例數(shù)

24、 各樣本的秩和 相等秩次的個數(shù) 校正H值 nk個樣本例數(shù)之和 (總例數(shù)) C 校正數(shù) 將各組數(shù)據(jù)混合由小到大排列，列出等級，相同數(shù)量則取平均等級。當(dāng)進行排序時無相同數(shù)據(jù)，則計算H值。如在排序時出現(xiàn)相同數(shù)據(jù)，則計算校正H值 ( ) 當(dāng)這些組的樣本不大時，即k3， 5時，可查附表2（H值與概率對照表）進行判斷。如果超出附表范圍，根據(jù) 不太小時H近似于自由度為k-1的 分布，可以查相應(yīng)的 界值表而得到近似的P值，其中k為處理組數(shù)。這種近似方法只是M-H平均得分統(tǒng)計量對一個層采用等級（rank）得分時的特殊情況。 6.2.2 例 6.2 50只小鼠隨機分配到5個不同飼料組，每組10只小鼠。在喂養(yǎng)一定時

25、間后測得鼠肝中鐵的含量（ /g）如表6-3所示（括號中為等級數(shù)）。試檢驗各組含量的差別有無統(tǒng)計意義。檢驗步驟：編秩表表6-3 不同飼料組小鼠肝臟中鐵（不同飼料組小鼠肝臟中鐵（g/g）含量及其等級）含量及其等級 (2).計算各組等級和 (3).計算各 3667.225+7924.225+17640+970.225+8037.225 =38238.900 (4).按公式計算H值 本例有相同值，所以應(yīng)該計算 值。 因為本例k5，自由度等于514，查 ，H P0.001， 差別有高度統(tǒng)計意義。 6.2.3 SAS程序中CMH的使用 對于成組比較的K-W非參數(shù)檢驗的方法，在SAS系統(tǒng)中可以通過以下的程序

26、得到解決，該程序是在求等級（rank）的基礎(chǔ)上是使用Mantel-Haenszel平均分差檢驗的方法。 程序6-2 data mouse; input group quality ; cards; 1 2.23 2 5.59 3 4.50 4 1.35 5 1.40 1 1.14 2 0.96 3 3.92 4 1.06 5 1.06 1 2.63 2 6.96 3 10.33 4 0.74 5 2.49 1 1.00 2 1.23 3 8.23 4 0.96 5 1.74 1 1.35 2 1.61 3 2.07 4 1.16 5 1.59 1 2.01 2 2.94 3 4.90 4 2.08 5 1.36 1 1.64 2 1.96 3 6.84 4 0.69 5 3.00 1 1.13 2 3.68 3 6.42 4 0.68 5 4.81 1 1.01 2 1.54 3 3.72 4 0.84 5 5.21 1 1.70 2 2.59 3 6.00 4 1.34 5 5.12 ； proce freq; tables group * weight/noprint CMH2 scores = rank; run; 經(jīng)SAS運行結(jié)果如下參看第二行row mean sc