Построение графической модели приводит к решетке
Я подгоняю модель используя gam
от mgcv
пакет. Я сохраняю результат в model
и до сих пор я смотрел на гладкие компоненты, используя plot(model)
, Я недавно начал использовать lattice
и нравится его вывод. Поэтому мне интересно, возможно ли построить эти графики с помощью lattice
?
Это мой набор данных: https://gist.github.com/plxsas/fcef4a228c18c772b4f3
m2<- gam(TotalInd ~ s(dayinyear, by=as.numeric(Site=="1"), bs="cr")
+s(dayinyear, by=as.numeric(Site=="2"), bs="cr") + s(dayinyear,
by=as.numeric(Site=="3"), bs="cr"), random=list(Replicate=~ 1), data=data)
Как я могу сделать сюжет этой модели в lattice
Пакет с тремя панелями, представляющими мои три сайта более гладкими, пожалуйста?
Вы также могли заметить, что я использовал dayinyear
вместо правильного формата месяца (первый столбец в данных). Это связано с тем, что обобщенные аддитивные модели не обрабатывают категориальные переменные. Тем не менее, я хотел бы представить время в моем графике с названиями месяцев (как в первом столбце). Кто-нибудь знает, как двигаться вперед в этом lattice
сюжет?
1 ответ
Вот общий способ сделать это, используя поддельные данные. Вам нужно настроить это, чтобы убедиться, что имена такие, как вам нравится,
library(reshape)
library(mgcv)
library(lattice)
X1<-rnorm(100) # Make some fake data
X2<-rnorm(100)
X3<-rnorm(100)
Y<-rnorm(100)
Mod<-gam(Y~s(X1,bs="cr")+s(X2,bs="cr")+s(X3, bs="cr")) # make a model
Z<- predict(Mod,type="terms", se.fit=T) #Z is the predicted value
#for each smooth term, se.fit give you SE
Z2<-melt(Z$fit) #Z was in wide form, Z2 is long form
Z2$XX<-c(X1,X2,X3) #add the original values for he predictors
Z2$SE<-melt(Z$se.fit)$value #add SE
Z2$UP<-Z2$value+2*Z2$SE #+2 SE
Z2$Low<-Z2$value-2*Z2$SE # - 2 SE
Z2<-Z2[order(Z2$XX),]
xyplot(value~XX|X2,data=Z2,type="l",col="black",as.table=T,
prepanel=function (x,y,...)list(ylim=c(min(Z2$Low),max(Z2$UP))),
panel=function(x,y,groups,subscripts,...){
panel.xyplot(x,y,...)
panel.lines(Z2$UP[subscripts]~Z2$XX[subscripts],lty=2, col="red")
panel.lines(Z2$Low[subscripts]~Z2$XX[subscripts],lty=2, col="red")
}
)
value
является прогнозируемым значением для каждого предиктора и X2
где переменная группировки (указывает, какие данные принадлежат каждому предиктору). Если вы много работаете с нами, вы должны переименовать вещи, чтобы прояснить ситуацию. order
часть просто избегает спагетти
Вы можете контролировать способ обозначения оси X, используя at
а также labels
аргументы для оси х в scales
аргумент. Подробнее см. ?xyplot
Обновление - вот версия, которая работает с этими данными
m2<- gam(TotalInd ~ s(dayinyear, by=as.numeric(Site=="1"), bs="cr")
+s(dayinyear, by=as.numeric(Site=="2"), bs="cr")
+ s(dayinyear, by=as.numeric(Site=="3"), bs="cr"),
random=list(Replicate=~ 1), data=Data)
Z<- predict(m2,type="terms",se.fit=T) #Z is the predicted value and SE
Z2<-melt(Z$fit) #Z was in wide form, Z2 is long form
Z2$dayinyear<-Data$dayinyear #add the original values for he predictors
Z2$SE<-melt(Z$se.fit)$value
Z2$UP<-Z2$value+2*Z2$SE
Z2$Low<-Z2$value-2*Z2$SE
Z2<-Z2[Z2$value!=0,] #gets rid of excess zeroes
Z2<-Z2[order(Z2$dayinyear),]
xyplot(value~dayinyear|X2,data=Z2,type="l",col="black",as.table=T,
prepanel=function (x,y,...)list(ylim=c(min(Z2$Low),max(Z2$UP))),
panel=function(x,y,groups,subscripts,...){
panel.xyplot(x,y,...)
panel.lines(Z2$UP[subscripts]~Z2$dayinyear[subscripts],lty=2, col="red")
panel.lines(Z2$Low[subscripts]~Z2$dayinyear[subscripts],lty=2, col="red")
}
)
Обратите внимание, что я изменил название стартового data.frame
от data
в Data
РЕДАКТИРОВАТЬ - я добавил две пунктирные линии, которые показывают + /- 2 SE для каждого графика