胡純嚴(yán) ，胡良平 ，2*

（1.軍事科學(xué)院研究生院，北京 100850；2.世界中醫(yī)藥學(xué)會聯(lián)合會臨床科研統(tǒng)計學(xué)專業(yè)委員會，北京 100029*通信作者：胡良平，E-mail：lphu927@163.com）

在橫斷面設(shè)計的二維列聯(lián)表資料或稱為R×C列聯(lián)表資料中，根據(jù)兩變量是否為有序變量，可以歸納出三種不同的二維列聯(lián)表資料：①雙向無序R×C列聯(lián)表資料；②結(jié)果變量為有序變量R×C列聯(lián)表資料；③雙向有序且屬性不同R×C列聯(lián)表資料。對其進(jìn)行統(tǒng)計分析的方法分別為“獨立性假設(shè)檢驗”“秩和檢驗”和“相關(guān)分析”，它們可以被概括成一種統(tǒng)計分析方法，即CMH χ2檢驗。本文介紹廣義CMH χ2檢驗統(tǒng)計量及其三種變形，并基于SAS軟件實現(xiàn)統(tǒng)計計算。

1 適合進(jìn)行CMH χ2檢驗的三種R×C列聯(lián)表資料的實例

1.1 雙向無序R×C列聯(lián)表資料的實例

【例1】文獻(xiàn)［1］中有一個雙向無序R×C表資料（說明：不考慮“時間”的有序性），見表1，試分析“條目”與“時間”之間是否存在關(guān)聯(lián)性。

表1 基層精防醫(yī)護(hù)人員K6評定結(jié)果

1.2 結(jié)果變量為有序變量R×C列聯(lián)表資料的實例

【例2】文獻(xiàn)［2］給出了如下臨床資料：比較三種對肥厚性鼻炎（HR）治療方法的效果。將162例HR患者分為三組，A組（n=56）行下鼻甲骨黏骨膜下切除術(shù)，B組（n=43）行下鼻甲部分切除術(shù)，C組（n=63）行下鼻甲黏膜下微波熱凝術(shù)。治療效果見表2。試對三組治療效果進(jìn)行比較。

表2 三組下鼻甲手術(shù)治療的臨床效果

1.3 雙向有序R×C列聯(lián)表資料的實例

【例3】在文獻(xiàn)［3］中，為探討新冠肺炎疫情期間居民急性應(yīng)激障礙（ASD）癥狀的檢出情況及影響因素，得到有效問卷16 048份，有效問卷回收率為72.83%。在該研究中，因素“媒體暴露情況（時間）”與結(jié)果變量（出現(xiàn)ASD癥狀的程度）之間的數(shù)量關(guān)系如表3所示。試分析“媒體暴露情況（時間）”與“ASD癥狀的程度”之間是否存在線性相關(guān)關(guān)系。

表3 媒體暴露情況（時間）與ASD癥狀檢出情況的調(diào)查結(jié)果

2 CMH χ2檢驗的檢驗統(tǒng)計量及其三種變形

2.1 “CMH”的含義

“CMH”是“Cochran-Mantel-Haenszel”的縮寫，即采用三位作者姓名的第一個字母來命名該分析方法 ，實際上，Cochran、Mantel、Haenszel、Mantel、Birch、Landis、Heyman和Koch都對該分析方法做出過貢獻(xiàn)［4］。

2.2 廣義CMH χ2檢驗統(tǒng)計量的定義

一般來說，CMHχ2檢驗統(tǒng)計量是為分析高維列聯(lián)表資料而構(gòu)造出來的。將高維列聯(lián)表資料按某一個或多個因素進(jìn)行分層，每層都應(yīng)該是規(guī)模相同的R×C列聯(lián)表資料。下面以表4形式表示第h層的R×C表，h=1、2、…、q。q為層數(shù)，R為行數(shù)，C為列數(shù)。

表4 第h層R×C列聯(lián)表的列表格式

廣義CMHχ2檢驗統(tǒng)計量［4］定義如下：

需注意的是，當(dāng)各層間效應(yīng)方向不一致時，CMHχ2檢驗統(tǒng)計量的檢驗功效很低。

2.3 廣義CMH χ2檢驗統(tǒng)計量的三種變形

2.3.1 概述

式（1）是以矩陣形式呈現(xiàn)的計算公式，很不直觀。由于它試圖達(dá)到高維列聯(lián)表資料的多種不同分析目的（例如獨立性分析、相關(guān)性分析、差異性分析），故其隱含諸多前提條件，包括“是否存在有序變量”“如何給有序變量賦值”。基于不同的前提條件，就會產(chǎn)生出不同的統(tǒng)計計算公式。為了直觀起見，現(xiàn)假定只有一層，即一個二維R×C列聯(lián)表資料，基于不同的前提條件，將式（1）轉(zhuǎn)變成三種具體的、直觀的檢驗統(tǒng)計量。

