4
$\begingroup$

I am interested in the behavioral response of floral visitors to a treatment, applied in a paired fashion within plants. That is, one stem on each plant receives the treatment, and another stem serves as a control. The response that I am interested in is total number of flowers visited per day. Each plant is only used for one day and we sampled for a total of three days.

The number of flowers on each stem varies within and between plants, and will affect visitation. However, I am not interested in this effect.

I do not know whether it is more appropriate to specify a glmm, with flower number nested within individual plant as a random effect, and date as a separate random effect, such as this:

mod1<-glmer(data=data, total_flr_vis ~ treatment + (1|date) + (1|individual/flr_num), family=poisson)

or whether it is more appropriate to specify a glm such as this:

mod2<-glm(data=data, total_flr_vis ~ treatment + individual + date + flr_num, family=quasipoisson)

(I've used quasipoisson here as the data are over-dispersed when the effects are all fixed)

I think that it is more philosophically sound to use the glmm, and would go with that if the results of each model were qualitatively similar. However, they are very different (output below). Specifically, the glmm tells me that there is a significant effect of treatment, whereas the glm does not. Thus, I would like to be extra careful that I am using the most appropriate model.

Any advice is most welcome! Many thanks!

> mod1<-glmer(data=newdata, total_flr_vis ~ treatment + (1|date) +    (1|individual/flr_num), family=poisson)
> mod1
Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: poisson  ( log )
Formula: total_flr_vis ~ treatment + (1 | date) + (1 | individual/flr_num)
   Data: newdata
      AIC       BIC    logLik  deviance  df.resid 
 457.3714  472.4908 -223.6857  447.3714       147 
Random effects:
 Groups             Name        Std.Dev.
 flr_num:individual (Intercept) 2.864   
 individual         (Intercept) 1.008   
 date               (Intercept) 0.000   
Number of obs: 152, groups:  flr_num:individual, 149; individual, 42; date, 3
 Fixed Effects:
(Intercept)   treatmentR  
 -1.712       -1.503  
> summary(mod1)
Generalized linear mixed model fit by maximum likelihood (Laplace  Approximation) ['glmerMod']
 Family: poisson  ( log )
Formula: total_flr_vis ~ treatment + (1 | date) + (1 | individual/flr_num)
   Data: newdata

     AIC      BIC   logLik deviance df.resid 
   457.4    472.5   -223.7    447.4      147 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.0196 -0.2783 -0.1751  0.1270  0.6132 

Random effects:
 Groups             Name        Variance Std.Dev.
 flr_num:individual (Intercept) 8.201    2.864   
 individual         (Intercept) 1.015    1.008   
 date               (Intercept) 0.000    0.000   
Number of obs: 152, groups:  flr_num:individual, 149; individual, 42; date, 3

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -1.7117     0.6530  -2.621 0.008757 ** 
treatmentR   -1.5032     0.3959  -3.797 0.000146 ***
 ---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
           (Intr)
treatmentR -0.165

Compared to:

> mod2 = glm(data=newdata, total_flr_vis ~ treatment + individual + date + flr_num, family=quasipoisson)
> mod2

Call:  glm(formula = total_flr_vis ~ treatment + individual + date + 
    flr_num, family = quasipoisson, data = newdata)

Coefficients:
(Intercept)   treatmentR   individual         date      flr_num  
  -0.8641906   -0.4262042   -0.0004565    0.0072634    0.0150606  

Degrees of Freedom: 151 Total (i.e. Null);  147 Residual
Null Deviance:      1011 
Residual Deviance: 963.4    AIC: NA
> summary(mod2)

Call:
glm(formula = total_flr_vis ~ treatment + individual + date + 
    flr_num, family = quasipoisson, data = newdata)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-3.4249  -1.9918  -1.6997  -0.5174  10.9899  

Coefficients:
              Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.8641906 20.9133728  -0.041    0.967
treatmentR  -0.4262042  0.3971152  -1.073    0.285
individual  -0.0004565  0.0179801  -0.025    0.980
date         0.0072634  0.1155358   0.063    0.950
flr_num      0.0150606  0.0100614   1.497    0.137

(Dispersion parameter for quasipoisson family taken to be 11.20508)

    Null deviance: 1011.19  on 151  degrees of freedom
Residual deviance:  963.39  on 147  degrees of freedom
AIC: NA

Number of Fisher Scoring iterations: 7
$\endgroup$
  • 1
    $\begingroup$ Would it be possible to post a sample of your data, e.g. first 10-20 rows of your dataframe? There are multiple things about your models that don't make a lot of sense to me. (1) Does date only have three levels, are the dates the same for all flowers? Then including date as a random effect is arguably suboptimal, it's too few levels for a random effect. (2) What is individual? Is it an id number of a flower? Including it as a continuous variable in your glm model does not make any sense at all. You should include it as a factor. (3) Similarly, flr_num as random is strange. Etc. $\endgroup$ – amoeba says Reinstate Monica Oct 17 '16 at 23:29
  • $\begingroup$ @amoeba (1) There are 14 levels of date; (2) Yes, individual is an ID (plant, not flower). I'm sorry that it isn't clear from my code, but it is defined as a factor; (3) I included flr_num as random as it represents the display size (number of open flowers at the time of sampling) of each individual on each day of sampling. It changes at each sampling date; the levels would not be repeated exactly were I to run the experiment again; the value is expected to have an influence on visitor behavior, but is not something that I am interested in testing. $\endgroup$ – JKO Mar 31 '17 at 20:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.