r語言 似然比檢驗置信區間_r語言 模型似然比檢驗 - CSDN

介紹

首先，請注意，圍繞多級模型的術語非常不一致。例如，多級模型本身可以稱為分級線性模型，隨機效應模型，多級模型，隨機截距模型，隨機斜率模型或匯集模型。根據學科，使用的軟體和學術文獻，許多這些術語可能指的是相同的一般建模策略。 

 

讀入數據

多級模型適用於特定類型的數據結構，其中單元嵌套在組內（通常為5個以上組），並且我們希望對數據的組結構進行建模。

## id extro open agree social class school## 1 1 63.69 43.43 38.03 75.06 d IV## 2 2 69.48 46.87 31.49 98.13 a VI## 3 3 79.74 32.27 40.21 116.34 d VI## 4 4 62.97 44.41 30.51 90.47 c IV## 5 5 64.25 36.86 37.44 98.52 d IV## 6 6 50.97 46.26 38.83 75.22 d I

在這裡，我們有關於嵌套在課堂內和學校內的科目外向性的數據。

在開始之前讓我們先了解一下數據的結構：

## 'data.frame': 1200 obs. of 7 variables:## $ id : int 1 2 3 4 5 6 7 8 9 10 ...## $ extro : num 63.7 69.5 79.7 63 64.2 ...## $ open : num 43.4 46.9 32.3 44.4 36.9 ...## $ agree : num 38 31.5 40.2 30.5 37.4 ...## $ social: num 75.1 98.1 116.3 90.5 98.5 ...## $ class : Factor w/ 4 levels "a","b","c","d": 4 1 4 3 4 4 4 4 1 2 ...## $ school: Factor w/ 6 levels "I","II","III",..: 4 6 6 4 4 1 3 4 3 1 ...

在這裡，我們看到我們有兩個可能的分組變量 - class和school。讓我們進一步探索它們：

## ## a b c d ## 300 300 300 300

## ## I II III IV V VI ## 200 200 200 200 200 200

## ## I II III IV V VI## a 50 50 50 50 50 50## b 50 50 50 50 50 50## c 50 50 50 50 50 50## d 50 50 50 50 50 50

這是一個完美平衡的數據集。您很可能沒有使用完美平衡的數據集，但我們將在未來探討其中的含義。現在，讓我們稍微繪製一下數據。使用包中的優秀xyplot功能lattice，我們可以跨變量探索學校和班級之間的關係。

 ​

 

在這裡，我們看到，班內有明顯的分層，我們也看到，social變量是從強不同open和agree變量。我們還看到，班級a和班級d在最低和最高頻段分別有更多的傳播。讓我們接下來繪製數據school。​

 

通過學校我們看到學生緊密分組，但學校I和學校的VI分散程度遠遠高於其他學校。我們的預測因子中的相同模式在學校之間就像在課堂之間一樣。​

 

在這裡我們可以看到，學校和階級似乎在密切區分我們的預測者和外向性之間的關係。

探索merMod對象的內部

在上一個教程中，我們為嵌套數據擬合了一系列隨機攔截模型。我們lmerMod將更深入地研究在擬合此模型時生成的對象，以便了解如何使用R中的混合效果模型。我們首先擬合下面按類分組的基本示例：

## [1] "lmerMod"## attr(,"package")## [1] "lme4"

首先，我們看到它MLexamp1現在是該類的R對象lmerMod。這個lmerMod對象是一個S4類，為了探索它的結構，我們使用slotNames：

## [1] "resp" "Gp" "call" "frame" "flist" "cnms" "lower" ## [8] "theta" "beta" "u" "devcomp" "pp" "optinfo"

在lmerMod對象中，我們看到了許多我們可能希望探索的對象。要查看其中的任何一個，我們只需鍵入MLexamp1@然後鍵入插槽名稱本身即可。例如：

## lmer(formula = extro ~ open + agree + social + (1 | school), ## data = lmm.data)

## [1] 59.116514 0.009751 0.027788 -0.002151

## [1] "data.frame"

## extro open agree social school## 1 63.69 43.43 38.03 75.06 IV## 2 69.48 46.87 31.49 98.13 VI## 3 79.74 32.27 40.21 116.34 VI## 4 62.97 44.41 30.51 90.47 IV## 5 64.25 36.86 37.44 98.52 IV## 6 50.97 46.26 38.83 75.22 I

