我正在跟进这个很好的答案。功能foo下面，取Name的列VarCorr(fit)输出并使它们成为列名summary(rePCA(fit)) call.
当我们输入时效果很好fm1, fm2，但我想知道为什么它失败了fm3？有解决办法吗？
library(lme4)
dat <- read.csv('https://raw.githubusercontent.com/rnorouzian/e/master/sng.csv')
fm1 <- lmer(diameter ~ 1 + (1|plate) + (1|sample), Penicillin)
fm2 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
fm3 <- lmer(y ~ A * B * C + (A + B | group) + (C|group), data = dat)
foo <- function(fit) {
obj <- summary(rePCA(fit))
model <- VarCorr(fit)
Map(function(x, z) {
colnames(x$importance) <- paste(z, unique(sapply(model, colnames)), sep = '_')
x
}, obj, names(obj))
}
#EXAMPLE OF USE:
foo(fm1) ### OK !
foo(fm2) ### OK !
foo(fm3) ### :-( Error in dimnames(x) <- dn
当长度为obj and model是不同的。这是一个让它工作的技巧fm3.
foo <- function(fit) {
obj <- summary(rePCA(fit))
model <- VarCorr(fit)
if(length(obj) == length(model)) {
obj <- Map(function(x, z) {
colnames(x$importance) <- paste(z, unique(sapply(model, colnames)), sep = '_')
x
}, obj, names(obj))
}
else if(length(obj) == 1) {
colnames(obj[[1]]$importance) <- unlist(mapply(paste, names(model), sapply(model, colnames), MoreArgs = list(sep = '_')))
}
return(obj)
}
这将返回以下输出：
foo(fm1)
#$plate
#Importance of components:
# plate_(Intercept)
#Standard deviation 1.54
#Proportion of Variance 1.00
#Cumulative Proportion 1.00
#$sample
#Importance of components:
# sample_(Intercept)
#Standard deviation 3.51
#Proportion of Variance 1.00
#Cumulative Proportion 1.00
foo(fm2)
#$Subject
#Importance of components:
# Subject_(Intercept) Subject_Days
#Standard deviation 0.967 0.2309
#Proportion of Variance 0.946 0.0539
#Cumulative Proportion 0.946 1.0000
foo(fm3)
#$group
#Importance of components:
# group_(Intercept) group_A group_B group.1_(Intercept) group.1_C
#Standard deviation 1.418 1.291 0.5129 0.4542 0.0000497
#Proportion of Variance 0.485 0.402 0.0634 0.0498 0.0000000
#Cumulative Proportion 0.485 0.887 0.9502 1.0000 1.0000000
本文内容由网友自发贡献，版权归原作者所有，本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容，请联系:hwhale#tublm.com(使用前将#替换为@)
