[R-sig-ME] glmmPQL: Temporal correlation structure for grouped data with observation gaps

Thierry Onkelinx th|erry@onke||nx @end|ng |rom |nbo@be
Wed Mar 10 12:21:15 CET 2021

Dear Jorge,

glmmPQL is not the only option if you need a correlation structure. glmmTMB
and INLA allow for correlated random effects. That is IMHO superior to
correlated residuals within the most detailed levels of the random effect
(what nlme does).

Furthermore you need the normalised residuals, as they take the correlation
structure into account.

Best regards,

ir. Thierry Onkelinx
Statisticus / Statistician

Vlaamse Overheid / Government of Flanders
Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
thierry.onkelinx using inbo.be
Havenlaan 88 bus 73, 1000 Brussel

To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to say
what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey


Op wo 10 mrt. 2021 om 05:58 schreef Jorge Garcia Molinos <garciamj using tcd.ie>:

> Dear list members,
> I have a question regarding the implementation of a temporal correlation
> structure in a GLMM context in R. This is my first time posting here so I
> hope the question is suitable for this forum.
> I have a data set that comprises summer (June-September) paired daily air
> and water temperatures measured at 83 monitoring sites over a period of 3
> years. The purpose is to predict water temperatures based on air
> temperatures and a suite of environmental variables related to the
> characteristics of the rivers and contributing catchment (land use, slope,
> elevation and so on). The number of observations among sites differ somehow
> (the series at some sites started at slightly different dates and some
> sites have missing values within the series).
> My approach so far has been as follows:
> 1. Fit a GLMM on water temperatures (WT) with air temperature (AT) and the
> different environmental variables as fixed factors and site as random
> factor (random intercept).
> model <- glmmPQL(WT ~ AT + H2OArea + avEleE + RipB60_fdec + RipB60_fev,
> random = ~1|Site, data = t_d, family = Gamma(link = "log"))
> The reason for using glmmPQL is because of the need to incorporate a
> temporal correlation structure (see point 2)
> 2. Given the nature of the data, temporal correlation is to be expected. A
> check of acf and pacf functions on the model residuals by site shows
> consistently a gradually decaying acf, and a sharp drop of the pacf below
> the significance threshold after lag 1 (for which the pacf equals approx..
> 0.8), suggesting an AR1 process is suitable. Considering the different
> observation gaps in the temperature series by station and the seasonal gaps
> between years (only summer data), I decided to try with a corCAR1 structure
> of the form:
> modelAR <- update(m5, correlation = corCAR1(0.8, form = ~days|Site))
> Where days correspond to the day number since the start of the series
> across all sites (the first day an observation was made across all series).
> That is:
> t_d$days <- difftime(t_d$date, min(t_d$date), units = "days")
> 3. However, the acf and pacf plots for the residuals in the updated model
> look nearly identical from those of the original model, which makes me
> think that I have not specified the correlation structure or the time
> variable correctly or that I’m missing something else?
> An extract of the data set for a few stations can be downloaded here:
> https://drive.google.com/file/d/1iMaXSSYvCJ1Xo09ZoJjbLdzmc0gjfJKy/view?usp=sharing
> Any help/suggestions will be greatly appreciated!
> Many thanks,
> Jorge
>         [[alternative HTML version deleted]]
> _______________________________________________
> R-sig-mixed-models using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

	[[alternative HTML version deleted]]

More information about the R-sig-mixed-models mailing list