該merMod對象有許多可用的方法 - 在這裡枚舉太多。但是，我們將在下面的列表中介紹一些更常見的內容：

## [1] anova.merMod* as.function.merMod* coef.merMod* ## [4] confint.merMod deviance.merMod* df.residual.merMod* ## [7] drop1.merMod* extractAIC.merMod* extractDIC.merMod* ## [10] family.merMod* fitted.merMod* fixef.merMod* ## [13] formula.merMod* isGLMM.merMod* isLMM.merMod* ## [16] isNLMM.merMod* isREML.merMod* logLik.merMod* ## [19] model.frame.merMod* model.matrix.merMod* ngrps.merMod* ## [22] nobs.merMod* plot.merMod* predict.merMod* ## [25] print.merMod* profile.merMod* qqmath.merMod* ## [28] ranef.merMod* refit.merMod* refitML.merMod* ## [31] residuals.merMod* sigma.merMod* simulate.merMod* ## [34] summary.merMod* terms.merMod* update.merMod* ## [37] VarCorr.merMod* vcov.merMod weights.merMod* ## ## Non-visible functions are asterisked

常見的需求是從merMod對象中提取固定效果。fixef提取固定效果的命名數字向量，這很方便。

## (Intercept) open agree social ## 59.116514 0.009751 0.027788 -0.002151

如果您想了解這些參數的p值或統計顯著性，請首先查看lme4幫助，?mcmcsamp了解各種方法的執行情況。內置於包中的一種便捷方式是：

## 0.5 % 99.5 %## .sig01 4.91840 23.88758## .sigma 2.53287 2.81456## (Intercept) 46.27751 71.95610## open -0.02465 0.04415## agree -0.01164 0.06721## social -0.01493 0.01063

從這裡我們可以首先看到我們的固定效果參數重疊0表示沒有效果的證據。我們還可以看到.sig01，這是我們對隨機效應變化的估計，是非常大且非常廣泛的定義。這表明我們的團隊之間可能缺乏精確性 - 要麼是因為群體之間的群體影響很小，要麼得到更精確的估計的群體太少，我們每個群體中的單位太少，或者所有群體的組合都是以上。

另一個常見的需求是提取殘餘標準誤差，這是計算效果大小所必需的。要獲得殘差標準誤的命名向量：

## [1] 2.671

例如，教育研究中的常見做法是通過將固定效應參數除以殘差標準誤差來將固定效果標準化為「效果大小」，這可以lme4很容易地實現：

## (Intercept) open agree social ## 22.1336705 0.0036508 0.0104042 -0.0008055

從這一點，我們可以看到，我們對開放性，宜人性和社交性的預測因素在預測外向性方面幾乎毫無用處 - 正如我們的情節所顯示的那樣。讓我們把注意力轉向下一個隨機效應。

探索組變化和隨機效果

您很可能適合混合效果模型，因為您直接對模型中的組級變化感興趣。目前還不清楚如何從結果中探索這種群體水平的變化summary.merMod。我們從這個輸出得到的是組效應的方差和標準偏差，但我們不會對個別組產生影響。這是ranef功能派上用場的地方。

## $school## (Intercept)## I -14.091## II -6.183## III -1.971## IV 1.966## V 6.331## VI 13.948

運行該ranef功能為我們提供了每個學校的攔截，但沒有太多額外的信息 - 例如這些估計的精確度。為此，我們需要一些額外的命令：

## [1] "ranef.mer"

## , , 1## ## [,1]## [1,] 0.03565## ## , , 2## ## [,1]## [1,] 0.03565## ## , , 3## ## [,1]## [1,] 0.03565## ## , , 4## ## [,1]## [1,] 0.03565## ## , , 5## ## [,1]## [1,] 0.03565## ## , , 6## ## [,1]## [1,] 0.03565

該ranef.mer對象是一個列表，其中包含每個組級別的data.frame。數據框包含每個組的隨機效果（這裡我們只對每個學校進行攔截）。當我們要求lme4隨機效應的條件方差時，它被存儲在attribute那些數據幀的一個中作為方差 - 協方差矩陣的列表。

這種結構確實很複雜，但它很強大，因為它允許嵌套，分組和跨級隨機效果。此外，創建者lme4已經為用戶提供了一些簡單的快捷方式，以便從ranef.mer對象中獲得他們真正感興趣的內容。

