統計模型總結 - GLM的理解

最近遇到一個問題，如果因變量為一個連續變量（如胰島素水平），主要考察的變量為分組變量（如正常血糖組，前糖尿病組，糖尿病組三組），現在的目的是想看調整多種變量（包括多個連續性變量和分類變量）后，胰島素水平是否一致。


一開始的思路想到的是采用GLM進行協方差分析來解決。

但是有覺得似乎不是很對，因為經典的協方差分析通常只考慮一個連續變量（協變量）和一個分組變量，同時協變量和分組變量只有不存在交互的時候（經典協方差分析的前提）才能使用協方差分析。
針對我目前的問題，如果想調整多個連續性變量和分類變量，這種方法能否再叫協方差分析？如果可以認為是協方差的思想，用不用檢驗協方差分析的前提（如協變量與分組變量之間的交互）？多個連續性變量和分類變量存在時，該前提應該怎么檢驗？

通過跟別人交流之后，有一句話非常受用：線性模型其實最重要的不在於用的方差分析還是協方差分析，而主要是檢驗殘差是否符合線性的幾個條件。

受到該啟發后，認真再復習GLM的相關資料，得到更加重要的總結如下（來自高惠璇SAS/STAT軟件使用手冊，實際是SAS8.2的User's guide的中文版，但是目前SAS 9.2，9.3的User's guide關於GLM模型的介紹中已經刪去了這么經典的總結，實在可惜，倒讓人看不到GLM的真正長處了）：
如果X1-X3，Y1-Y2為連續性變量，Y3為分類變量，a-c為分類變量，time為時間變量，目前我們熟悉的模型可以簡單概括如下：
(1) y1= x1    簡單回歸
(2) y1= x1 x2 x3  多重回歸（multiple regression）
(3) y1 y2=x1 x2   多元回歸（multivariate regression）
(4) y1= a   單因素方差分析
(5) y1= a b   （析因設計的）主效應分析
(6) y1= a b  a*b （析因設計的）主效應加交互項分析
(7) y1= a x1   協方差分析
(8) y3= a   單因素logistic回歸
(9) y3= a b c x1 x2 x3   多因素logistic回歸
(10) y3(time) =a   單因素cox回歸
(11) y3(time) = a b c x1 x2 x3  多因素cox回歸

1-7采用SAS的一般線性模型GLM都能實現，而1-9采用SAS的廣義線性模型GENMOD都能實現，具體驗證詳見后面舉例

再次回到開始的問題：掌握上述的基本思路后，因為因變量為連續變量，所以采用線性模型肯定是對的。如果因變量可以認為是正態的，那采用一般線性模型是合理的。所以現在的關鍵問題是：如果調整多個變量（包含分類和連續變量）后看不同分組間因變量（連續變量）是否仍有差異時，能否再稱為協方差分析？我目前認為應該是可以的，但是事實上我們遇到這種情況后，並不再去強調它是協方差的思想，而只是回到線性模型分析的最初的起點，也即是檢查殘差是否符合線性的基本條件即可。

但是現實中，我們在使用GLM解決前面遇到的類似問題時，只是簡單地用了，而很多時候我們都沒有認真去檢驗殘差是否符合這個條件，這可能是我們濫用GLM的表現之一，因為我們更多只關注模型的參數是否有意義，而不去關心對結果“無關緊要”的前提條件。

再次思考一個問題，上述列舉的1-7模型，在GLM中並沒有特定的選項指定是哪一種模型，而采用一種表達方式。由此，可以進一步深入概括一句話，GLM模型，對於上述列舉的1-7模型並沒有本質區別，唯一的區別只是模型中自變量的屬性和數量不同。但是我們對1-7模型的叫法卻不相同。而其原因是我們對事物的認識是一個由淺到深的過程，之前我們認為他們是不同的7件事情，隨着認識的加深，發現原來這些問題可以用一個方式表達出來。而SAS的GENMOD則更能說明這一問題。現在還沒有一個模型能把上述模型1-11用一種表達方式表示，但是COX回歸在拋開基線生存函數之后剩下的部分也是線性模型，所以說不定哪天真的能夠把上述所有模型用一種表達方式表示出來。到時候更應該相信人們對事物的認識絕對是一個由淺到深的過程啦。

data drugtest; 
    input Drug $ PreTreatment PostTreatment @@; 
    datalines; 
   A 11  6   A  8  0   A  5  2   A 14  8   A 19 11 
   A  6  4   A 10 13   A  6  1   A 11  8   A  3  0 
   D  6  0   D  6  2   D  7  3   D  8  1   D 18 18 
   D  8  4   D 19 14   D  8  9   D  5  1   D 15  9 
   F 16 13   F 13 10   F 11 18   F  9  5   F 21 23 
   F 16 12   F 12  5   F 12 16   F  7  1   F 12 20 
; 

proc glm  data=drugtest; 
    class Drug; 
    model PostTreatment = Drug PreTreatment / solution; 
    lsmeans Drug / stderr tdiff cov out=adjmeans; 
run; 

proc genmod data=drugtest; 
    class Drug; 
    model PostTreatment = Drug PreTreatment / dist=NOR link=ID obstats type1; 
    lsmeans Drug ; 
run; 

 

 參考：http://bbs.pinggu.org/thread-2736449-1-1.html  bymoonstone君