2.3.2 與備擇假設(shè)為“非零相關(guān)”對應(yīng)的檢驗統(tǒng)計量

若擬分析的資料如本文表3（即雙向有序R×C列聯(lián)表資料）所示，此時式（1）就簡化成下式：

2.3.3 與備擇假設(shè)為“行平均評分不同”對應(yīng)的檢驗統(tǒng)計量

若擬分析的資料如本文表2（即結(jié)果變量為有序變量R×C列聯(lián)表資料）所示，此時式（1）就簡化成下面兩個式子之一：

在式（3）中，“∝”代表“呈正比”；當(dāng)給R×C列聯(lián)表的列變量（即結(jié)果變量）賦值或評分時，若按“表評分”方式評分，就采用式（3）（即單因素R水平設(shè)計一元定量資料方差分析）計算［4，7］；否則［包括“rank（秩）評分”“Ridit評分”和“修正的Ridit評分”］，就采用式（4）（即單因素R水平設(shè)計一元定量資料Kruskal-Wallis’s H秩和檢驗）計算［4，7］。在式（4）中，ni與Mi分別代表第i組中的“樣本含量”與“秩和”。

在式（3）中，可對F統(tǒng)計量進(jìn)行變換，見下式：

通過對式（5）變形，可以得到下式［8-9］：

式（6）的含義是：采用單因素R水平設(shè)計一元定量資料方差分析得到的檢驗統(tǒng)計量F值后，乘以組間自由度df組間所得到的結(jié)果，近似服從自由度為df組間的χ2分布。

2.3.4 與備擇假設(shè)為“一般關(guān)聯(lián)性”對應(yīng)的檢驗統(tǒng)計量

若擬分析的資料如本文表1（即雙向無序R×C列聯(lián)表資料）所示，此時式（1）就簡化成下式：

3 CMH χ2檢驗的SAS實現(xiàn)

【例4】沿用例3的資料，試將其分別視為“雙向有序R×C列聯(lián)表資料”“結(jié)果變量為有序變量R×C列聯(lián)表資料”和“雙向無序R×C列聯(lián)表資料”，同時對其進(jìn)行CMHχ2檢驗?；凇氨碓u分”算法所需要的SAS程序如下：

以上是基于“表評分”算法所得CMH χ2檢驗的三種計算結(jié)果，分別對應(yīng)三種不同的備擇假設(shè)：第一種，“非零相關(guān)（即行、列兩有序變量之間具有非零的相關(guān)關(guān)系）”，=6.635，df=1，P=0.0100；第二種，“行評分均值不同（即各行上有序結(jié)果變量的評分值的平均數(shù)不等）”，=33.7083，d f=4，P＜0.0001；第三種，“常規(guī)關(guān)聯(lián)（即行、列兩無序變量之間存在關(guān)聯(lián)性）”，=36.6435，df=8，P＜0.0001。

以上三種檢驗都得出P＜0.01，故接受備擇假設(shè)。對本例而言，可以認(rèn)為：①“媒體暴露情況（時間）”與“ASD程度”之間存在“非零相關(guān)”，具體地說，隨著媒體暴露時間延長，出現(xiàn)重度ASD癥狀的比例呈非線性上升趨勢；②五種不同媒體暴露時間所對應(yīng)的平均ASD值是不相等的，除第2個時間段之外，其他4個時間段都隨著“暴露時間延長”，平均ASD值變大（即ASD癥狀程度逐漸加重）；③媒體暴露情況（時間）與ASD程度之間存在關(guān)聯(lián)性（前面的“①”和“②”可以清楚體現(xiàn)出兩變量之間是如何關(guān)聯(lián)的）。

4 給有序變量評分的方法

4.1 表評分方法

在前面的SAS程序中，在“tables語句”中出現(xiàn)了選項“SCORES=table”，其含義是用表格里規(guī)定的方式給有序變量評分。這里的“表格”，其實是指SAS程序中如何給“行變量”或“列變量”賦值。“DO A=1 TO 5；”就是讓“行變量 A”依次取 1、2、3、4、5分；同理，“DO B=1 TO 3；”就是讓“列變量B”依次取1、2、3分。當(dāng)然，也可以給變量A或B賦任意的幾個由小到大的數(shù)值，例如“DO A=0.1，1.2，13.5，21.8，34.6；”。

4.2 其他評分方法

在SAS/STAT的FREQ過程中，還有另外三種評分方法，分別為“rank評分法”“ridit評分法”和“修正的ridit評分法”，對應(yīng)的選項依次為“SCORES=rank”“SCORES=ridit”和“SCORES=modridit”。

