anovaOD: Likelihood-Ratio Test Corrected for Overdispersion

Description

Compute likelihood-ratio test between a given model and a simpler model.

Usage

anovaOD(mod.simple, mod.complex, c.hat = 1, 
        nobs = NULL, ...)
# S3 method for glm
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitOccu
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitColExt
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitOccuRN
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitPCount
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitPCO
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitDS
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitGDS
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitOccuFP
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitMPois
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitGMM
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitGPC
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitOccuMS
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitOccuTTD
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitMMO
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for unmarkedFitDSO
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for glmerMod
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for maxlikeFit
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for multinom
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)
# S3 method for vglm
anovaOD(mod.simple, mod.complex, c.hat = 1,  
        nobs = NULL, ...)

Value

anovaOD returns an object of class anovaOD as a list with the following components:

form.simple: a character string of the parameters estimated in mod.simple.
form.complex: a character string of the parameters estimated in mod.complex.
c.hat: the c.hat estimate used to adjust the likelihood-ratio test.
devMat: a matrix storing as columns the number of parameters estimated (K), the log-likelihood of each model logLik, the difference in estimated parameters between the two models (Kdiff), minus twice the difference in log-likelihoods between the models (-2LL), the test statistic, and the associated P-value.

Arguments

mod.simple: an object of class glm, glmmTMB, maxlikeFit, mer, merMod, multinom, vglm, and various unmarkedFit classes containing the output of a model. This model should be a simpler version of mod.complex resulting from a deletion of certain terms (i.e., nested model).
mod.complex: an object of the same class as mod.simple.
c.hat: value of overdispersion parameter (i.e., variance inflation factor) such as that obtained from c_hat, mb.gof.test, or Nmix.gof.test. Typically, this value should be computed for the most complex model and applied to simpler models.
nobs: the number of observations used in the analysis. If nobs = NULL, the total number of rows are used as the sample size to compute the residual degrees of freedom as $nobs - K$, where K is the number of estimated parameters. This is relevant only for mixed models or various models of unmarkedFit classes where sample size is not straightforward. In such cases, one might use total number of observations or number of independent clusters (e.g., sites) as the value of nobs.
...: additional arguments passed to the function.

Author

Marc J. Mazerolle

Details

This function applies a correction for overdispersion on the likelihood-ratio test between a model and its simpler counterpart. The simpler model must be nested within the more complex model, typically as the result of deleting terms. You should supply the c.hat value of the most complex of the two models you are comparing.

When $1 < c.hat < 4$, the likelihood-ratio test is computed as:

$$LR = \frac{-2 * (LL.simple - LL.complex)}{(K.complex - K.simple) * c.hat}$$

where LL.simple and LL.complex are the log-likelihoods of the simple and complex models, respectively, and where K.complex and K.simple are the number of estimated parameters in each model. The test statistic is approximately distributed as $F_{K.complex - K.simple, n - K.complex}$, where n is the number of observations (i.e., nobs) used in the analysis (Venables and Ripley 2002).

When nobs = NULL, the number of observations is based on the number of rows of the data frame used in the analysis. For mixed models or various models of unmarkedFit, sample size is less straightforward, and nobs could be based on the total number of observations or on the number of independent clusters (e.g., sites), among other choices.

When c.hat = 1, the likelihood-ratio test simplifies to:

$$LR = -2 * (LL.simple - LL.complex)$$

where in this case the test statistic is distributed as a $\chi^2_{K.complex - K.simple}$ (McCullagh and Nelder 1989).

The function supports different model types such as Poisson GLM's and GLMM's, single-season and dynamic occupancy models (MacKenzie et al. 2002, 2003), and various N-mixture models (Royle 2004, Dail and Madsen 2011).

References

Dail, D., Madsen, L. (2011) Models for estimating abundance from repeated counts of an open population. Biometrics 67, 577--587.

MacKenzie, D. I., Nichols, J. D., Lachman, G. B., Droege, S., Royle, J. A., Langtimm, C. A. (2002) Estimating site occupancy rates when detection probabilities are less than one. Ecology 83, 2248--2255.

MacKenzie, D. I., Nichols, J. D., Hines, J. E., Knutson, M. G., Franklin, A. B. (2003) Estimating site occupancy, colonization, and local extinction when a species is detected imperfectly. Ecology 84, 2200--2207.

Mazerolle, M. J. (2006) Improving data analysis in herpetology: using Akaike's Information Criterion (AIC) to assess the strength of biological hypotheses. Amphibia-Reptilia 27, 169--180.

McCullagh, P., Nelder, J. A. (1989) Generalized Linear Models. Second edition. Chapman and Hall: New York.

Royle, J. A. (2004) N-mixture models for estimating population size from spatially replicated counts. Biometrics 60, 108--115.

Venables, W. N., Ripley, B. D. (2002) Modern Applied Statistics with S. Second edition. Springer-Verlag: New York.

Examples

Run this code

##anuran larvae example from Mazerolle (2006)
data(min.trap)
##assign "UPLAND" as the reference level as in Mazerolle (2006)          
min.trap$Type <- relevel(min.trap$Type, ref = "UPLAND") 

##run model
m1 <- glm(Num_anura ~ Type + log.Perimeter + Num_ranatra,
          family = poisson, offset = log(Effort),
          data = min.trap) 
##null model
m0 <- glm(Num_anura ~ 1,
          family = poisson, offset = log(Effort),
          data = min.trap) 

##check c-hat for global model
c_hat(m1) #uses Pearson's chi-square/df

##likelihood ratio test corrected for overdispersion
anovaOD(mod.simple = m0, mod.complex = m1, c.hat = c_hat(m1))
##compare without overdispersion correction
anovaOD(mod.simple = m0, mod.complex = m1)


##example with occupancy model
if (FALSE) {
##load unmarked package
if(require(unmarked)){
   
   data(bullfrog)
     
   ##detection data
   detections <- bullfrog[, 3:9]

   ##assemble in unmarkedFrameOccu
   bfrog <- unmarkedFrameOccu(y = detections)
     
   ##run model
   fm <- occu(~ 1 ~ Reed.presence, data = bfrog)
   ##null model
   fm0 <- occu(~ 1 ~ 1, data = bfrog)

   ##check GOF
   ##GOF <- mb.gof.test(fm, nsim = 1000)
   ##estimate of c-hat:  1.89

   ##display results after overdispersion adjustment
   anovaOD(fm0, fm, c.hat = 1.89)

   detach(package:unmarked)
}
}

Run the code above in your browser using DataLab