Dear R community,
I've got several sites with 2 sensors measuring sapflow on tree level. Each
sensors contains missing data. To fill these gaps, I'm trying to fit a
linear (mixed) model of 1 sensor in dependence of the other one across all
sites.
I tried to do this both with lme (nlme) and gamm (mgcv), but somehow I seem
to miss the correct model definition. There are strange patterns within the
residuals which I don't know how to get rid of (see
*https://webdisk.ads.mwn.de/Handlers/AnonymousDownload.ashx?folder=1a7cbaa4&path=xylemRohWeekXnn2011.jpg
<https://webdisk.ads.mwn.de/Handlers/AnonymousDownload.ashx?folder=1a7cbaa4&path=xylemRohWeekXnn2011.jpg>*
).
I'm not quite sure, but it seems that I'm not quire getting rid of the
autocorrelation, though playing with corARMA and different p & q values
between 0 and 3 didn't improve it...
I'm quite clueless how to continue and I'm grateful for any hint and help!
Thanks,
Katharina
<snip>
library(RCurl)
library(mgcv)
library(nlme)
#retrieve xylemRohWeekXnn2011 test data frame
eval( expr = parse( text =
getURL("
https://webdisk.ads.mwn.de/Handlers/AnonymousDownload.ashx?folder=1a7cbaa4&path=xylemRohWeekXnn2011.R
")
))
#plot data
ggplot() + geom_point(aes(x=RecordTimestamp, y=sensor2, color=site),
size=0.4, data=xylemRohWeekXnn2011)
ggplot() + geom_point(aes(x=RecordTimestamp, y=sensor1, color=site),
size=0.4, data=xylemRohWeekXnn2011)
ggplot() + geom_point(aes(x=sensor2, y=sensor1, color=site), size=0.4,
data=xylemRohWeekXnn2011)
#model lme
xylemRohWeekXnn.fit.lme <- lme(sensor1 ~ sensor2 ,random=~sensor2|site,
data=xylemRohWeekXnn2011,
correlation=corAR1(),na.action=na.omit,control=lmeControl(opt = "optim"))
plot(residuals(xylemRohWeekXnn.fit.lme) ~ fitted(xylemRohWeekXnn.fit.lme),
col=xylemRohWeekXnn2011$site, main="xylemRohWeekXnn.fit.lme")
#model gamm
xylemRohWeekXnn2011$dum <- rep(1,times=nrow(xylemRohWeekXnn2011))
xylemRohWeekXnn.fit.gamm <- gamm(sensor1 ~ s(sensor2) +s(site,
bs="re",by=dum) + s(site, sensor2, bs="re",by=dum),
data=xylemRohWeekXnn2011, correlation=corAR1(), na.action=na.omit)
gam.check( xylemRohWeekXnn.fit.gamm$gam)
</snap>
[[alternative HTML version deleted]]
_______________________________________________
R-sig-ecology mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