4.3 評分方法與檢驗統(tǒng)計量計算方法之間的關(guān)系

對“常規(guī)關(guān)聯(lián)”的計算結(jié)果而言，四種評分方法所得計算結(jié)果是完全相同的，因為計算公式中不涉及“評分”；對“非零相關(guān)”的計算結(jié)果而言，“表評分法”（采用Pearson’s相關(guān)分析）與另外三種評分法（其結(jié)果是一樣的，采用Spearman’s秩相關(guān)分析）計算結(jié)果可能相差很大、甚至結(jié)論相反；對“行評分均值不同”的計算結(jié)果而言，“表評分法”（采用單因素R水平設(shè)計一元定量資料方差分析）與另外三種評分法（其結(jié)果是一樣的，采用單因素R水平設(shè)計一元定量資料Kruskal-Wallis’s H秩和檢驗）計算結(jié)果稍有差別。例如，例4“基于秩評分”的結(jié)果如下：

【說明】“ridit評分法”和“修正的ridit評分法”輸出結(jié)果與“基于秩評分”方法輸出結(jié)果相同，此處從略。

注意：“非零相關(guān)”的結(jié)果與前面“基于表評分”的結(jié)果相反；四種評分法對應(yīng)的“常規(guī)關(guān)聯(lián)”結(jié)果都是χ2=36.6435、P＜0.0001；四種評分法對應(yīng)的“行評分均值不同”有兩個不同的計算結(jié)果，分別為“χ2=33.7083［采用方差分析計算再基于式（6）轉(zhuǎn)換的結(jié)果］”與“χ2=31.2716（采用 Kruskal-Wallis’s H秩和檢驗計算的結(jié)果）”，但P值都小于0.0001。

5 討論與小結(jié)

5.1 討論

在給“scores=”選項指定“table”或其他三種非參數(shù)評分法（即rank、ridit和modridit）時，為什么“非零相關(guān)”有兩個完全相反的結(jié)果？原因在于：前者采用的是“Pearson’s相關(guān)分析”，而后者采用的是“Spearman’s秩相關(guān)分析”。就例3資料而言，分別基于“Pearson’s相關(guān)分析”與“Spearman’s秩相關(guān)分析”得到“r=0.02033”與“rs=0.00568”。依據(jù)式（2）可算出對應(yīng)的CMH χ2值分別為6.632368（與前面相應(yīng)的計算結(jié)果6.6350很接近）和0.517715（與前面相應(yīng)的計算結(jié)果0.5174很接近）。

在“常規(guī)關(guān)聯(lián)”的分析中，將有序變量視為無序變量或名義變量，故給“scores=”選項指定四種內(nèi)容（即table、rank、ridit和modridit）時，所得結(jié)果完全相同。就例3資料而言，。它在本質(zhì)上就是=36.6458（即Pearson’s擬合優(yōu)度χ2值），由式（7）可知：

在備擇假設(shè)為“行評分均值不同”時，SAS/STAT中FREQ過程將此時的資料視為“單因素R水平設(shè)計一元定量資料”。若采用“表評分（即scores=table）”，則SAS/STAT中FREQ過程采用“單因素方差分析”；若采用其他三種評分，SAS/STAT中FREQ過程采用Kruskal-Wallis’s秩和檢驗。事實上，就是把結(jié)果變量B視為“定量變量”，因此，對結(jié)果變量B的賦值方法不同，將直接影響R×C列聯(lián)表資料各行上結(jié)果變量平均值的計算結(jié)果。四種評分方法對應(yīng)結(jié)果變量B的3個評分值如下：

若直接采用SAS/STAT中ANOVA過程計算，得到F=8.44、df=4，基于式（6）可算得≈ 8.44 × 4=33.76，與“scores=table”對應(yīng)的“行評分均值不同”時的“=33.7083”比較接近；若直接采用SAS/STAT中NPAR1WAY過程計算，得到H=31.2716，與“scores=rank”對應(yīng)的“行評分均值不同”時的“=31.2716”完全相同。

值得一提的是，在分析二維列聯(lián)表資料時，三種非參數(shù)評分法輸出的三種檢驗結(jié)果相同；但若是分析高維列聯(lián)表資料，情況就不盡相同了。因篇幅所限，此處從略。

5.2 小結(jié)

本文呈現(xiàn)了三種R×C列聯(lián)表資料的實例及其三種相應(yīng)的統(tǒng)計分析方法，這些方法可概括成廣義CMH χ2檢驗；基于SAS/STAT中FREQ過程對R×C列聯(lián)表資料進(jìn)行分析時，可同時輸出三種統(tǒng)計分析結(jié)果；本文還詳細(xì)揭示了選取不同的“評分”方法，將會影響“非零相關(guān)”與“行評分均值不同”計算結(jié)果的真實原因。