2017-04-26 12 views
0

sjpPlot, sjp.int 함수를 사용하여 lme의 상호 작용을 계획합니다. 운영자 값의 옵션은 평균 +/- sd, 4 분위수, 모두, 최대/최소입니다. 평균 +/- 2sd를 플롯 할 수있는 방법이 있습니까?sjPlot 상호 작용에서 중재자 옵션 조정

일반적으로는 다음과 같이 될 것이다 :

재현 예를

model <- lme(outcome ~ var1+var2*time, random=~1|ID, data=mydata, na.action="na.omit") 
sjp.int(model, show.ci=T, mdrt.values="meansd") 

많은 감사 :

#create data 
mydata <- data.frame(SID=sample(1:150,400,replace=TRUE),age=sample(50:70,400,replace=TRUE), sex=sample(c("Male","Female"),200, replace=TRUE),time= seq(0.7, 6.2, length.out=400), Vol =rnorm(400),HCD =rnorm(400)) 
mydata$time <- as.numeric(mydata$time) 

#insert random NAs 
    NAins <- NAinsert <- function(df, prop = .1){ 
n <- nrow(df) 
m <- ncol(df) 
num.to.na <- ceiling(prop*n*m) 
id <- sample(0:(m*n-1), num.to.na, replace = FALSE) 
rows <- id %/% m + 1 
cols <- id %% m + 1 
sapply(seq(num.to.na), function(x){ 
    df[rows[x], cols[x]] <<- NA 
} 
) 
return(df) 
} 


mydata2 <- NAins(mydata,0.1) 

#run the lme which gives error message 
model = lme(Vol ~ age+sex*time+time* HCD, random=~time|SID,na.action="na.omit",data=mydata2);summary(model) 

mydf <- ggpredict(model, terms=c("time","HCD [-2.5, -0.5, 2.0]")) 

#lmer works 
model2 = lmer(Vol ~ age+sex*time+time* HCD+(time|SID),control=lmerControl(check.nobs.vs.nlev = "ignore",check.nobs.vs.rankZ = "ignore", check.nobs.vs.nRE="ignore"), na.action="na.omit",data=mydata2);summary(model) 
mydf <- ggpredict(model2, terms=c("time","HCD [-2.5, -0.5, 2.0]")) 

#plotting gives problems (jittered lines) 
plot(mydf) 

enter image description here

답변

1

sjPlot , 그것은 현재 수 없습니다. 그러나 나는 한계 효과를 계산하고 플로팅하기위한 전용 패키지를 작성했습니다 : ggeffects. 이 패키지는 약간의 융통성이 있습니다 (가장자리 효과 플롯의 경우).

ggeffects -package에는 특정 값에서 한계 효과를 계산할 수있는 ggpredict() -function이 있습니다. 당신은 문제의 모델 기간의 SD를 알게되면, 당신은 당신의 상호 작용을 플롯 함수 호출이 값을 지정할 수 있습니다

package-vignette에 몇 가지 예제가 있습니다
library(ggeffects) 
# plot interaction for time and var2, for values 
# 10, 30 and 50 of var2 
mydf <- ggpredict(model, terms = c("time", "var2 [10,30,50]")) 
plot(mydf) 

, 특히 this section를 참조하십시오. 여기에 편집

은 재현 예에 따라 결과이다 (GitHub의-버전이 현재 필요하다고 참고!) :

# requires at least the GitHub-Versiob 0.1.0.9000! 
library(ggeffects) 
library(nlme) 
library(lme4) 
library(glmmTMB) 

#create data 
mydata <- 
    data.frame(
    SID = sample(1:150, 400, replace = TRUE), 
    age = sample(50:70, 400, replace = TRUE), 
    sex = sample(c("Male", "Female"), 200, replace = TRUE), 
    time = seq(0.7, 6.2, length.out = 400), 
    Vol = rnorm(400), 
    HCD = rnorm(400) 
) 
mydata$time <- as.numeric(mydata$time) 

#insert random NAs 
NAins <- NAinsert <- function(df, prop = .1) { 
    n <- nrow(df) 
    m <- ncol(df) 
    num.to.na <- ceiling(prop * n * m) 
    id <- sample(0:(m * n - 1), num.to.na, replace = FALSE) 
    rows <- id %/% m + 1 
    cols <- id %% m + 1 
    sapply(seq(num.to.na), function(x) { 
    df[rows[x], cols[x]] <<- NA 
    }) 
    return(df) 
} 

mydata2 <- NAins(mydata, 0.1) 

# run the lme, works now 
model = lme(
    Vol ~ age + sex * time + time * HCD, 
    random = ~ time | 
    SID, 
    na.action = "na.omit", 
    data = mydata2 
) 
summary(model) 

mydf <- ggpredict(model, terms = c("time", "HCD [-2.5, -0.5, 2.0]")) 
plot(mydf) 

LME 플롯

enter image description here

# lmer also works 
model2 <- lmer(
    Vol ~ age + sex * time + time * HCD + (time | 
              SID), 
    control = lmerControl(
    check.nobs.vs.nlev = "ignore", 
    check.nobs.vs.rankZ = "ignore", 
    check.nobs.vs.nRE = "ignore" 
), 
    na.action = "na.omit", 
    data = mydata2 
) 
summary(model) 
mydf <- ggpredict(model2, terms = c("time", "HCD [-2.5, -0.5, 2.0]"), ci.lvl = NA) 

# plotting works, but only w/o CI 
plot(mydf) 

lmer 플롯

enter image description here

# lmer also works 
model3 <- glmmTMB(
    Vol ~ age + sex * time + time * HCD + (time | SID), 
    data = mydata2 
) 
summary(model) 
mydf <- ggpredict(model3, terms = c("time", "HCD [-2.5, -0.5, 2.0]")) 
plot(mydf) 
plot(mydf, facets = T) 

glmmTMB-플롯

enter image description here

enter image description here

,
+0

안녕 다니엘, 고마워! 나는 ggeffects로 놀고있다. lme 함수를 사용하여 플롯 할 수 없었지만 lmer가 작동하는 것 같습니다. lme의 오류 메시지 : "predict.lme의 오류 (모델, newdata = fitfram, type ="응답 ", ...) : 'newdata'에서 원하는 수준의 그룹을 평가할 수 없습니다." 그게 lme로 끝났어?(random intercept/random slope) 그러나 mydf의 출력을 plot 또는 ggplot으로 플로팅 할 때 기울기는 직선이 아니며 지터가 있습니다. 왜 그랬을까요? – user6121484

+0

재현 가능한 예가 있습니까? 원하는 경우 전자 메일로 보낼 수 있습니다. – Daniel

+0

현재 [GitHub-version of ggeffects] (https://github.com/strengejacke/ggeffects)는 lme-models로 문제를 해결해야합니다. – Daniel