## $school## (Intercept)## I -14.091## II -6.183## III -1.971## IV 1.966## V 6.331## VI 13.948## ## with conditional variances for "school"

## $school

 

此圖顯示了dotplot隨機效應。

使用模擬和圖來探索隨機效應

一種常見的計量經濟學方法是創建所謂的集團級術語的經驗貝葉斯估計。不幸的是，關於什麼構成隨機效應項的適當標準誤差甚至如何一致地定義經驗貝葉斯估計，沒有太多的一致意見。

## X1 X2 mean median sd## 1 I (Intercept) -14.053 -14.093 3.990## 2 II (Intercept) -6.141 -6.122 3.988## 3 III (Intercept) -1.934 -1.987 3.981## 4 IV (Intercept) 2.016 2.051 3.993## 5 V (Intercept) 6.378 6.326 3.981## 6 VI (Intercept) 13.992 14.011 3.987

該REsim函數為每個學校返回級別名稱X1，估計名稱，X2估計值的平均值，中位數和估計的標準差。

另一個便利功能可以幫助我們繪製這些結果，看看他們如何與以下結果進行比較dotplot：

 

這提供了對隨機效應分量之間的變化的更保守的觀點。根據您的數據收集方式和研究問題，可以採用其他方法來估算這些影響大小。但是，請謹慎行事。

作者推薦的另一種方法lme4涉及RLRsim包。使用該軟體包，我們可以測試隨機效應的包含是否改善了模型擬合，我們可以使用基於模擬的似然比檢驗來評估其他隨機效應項的p值。

## ## simulated finite sample distribution of LRT. (p-value based on## 10000 simulated values)## ## data: ## LRT = 2958, p-value < 2.2e-16

這裡exactLRT發出警告，因為我們最初使用REML而不是完全最大可能性來擬合模型。幸運的是，該refitML功能lme4允許我們使用完全最大可能性輕鬆地重新調整我們的模型，以便輕鬆地進行精確測試。

## ## simulated finite sample distribution of LRT. (p-value based on## 10000 simulated values)## ## data: ## LRT = 2958, p-value < 2.2e-16

在這裡，我們可以看到包含我們的分組變量是重要的，即使每個單獨的組的影響可能實質上很小和/或不精確地測量。這對於理解模型的正確規範很重要。我們的下一個教程將更詳細地介紹這樣的規範測試。

隨機效應至關重要？

如何解釋我們隨機效應的實質性影響？這在嘗試使用多級結構來理解分組可能對個體觀察產生的影響的觀察工作中通常是至關重要的。為此，我們選擇了12個隨機病例，然後我們模擬了extro它們在6所學校中的每一所學校的預測值。注意，這是一個非常簡單的模擬，僅使用固定效應的平均值和隨機效應的條件模式，而不是複製或採樣以獲得可變性的感覺。這將留給讀者和/或未來的教程練習！

現在我們已經建立了一個模擬數據幀，​​

 

這個圖表告訴我們，在每個情節中，代表一個案例，學校有很大的變化。因此，將每個學生轉移到另一所學校對外向學分數有很大影響。但是，每個學校的每個案例都有所不同嗎？

 

在這裡我們可以清楚地看到，在每個學校中，案例相對相同，表明群體效應大於個體效應。

這些圖可用於以實質方式證明群體和個體效果的相對重要性。可以做更多的事情來使圖表更具信息性，例如放置對結果的總可變性的參考，並且還觀察距離，移動組將每個觀察值從其真實值移開。

結論

lme4提供了一個非常強大的面向對象的工具集，用於處理R中的混合效果模型。理解lme4對象的模型擬合和置信區間需要一些勤奮的研究和使用各種函數和擴展lme4本身。在下一個教程中，我們將探索如何lme4為難以指定的模型確定隨機效應模型的適當規範和框架的貝葉斯擴展。我們還將探討廣義線性模型框架和glmer多級廣義線性建模的功能。

還有問題嗎？聯繫我們！

 

大數據部落 -中國專業的第三方數據服務提供商，提供定製化的一站式數據挖掘和統計分析諮詢服務

統計分析和數據挖掘諮詢服務：y0.cn/teradat（諮詢服務請聯繫官網客服）

​QQ：3025393450

​

【服務場景】  

科研項目; 公司項目外包;線上線下一對一培訓;數據採集;學術研究;報告撰寫;市場調查。

【大數據部落】提供定製化的一站式數據挖掘和統計分析諮詢服務

​