# Count data with mixed effects

EDIT: Thanks for the help so far. I have updated the question based on further work.

I am interested in finding any difference in response of counts of two competing species to dryness and elevation. Both are caught in the same type of trap, but trap placement (bushes vs buildings) probably biases toward one species or the other, so each trap placement is rated according to indoors/outdoors on a scale of 1-3. A trap does not move once it is placed. Traps are nested in sites, dryness is measured at each site/visit combination, and elevation is of course constant at each site.

My data includes the number of each species caught in 4-12 traps, at 8-12 sites during 6-8 visits to each site. The visits and sites are far enough spaced that I assume they are independent. Most data are zeros, especially for Sp_A which is rare. Below is the structure of the data.

R> str(mydata)
'data.frame':   300 obs. of  8 variables:
$count : num 0 5 1 0 1 1 0 0 0 0 ...$ species  : Factor w/ 2 levels "a","b": 1 1 1 1 1 1 1 1 1 1 ...
$elevation: int 1 1 1 1 1 1 1 1 1 1 ...$ dryness  : num  0.179 0.179 0.179 0.179 0.179 ...
$site : Factor w/ 5 levels "a","b","c","d",..: 1 1 1 1 1 1 1 1 1 1 ...$ trap     : Factor w/ 50 levels "MT10a","MT10b",..: 6 11 16 21 26 31 36 41 46 1 ...
$visit : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...$ in_out   : num  1 1 1 1 1 1 1 1 1 1 ...


As an initial attempt, I had no trouble using R function lmer to model Sp_A presence/absence with Sp_B as a covariate: (SpA>0) ~ Sp_B + elevation + (1|site/trap) + (1|visit)

However on suggestion from R-sig-mixed-models list, and also "ils" from this site, I have reshaped the data as seen above, and tried to model count with species as a factor in the model. I have used two different packages (lme4 and glmmADMB). The functions handle the data ok in simple models, until I add the factor "species." Below are the R code and error messages. These same errors also happen when I use simpler dummy data. Any ideas on syntax or other packages?

require(lme4)
glmera1 <- glmer(count~elevation*species*in_out
+ (1|visit) + (1|site/trap),
data=mydata, family="poisson")
Error in asMethod(object) : matrix is not symmetric [1,2]
In addition: Warning messages:
1: In mer_finalize(ans) :
Cholmod warning 'not positive definite' at file:../Cholesky/t_cholmod_rowfac.c, line 432
2: In mer_finalize(ans) :
Cholmod warning 'not positive definite' at file:../Cholesky/t_cholmod_rowfac.c, line 432
3: In mer_finalize(ans) : singular convergence (7)

group="site", data=mydata, family="nbinom")
Error in glmm.admb(count ~ elevation * species * in_out, random = ~visit,  :
The function maximizer failed
In addition: Warning message:
running command './nbmm -maxfn 500 ' had status 1

• I'm not aware of any question that covers "exactly the same ground." While these criteria are sometimes negative, my question is, have I done something wrong? – J. Win. Apr 7 '11 at 22:52
• @J. Winchester I'm afraid your question still appears to focus on AIC, and responses to the related question on AIC apply as well here. However, I feel you are also wondering whether your lmer syntax allows you to test the hypothesis formulated in your first paragraph. If this is the case, you can simply edit the title and your question and focus on this issue, rather than AIC. You might also consider adding the output of summary(lmera.1). Finally, the outcome you consider is a ratio of counts (% Sp_A): Is it usual in your field to treat this as a continuous response variable? – chl Apr 8 '11 at 19:30
• @chl. Thank you. I have rewritten the question and title; regarding your last sentence, I cannot say that it is usual, but I do not have a better idea. – J. Win. Apr 12 '11 at 8:46
• @J. Winchester: It might be a crazy idea, but why not treat abundace as the dependent variable and add species_identity as a factor into the model? Then, when you have estimated the magnitude of species_identity you will be able to compare the abundance of both species in relative terms. This way you won't need the model comparison stage and things will be much more simple. I also don't see the purpose of visit_day as an estimable parameter; the species abundance will clearly depend on it as an exposure parameter, but then how one will interpret the dryness * visit_day interaction? – ils Apr 12 '11 at 11:34
• @ils Thank you, I have tried that and updated the question. I wonder what is the difference between your suggestion, and using lmer(cbind(sp_a, sp_b)~elevation*dryness ... ) which is a format used at the bottom of the lmer help page. – J. Win. Apr 20 '11 at 23:38