จะพล็อตจุดตัดและความชันแบบสุ่มในแบบจำลองผสมกับตัวทำนายหลายตัวได้อย่างไร

เป็นไปได้ไหมที่จะพล็อตจุดตัดสุ่มหรือความชันของแบบจำลองผสมเมื่อมีตัวทำนายมากกว่าหนึ่งตัว

ด้วยตัวทำนายหนึ่งตัวฉันจะทำสิ่งนี้:

#generate one response, two predictors and one factor (random effect)
resp<-runif(100,1, 100)
pred1<-c(resp[1:50]+rnorm(50, -10, 10),resp[1:50]+rnorm(50, 20, 5))
pred2<-resp+rnorm(100, -10, 10)
RF1<-gl(2, 50)

#gamm
library(mgcv)
mod<-gamm(resp ~ pred1, random=list(RF1=~1))
plot(pred1, resp, type="n")
for (i in ranef(mod$lme)[[1]]) {
abline(fixef(mod$lme)[1]+i, fixef(mod$lme)[2])
}

#lmer
library(lme4)
mod<-lmer(resp ~ pred1 + (1|RF1))
plot(pred1, resp, type="n")
for (i in ranef(mod)[[1]][,1]) {
abline(fixef(mod)[1]+i, fixef(mod)[2])
}

แต่ถ้าฉันมีโมเดลแบบนี้แทนล่ะ:

mod<-gamm(resp ~ pred1 + pred2, random=list(RF1=~1))

หรือกับแอลเมอร์

mod<-lmer(resp ~ pred1 + pred2 + (1|RF1))

ฉันควรพิจารณาสัมประสิทธิ์ทั้งหมดหรือเฉพาะตัวแปรที่ฉันกำลังพล็อตเท่านั้น

ขอบคุณ


r plot mixed-models random-effects lmer
person Oritteropus    schedule 14.07.2013    source แหล่งที่มา
comment
โดยพื้นฐานแล้ว คุณต้องตัดสินใจว่าคุณต้องการทำอะไรกับตัวแปรอื่นๆ ขั้นตอนที่พบได้บ่อยที่สุดคือการเลือกค่าอ้างอิงสำหรับตัวแปรหนึ่งตัว (เช่น pred2 เท่ากับค่าเฉลี่ย) และวาดจุดความชันเทียบกับ pred1 สำหรับค่านั้น หรือคุณสามารถเลือกค่า pred2 หลายค่าและพล็อตบรรทัด (ชุด) สำหรับแต่ละค่า อาจเป็นในแผนย่อยแยกกัน หรือ (น่าเกลียดที่สุด) ทำพล็อต 3 มิติและพล็อตระนาบ resp~f(pred1,pred2) แทน   -  person Ben Bolker    schedule 15.07.2013
comment
ขอบคุณเบน ขออภัย แต่ฉันไม่แน่ใจว่าจะติดตามคุณ คุณหมายความว่าอย่างไรสำหรับการเลือกค่าอ้างอิงสำหรับตัวแปรตัวเดียว คุณจะปฏิบัติอย่างไร?   -  person Oritteropus    schedule 15.07.2013

คำตอบ (1)


arrow_upward
5
arrow_downward 

## generate one response, two predictors and one factor (random effect)
set.seed(101)
resp <- runif(100,1,100)
pred1<- rnorm(100, 
           mean=rep(resp[1:50],2)+rep(c(-10,20),each=50),
           sd=rep(c(10,5),each=50))
pred2<- rnorm(100, resp-10, 10)

หมายเหตุ ว่าคุณควรจะ ไม่ กำลังพยายามปรับเอฟเฟกต์สุ่มให้พอดีสำหรับตัวแปรการจัดกลุ่มที่มีเพียงสองระดับ ซึ่งเกือบจะคงที่ส่งผลให้เกิดความแปรปรวนของเอฟเฟกต์สุ่มโดยประมาณของ ศูนย์ ซึ่งจะทำให้เส้นที่คุณคาดการณ์ไว้ทับกัน -- ฉันกำลังเปลี่ยนจาก gl(2,50) เป็น gl(10,10) ...

RF1<-gl(10,10)
d <- data.frame(resp,pred1,pred2,RF1)

#lmer
library(lme4)
mod <- lmer(resp ~ pred1 + pred2 + (1|RF1),data=d)

เวอร์ชันพัฒนาของ lme4 มีฟังก์ชัน predict() ที่ทำให้ง่ายขึ้นเล็กน้อย ...

  • คาดการณ์ช่วง pred1 โดยที่ pred2 เท่ากับค่าเฉลี่ย และในทางกลับกัน ทั้งหมดนี้ฉลาดกว่าที่ควรจะเป็นเล็กน้อย เนื่องจากมันสร้างค่าทั้งหมดสำหรับตัวทำนายโฟกัสและพล็อตพวกมันด้วย ggplot ในคราวเดียว ...

()

nd <- with(d,
           rbind(data.frame(expand.grid(RF1=levels(RF1),
                      pred1=seq(min(pred1),max(pred1),length=51)),
                      pred2=mean(pred2),focus="pred1"),
                 data.frame(expand.grid(RF1=levels(RF1),
                      pred2=seq(min(pred2),max(pred2),length=51)),
                      pred1=mean(pred1),focus="pred2")))
nd$val <- with(nd,pred1[focus=="pred1"],pred2[focus=="pred2"])
pframe <- data.frame(nd,resp=predict(mod,newdata=nd))
library(ggplot2)
ggplot(pframe,aes(x=val,y=resp,colour=RF1))+geom_line()+
         facet_wrap(~focus,scale="free")
  • อีกทางหนึ่ง มุ่งเน้นไปที่ pred1 และสร้างการคาดการณ์สำหรับช่วง (เล็ก/ไม่ต่อเนื่อง) ของค่า pred2 ...

()

nd <- with(d,
           data.frame(expand.grid(RF1=levels(RF1),
                      pred1=seq(min(pred1),max(pred1),length=51),
                      pred2=seq(-20,100,by=40))))
pframe <- data.frame(nd,resp=predict(mod,newdata=nd))
ggplot(pframe,aes(x=pred1,y=resp,colour=RF1))+geom_line()+
         facet_wrap(~pred2,nrow=1)

คุณอาจต้องการตั้งค่า scale="free" ใน facet_wrap() สุดท้าย ... หรือใช้ facet_grid(~pred2,labeller=label_both)

สำหรับการนำเสนอ คุณอาจต้องการแทนที่สุนทรียศาสตร์ colour ด้วย group หากคุณต้องการแยกความแตกต่างระหว่างกลุ่มต่างๆ (เช่น ลงจุดบรรทัดแยกกัน) แทนที่จะระบุ ...

person Ben Bolker    schedule 15.07.2013
comment
มีประโยชน์มาก! ขอบคุณ - person Oritteropus; 15.07.2013