2013-04-26 28 views
5

當您有一個具有許多因素和交互作用的多級模型時,固定效果矩陣的相關大小會變得相當大且不清楚。如何提取lmer輸出的固定效果部分的相關性

我可以使用symbolic.cor=T參數在打印方法,使摘要的更清晰的打印象下面這樣:

ratbrain <- 
within(read.delim("http://www-personal.umich.edu/~bwest/rat_brain.dat"), 
{ 
treatment <- factor(treatment, 
labels = c("Basal", "Carbachol")) 
region <- factor(region, 
labels = c("BST", "LS", "VDB")) 
}) 

print(mod<-lmer(activate ~ region * treatment + (0 + treatment | animal),ratbrain),symbolic.cor=T) 

此圖表用於大矩陣稍微更清晰的相關矩陣。儘管這個例子的矩陣並沒有那麼大。 但是如果我可以繪製相關性的熱圖,那將會很好。
如何提取固定效果的相關性,以便製作此熱圖?

編輯:

這裏是我創建感謝答案的功能。

fixeff.plotcorr<-function(mod,...) 
{ 
    #require(GGally) # contains another correlation plot using ggplot2 
    require(lme4) 

    fixNames<-names(fixef(mod)) 

    # Simon O'Hanlon's answer: 
    # so <- summary(mod) 
    # df<-as.matrix([email protected]@factors$correlation) for version lme4<1.0 
    # df<-as.matrix([email protected]$correlation) # lme4 >= 1.0 

    df<-as.matrix(cov2cor(vcov(mod))) #Ben Bolker's solution 

    rownames(df)<-fixNames 
    colnames(df)<-abbreviate(fixNames, minlength = 11) 

    colsc=c(rgb(241, 54, 23, maxColorValue=255), 'white', rgb(0, 61, 104, maxColorValue=255)) 
    colramp = colorRampPalette(colsc, space='Lab') 
    colors = colramp(100) 
    cols=colors[((df + 1)/2) * 100] 
    # I'm using function my.plotcorr which you can download here: 
    # http://hlplab.wordpress.com/2012/03/20/correlation-plot-matrices-using-the-ellipse-library/ 
    my.plotcorr(df, col=cols, diag='none', upper.panel="number", mar=c(0,0.1,0,0),...) 

    # Another possibility is the corrplot package: 
    # cols <- colorRampPalette(c("#67001F", "#B2182B", "#D6604D", "#F4A582", "#FDDBC7", 
    #        "#FFFFFF", "#D1E5F0", "#92C5DE", "#4393C3", "#2166AC", "#053061")) 
    # require(corrplot,quiet=T) 
    # corrplot(df, type="upper", method="number", tl.pos='tl', tl.col='black', tl.cex=0.8, cl.pos='n', col=cols(50)) 
    # corrplot(df,add=TRUE, method='ellipse', type='lower', tl.pos='n', tl.col='black', cl.pos='n', col=cols(50), diag=FALSE) 
} 

您必須從here下載my.plotcorr功能。 產生的上述使用命令fixeff.plotcorr(mod)現在的例子中的情節是這樣的: enter image description here

回答

4

使用S4方法列表功能,我們可以返回時print被稱爲"mer"類的一個對象在其上調度的功能:

selectMethod(print , "mer") 

望着那被返回,我們可以找到適用的線條給你的源代碼:

if (correlation) { 
      corF <- [email protected]@factors$correlation 

so被定義爲你的對象的總結,所以你的情況,你應該只需要提取:

so <- summary(mod) 
[email protected]@factors$correlation 
5

我不知道直接的方法。但這是解決方法。

diag(diag(1/sqrt(vcov(mod)))) %*% vcov(mod) %*% diag(diag(1/sqrt(vcov(mod)))) 
+1

+1我認爲這不是一種解決方法,而是從v-cov矩陣的相關性的數學推導! :-) – 2013-04-26 13:13:46

5

如何使用內置

cov2cor(vcov(mod)) 

+0

不錯的一個!我不知道這個功能。 – Robert 2014-04-04 12:27:22