[R-sig-ME] glmmPQL crashes on inclusion of corSpatial object

Patrick Johann Schratz patrick.schratz at gmail.com
Mon Jul 25 15:38:08 CEST 2016


Dear Thierry,

thanks for your answer.

I also recognized that the spatial correlation structure is causing the 
problems as the random effects seem to work without any problems when 
omitting the spatial correlations part.
However, trying the formula with a random effect of 2 levels 
(=evaluation) instead of 4(=date) works:

/fit <- glmmPQL(fo, random = ~1 | evaluation, data = d, //
//                correlation = correl, family = binomial)/

Which gives me troubles tracking down the problem. glmmPQL should be 
able to handle 4 level random effects in combination with spatial 
correlation structures?

Setting "fixed = TRUE" in the corSpatial setup still causes the fitting 
of glmmPQL() to break in my case (with "date" as random effect). Did you 
also change other values of the corSpatial setup?

Thanks for your link to the R-INLA package. Looks very promising. 
However it will take me quite some time to adapt my code to R-INLA´s 
syntax and I doubt I have the time right now.

Best regards,
Patrick



Am 25.07.16 um 13:36 schrieb Thierry Onkelinx:
> Dear Patrick,
>
> It seems like the correlation structure makes the model unstable. I 
> get convergence when setting fixed = TRUE, but the estimate are very 
> unstable.
>
> I'd suggest to have a look at the INLA package which allows to fit 
> spatially correlated random intercepts. It's described in "Spatial and 
> Spatio-Temporal Bayesian models with R-INLA" (Blangiardo & Cameletti, 
> 2015)
>
> Best regards,
>
>
> ir. Thierry Onkelinx
> Instituut voor natuur- en bosonderzoek / Research Institute for Nature 
> and Forest
> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
> Kliniekstraat 25
> 1070 Anderlecht
> Belgium
>
> 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
>
> 2016-07-25 12:38 GMT+02:00 Patrick Johann Schratz 
> <patrick.schratz at gmail.com <mailto:patrick.schratz at gmail.com>>:
>
>     MacbookPro:
>
>     > sessionInfo()
>     R version 3.3.1 (2016-06-21)
>     Platform: x86_64-apple-darwin15.5.0 (64-bit)
>     Running under: OS X 10.11.6 (El Capitan)
>
>     locale:
>     [1] de_DE.UTF-8/de_DE.UTF-8/de_DE.UTF-8/C/de_DE.UTF-8/de_DE.UTF-8
>
>     attached base packages:
>     [1] stats     graphics  grDevices utils     datasets methods base
>
>     other attached packages:
>     [1] lme4_1.1-12  Matrix_1.2-6 tibble_1.1   nlme_3.1-128
>     MASS_7.3-45 gstat_1.1-3  sp_1.2-3
>
>     loaded via a namespace (and not attached):
>      [1] Rcpp_0.12.6      lattice_0.20-33  intervals_0.15.1 FNN_1.1   
>           spacetime_1.1-5  zoo_1.7-13  assertthat_0.1 grid_3.3.1     
>      pacman_0.4.1  minqa_1.2.4      nloptr_1.0.4
>     [12] xts_0.9-7        splines_3.3.1    tools_3.3.1
>
>
>
>     Am 25.07.16 um 12:34 schrieb Phillip Alday:
>
>         Hi Patrick,
>
>         can you send us your sessionInfo()? Knowing the R version is
>         important, but a few other details matter for debugging this
>         type of problem!
>
>         Best,
>         Phillip
>
>
>             On 25 Jul 2016, at 19:55, Patrick Johann Schratz
>             <patrick.schratz at gmail.com
>             <mailto:patrick.schratz at gmail.com>> wrote:
>
>             Link to data:
>             <https://www.dropbox.com/s/yi3vf0bmqvydr8h/data.Rd?dl=0>
>             (1170 obs, 9 variables, .Rd file) [plain link in case sth
>             goes wrong
>             with the hyperlink:
>             https://www.dropbox.com/s/yi3vf0bmqvydr8h/data.Rd?dl=0]
>
>             Simply read it in using `readRDS(file)`.
>
>             I´m trying to setup a GLMM using the `glmmPQL` function
>             from the `MASS`
>             package including a random effects part and accounting for
>             spatial
>             autocorrelation. However, R (Version: 3.3.1) crashes upon
>             execution.
>
>                  library(nlme)
>
>                  # setup model formula
>                  fo <- hail ~ prec_nov_apr + t_min_nov_apr +
>             srad_nov_apr + age
>
>                  # setup corSpatial object
>                  correl = corSpatial(value = c(10000, 0.1), form = ~ry
>             + rx, nugget
>             = TRUE,
>                                      fixed = FALSE, type = "exponential")
>                  correl = Initialize(correl, data = d)
>
>                  # fit model
>                  fit5 <- MASS::glmmPQL(fo, random = ~1 | date, data = d,
>                                  correlation = correl, family = binomial)
>
>             What I tried so far:
>
>             - reduce number of observation
>             - play with `corSpatial` parameters (range and nugget)
>             - reduce number of fixed predictors
>             - execute code on Windows, Linux (Debian) and Mac R
>             installations
>
>
>             While I get no error message on my local pc (RStudio just
>             crashes),
>             running the script on a server returns the following error
>             message:
>
>             `R: malloc.c:3540: _int_malloc: Assertion (fwd->size &
>             0x4) == 0'
>             failed. Aborted`
>
>             Debugging leads me to a "glibc" c++ library problem. I
>             also run valgrind
>             on it. If you need the output, just ask!
>
>             Help ist highly appreciated!
>             Cheers, Patrick
>
>                     [[alternative HTML version deleted]]
>
>             _______________________________________________
>             R-sig-mixed-models at r-project.org
>             <mailto:R-sig-mixed-models at r-project.org> mailing list
>             https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
>     _______________________________________________
>     R-sig-mixed-models at r-project.org
>     <mailto:R-sig-mixed-models at 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