lrtest: Likelihood Ratio Test of Nested Models

Description Usage Arguments Details Value Warning Note See Also Examples

View source: R/lrwaldtest.R

Description

lrtest is a generic function for carrying out likelihood ratio tests. The default method can be employed for comparing nested VGLMs (see details below).

Usage

1
2
3
 lrtest(object, ...)

 lrtest_vglm(object, ..., no.warning = FALSE, name = NULL)

Arguments

object

a vglm object. See below for details.

...

further object specifications passed to methods. See below for details.

no.warning

logical; if TRUE then no warning is issued. For example, setting TRUE might be a good idea when testing for linearity of a variable for a "pvgam" object.

name

a function for extracting a suitable name/description from a fitted model object. By default the name is queried by calling formula.

Details

lrtest is intended to be a generic function for comparisons of models via asymptotic likelihood ratio tests. The default method consecutively compares the fitted model object object with the models passed in .... Instead of passing the fitted model objects in ..., several other specifications are possible. The updating mechanism is the same as for waldtest() in lmtest: the models in ... can be specified as integers, characters (both for terms that should be eliminated from the previous model), update formulas or fitted model objects. Except for the last case, the existence of an update method is assumed. See waldtest() in lmtest for details.

Subsequently, an asymptotic likelihood ratio test for each two consecutive models is carried out: Twice the difference in log-likelihoods (as derived by the logLik methods) is compared with a Chi-squared distribution.

Value

An object of class "VGAManova" which contains a slot with the log-likelihood, degrees of freedom, the difference in degrees of freedom, likelihood ratio Chi-squared statistic and corresponding p value. These are printed by stats:::print.anova(); see anova.

Warning

Several VGAM family functions implement distributions which do not satisfying the usual regularity conditions needed for the LRT to work. No checking or warning is given for these.

Note

The code was adapted directly from lmtest (written by T. Hothorn, A. Zeileis, G. Millo, D. Mitchell) and made to work for VGLMs and S4. This help file also was adapted from lmtest.

Approximate LRTs might be applied to VGAMs, as produced by vgam, but it is probably better in inference to use vglm with regression splines (bs and ns). This methods function should not be applied to other models such as those produced by rrvglm, by cqo, by cao.

See Also

lmtest, vglm, lrt.stat.vlm, score.stat.vlm, wald.stat.vlm, anova.vglm.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
set.seed(1)
pneumo <- transform(pneumo, let = log(exposure.time), x3 = runif(nrow(pneumo)))
fit1 <- vglm(cbind(normal, mild, severe) ~ let     , propodds, data = pneumo)
fit2 <- vglm(cbind(normal, mild, severe) ~ let + x3, propodds, data = pneumo)
fit3 <- vglm(cbind(normal, mild, severe) ~ let     , cumulative, data = pneumo)
# Various equivalent specifications of the LR test for testing x3
(ans1 <- lrtest(fit2, fit1))
ans2 <- lrtest(fit2, 2)
ans3 <- lrtest(fit2, "x3")
ans4 <- lrtest(fit2, . ~ . - x3)
c(all.equal(ans1, ans2), all.equal(ans1, ans3), all.equal(ans1, ans4))

# Doing it manually
(testStatistic <- 2 * (logLik(fit2) - logLik(fit1)))
(mypval <- pchisq(testStatistic, df = length(coef(fit2)) - length(coef(fit1)),
                  lower.tail = FALSE))

(ans4 <- lrtest(fit3, fit1))  # Test proportional odds (parallelism) assumption

Example output

Loading required package: stats4
Loading required package: splines
Likelihood ratio test

Model 1: cbind(normal, mild, severe) ~ let + x3
Model 2: cbind(normal, mild, severe) ~ let
  #Df  LogLik Df  Chisq Pr(>Chisq)
1  12 -25.087                     
2  13 -25.090  1 0.0076     0.9304
[1] TRUE TRUE TRUE
[1] 0.007619612
[1] 0.9304407
Likelihood ratio test

Model 1: cbind(normal, mild, severe) ~ let
Model 2: cbind(normal, mild, severe) ~ let
  #Df  LogLik Df  Chisq Pr(>Chisq)
1  12 -25.019                     
2  13 -25.090  1 0.1424     0.7059

VGAM documentation built on Jan. 16, 2021, 5:21 p.m.