r - Error when Fitting a glmer with poisson error structure -
i hope can me. i'm trying conduct analysis examines number of samples of hymenoptera caught on elevational gradient. want examine possibility of uni-modal distribution in relation elevation, linear distribution. hence including i(altitude^2)
explanatory variable in analysis.
i trying run following model includes poisson error structure (as dealing count data) , date , trap type (trap
) random effects.
model7 <- glmer(no.specimens~altitude+i(altitude^2)+(1|date)+(1|trap), family="poisson",data=santa.lucia,na.action=na.omit)
however keep receiving following error message:
error: (maxstephalfit) pirls step-halvings failed reduce deviance in pwrssupdate in addition: warning messages: 1: predictor variables on different scales: consider rescaling 2: in pwrssupdate(pp, resp, tolpwrss, gqmat, compdev, fac, verbose) : cholmod warning 'not positive definite' @ file:../cholesky/t_cholmod_rowfac.c, line 431 3: in pwrssupdate(pp, resp, tolpwrss, gqmat, compdev, fac, verbose) : cholmod warning 'not positive definite' @ file:../cholesky/t_cholmod_rowfac.c, line 431
clearly making big mistakes. can me figure out going wrong?
here structure of dataframe:
str(santa.lucia) 'data.frame': 97 obs. of 6 variables: $ date : factor w/ 8 levels "01-sep-2014",..: 6 6 6 6 6 6 6 6 6 6 ... $ trap.no : factor w/ 85 levels "n1","n10","n11",..: 23 48 51 14 17 20 24 27 30 33 ... $ altitude : int 1558 1635 1703 1771 1840 1929 1990 2047 2112 2193 ... $ trail : factor w/ 3 levels "cascadas","limones",..: 1 1 1 1 1 3 3 3 3 3 ... $ no.specimens: int 1 0 2 2 3 4 5 0 1 1 ... $ trap : factor w/ 2 levels "net","pan": 2 2 2 2 2 2 2 2 2 2 ...
and here complete data.set (these preliminary analyses)
date trap.no altitude trail no.specimens trap 1 28-aug-2014 w2 1558 cascadas 1 pan 2 28-aug-2014 w5 1635 cascadas 0 pan 3 28-aug-2014 w8 1703 cascadas 2 pan 4 28-aug-2014 w11 1771 cascadas 2 pan 5 28-aug-2014 w14 1840 cascadas 3 pan 6 28-aug-2014 w17 1929 tower 4 pan 7 28-aug-2014 w20 1990 tower 5 pan 8 28-aug-2014 w23 2047 tower 0 pan 9 28-aug-2014 w26 2112 tower 1 pan 10 28-aug-2014 w29 2193 tower 1 pan 11 28-aug-2014 w32 2255 tower 0 pan 12 30-aug-2014 n1 1562 cascadas 5 net 13 30-aug-2014 n2 1635 cascadas 0 net 14 30-aug-2014 n3 1723 cascadas 2 net 15 30-aug-2014 n4 1779 cascadas 0 net 16 30-aug-2014 n5 1842 cascadas 3 net 17 30-aug-2014 n6 1924 tower 2 net 18 30-aug-2014 n7 1979 tower 2 net 19 30-aug-2014 n8 2046 tower 0 net 20 30-aug-2014 n9 2110 tower 0 net 21 30-aug-2014 n10 2185 tower 0 net 22 30-aug-2014 n11 2241 tower 0 net 23 31-aug-2014 n1 1562 cascadas 1 net 24 31-aug-2014 n2 1635 cascadas 1 net 25 31-aug-2014 n3 1723 cascadas 0 net 26 31-aug-2014 n4 1779 cascadas 0 net 27 31-aug-2014 n5 1842 cascadas 0 net 28 31-aug-2014 n6 1924 tower 0 net 29 31-aug-2014 n7 1979 tower 7 net 30 31-aug-2014 n8 2046 tower 4 net 31 31-aug-2014 n9 2110 tower 6 net 32 31-aug-2014 n10 2185 tower 1 net 33 31-aug-2014 n11 2241 tower 1 net 34 01-sep-2014 w1 1539 cascadas 0 pan 35 01-sep-2014 w2 1558 cascadas 0 pan 36 01-sep-2014 w3 1585 cascadas 2 pan 37 01-sep-2014 w4 1604 cascadas 0 pan 38 01-sep-2014 w5 1623 cascadas 1 pan 39 01-sep-2014 w6 1666 cascadas 4 pan 40 01-sep-2014 w7 1699 cascadas 0 pan 41 01-sep-2014 w8 1703 cascadas 0 pan 42 01-sep-2014 w9 1746 cascadas 1 pan 43 01-sep-2014 w10 1762 cascadas 0 pan 44 01-sep-2014 w11 1771 cascadas 0 pan 45 01-sep-2014 w12 1796 cascadas 1 pan 46 01-sep-2014 w13 1825 cascadas 0 pan 47 01-sep-2014 w14 1840 tower 4 pan 48 01-sep-2014 w15 1859 tower 2 pan 49 01-sep-2014 w16 1889 tower 2 pan 50 01-sep-2014 w17 1929 tower 0 pan 51 01-sep-2014 w18 1956 tower 0 pan 52 01-sep-2014 w19 1990 tower 1 pan 53 01-sep-2014 w20 2002 tower 3 pan 54 01-sep-2014 w21 2023 tower 2 pan 55 01-sep-2014 w22 2047 tower 0 pan 56 01-sep-2014 w23 2068 tower 1 pan 57 01-sep-2014 w24 2084 tower 0 pan 58 01-sep-2014 w25 2112 tower 1 pan 59 01-sep-2014 w26 2136 tower 0 pan 60 01-sep-2014 w27 2150 tower 1 pan 61 01-sep-2014 w28 2193 tower 1 pan 62 01-sep-2014 w29 2219 tower 0 pan 63 01-sep-2014 w30 2227 tower 1 pan 64 01-sep-2014 w31 2255 tower 0 pan 85 03/06/2015 wt47 1901 tower 2 pan 86 03/06/2015 wt48 1938 tower 2 pan 87 03/06/2015 wt49 1963 tower 2 pan 88 03/06/2015 wt50 1986 tower 0 pan 89 03/06/2015 wt51 2012 tower 9 pan 90 03/06/2015 wt52 2033 tower 0 pan 91 03/06/2015 wt53 2050 tower 4 pan 92 03/06/2015 wt54 2081 tower 2 pan 93 03/06/2015 wt55 2107 tower 1 pan 94 03/06/2015 wt56 2128 tower 4 pan 95 03/06/2015 wt57 2155 tower 0 pan 96 03/06/2015 wt58 2179 tower 2 pan 97 03/06/2015 wt59 2214 tower 0 pan 98 03/06/2015 wt60 2233 tower 0 pan 99 03/06/2015 wt61 2261 tower 0 pan 100 03/06/2015 wt62 2278 tower 0 pan 101 03/06/2015 wt63 2300 tower 0 pan 102 04/06/2015 wt31 1497 cascadas 0 pan 103 04/06/2015 wt32 1544 cascadas 1 pan 104 04/06/2015 wt33 1568 cascadas 1 pan 105 04/06/2015 wt34 1574 cascadas 0 pan 106 04/06/2015 wt35 1608 cascadas 5 pan 107 04/06/2015 wt36 1630 cascadas 3 pan 108 04/06/2015 wt37 1642 cascadas 0 pan 109 04/06/2015 wt38 1672 cascadas 5 pan 110 04/06/2015 wt39 1685 cascadas 6 pan 111 04/06/2015 wt40 1723 cascadas 3 pan 112 04/06/2015 wt41 1744 cascadas 2 pan 113 04/06/2015 wt42 1781 cascadas 1 pan 114 04/06/2015 wt43 1794 cascadas 2 pan 115 04/06/2015 wt44 1833 cascadas 0 pan 116 04/06/2015 wt45 1855 cascadas 4 pan 117 04/06/2015 wt46 1876 cascadas 2 pan
you're there. @bondeddust suggests, it's not practical use two-level factor (trap
) random effect; in fact, doesn't seem right in principle either (the levels of trap
not arbitrary/randomly chosen/exchangeable). when tried model quadratic altitude, fixed effect of trap, , random effect of date
, warned might want rescale parameter:
some predictor variables on different scales: consider rescaling
(you saw warning mixed in error messages). continuous (and hence worth rescaling) predictor altitude
, centered , scaled scale()
(the disadvantage changes quantitative interpretation of coefficients, model practically identical). added observation-level random effect allow overdispersion.
the results seem ok, , agree picture.
library(lme4) santa.lucia <- transform(santa.lucia, scalt=scale(altitude), obs=factor(seq(nrow(santa.lucia)))) model7 <- glmer(no.specimens~scalt+i(scalt^2)+trap+(1|date)+(1|obs), family="poisson",data=santa.lucia,na.action=na.omit) summary(model7) ## random effects: ## groups name variance std.dev. ## obs (intercept) 0.64712 0.8044 ## date (intercept) 0.02029 0.1425 ## number of obs: 97, groups: obs, 97; date, 6 ## ## fixed effects: ## estimate std. error z value pr(>|z|) ## (intercept) 0.53166 0.31556 1.685 0.09202 . ## scalt -0.22867 0.14898 -1.535 0.12480 ## i(scalt^2) -0.52840 0.16355 -3.231 0.00123 ** ## trappan -0.01853 0.32487 -0.057 0.95451
test quadratic term comparing model lacks ...
model7r <- update(model7, . ~ . - i(scalt^2)) ## convergence warning, ok ... anova(model7,model7r)
on principle might worth looking @ interaction between quadratic altitude model , trap (allowing different altitude trends trap type), picture suggests won't ...
library(ggplot2); theme_set(theme_bw()) ggplot(santa.lucia,aes(altitude,no.specimens,colour=trap))+ stat_sum(aes(size=factor(..n..)))+ scale_size_discrete(range=c(2,4))+ geom_line(aes(group=date),colour="gray",alpha=0.3)+ geom_smooth(method="gam",family="quasipoisson", formula=y~poly(x,2))+ geom_smooth(method="gam",family="quasipoisson", formula=y~poly(x,2),se=false, aes(group=1),colour="black")
Comments
Post a Comment