I wasn't able to find any satisfying answer about that topic. I hope someone who understand correctly the subject could enlighten this shadow.

This is not very important, just for the sake of notations: In the following we work with exponential families: $$ \exp \{ \frac {y \theta - b( \theta ) } { a(\phi)} + c(y,\phi) \} $$ the deviance is defined as the difference in likelihood between the model and the saturated model: $$D = 2 \phi \left ( l(y, \phi , y ) - l ( \hat{\mu}, \phi, y ) \right ) $$ the scaled deviance as : $$ D^* = \frac D {\phi} $$

When one compares two nested models, in particular we are trying to check the goodness of fit, and we are trying to test the null hypothesis (smaller model is doing a better approximation), we can compute $D_0^* - D_1^* $ where $D_0^*$ is the scaled deviance of the smaller model. This, should be approximately, be following a $ \chi^2 $ distribution, of the difference in number of parameters.

My question is the following, how to test whether that approximation is good ? My professor is not very precise on that topic, he mentions doing "simulations". What does it mean please?


That is a large topic, and the quality of the approximation must be studied, generally, on a case-by-case basis. That's why simulation is a useful approach. So the best way to answer your question is by an example, which you can adapt for your data and models.

So I will simulate, in R, some data for a logistic regression, and I will simulate binomial responses from the null model, that is, all regression parameters (except the intercept) assumed to be zero. Then I fit the logistic model, calculate the deviance, repeat many times, and plot the histogram, with the chi-square approximation overlaid. For my case the result is (based on 10000 replications)

histogram of simulated deviances

which looks rather good. A qqplot could also be useful:

qqplot against chisquare distribution

This looks like a very good fit. Maybe some clumpiness (discreteness), but no systematic deviations. But also, the far tail is poorly represented, to check quality of fit there, maybe a larger simulation.

You can modify the code below for your experiments:

set.seed(7*11*13)  # My public seed
X <- data.frame(x1=rnorm(200, 20, sd=3),
                x2=factor(rep(1:3, c(100, 50, 50))),
                y=rbinom(200, 1, 0.25))

N <- 10000 # Simulating the null distribution

sim_dev <- function(N) {
    res <- numeric(N)
    for (i in seq_along(res)) {
        X$y <- rbinom(200, 1, 0.25)
        mod <- glm(y  ~ x1+x2, data=X, family=binomial)
        res[i] <- with(mod,  null.deviance-deviance)

deviance <- sim_dev(N)

mean(deviance); var(deviance)
 6.090653  # Rather close to theoretical values 3, 6

hist(deviance, prob=TRUE, breaks="FD")
plot( function(x) dchisq(x, df=3), col="red", add=TRUE, from=0, to=22) 
qqplot(qchisq(ppoints(N), df=3), deviance)
qqline(deviance, distribution=function(x)qchisq(x, df=3))
title("qqplot against chisquare distribution (df=3)") 
  • 1
    $\begingroup$ thank you ! So you would plot the distribution of the deviance for the model right? Which in your example, is the null model. If it looks chi square, then we can confidently say that it is a good approximation. Is it correct? $\endgroup$ – Marine Galantin May 20 '20 at 12:35
  • 1
    $\begingroup$ Yes, more or less. One could also do a more formal test, or a qqplot against the chisquare distribution could be a better description of eventual differences. But beware to do the simulations for cases very close to your actual data. $\endgroup$ – kjetil b halvorsen May 20 '20 at 17:18
  • $\begingroup$ Quick question (sorry if it s obvious, I am still learning). The CHI^2 is of 3DF because, there is one continuous parameter and then 3 categorical variables, adding 2 regressors to the model (hence it adds 2DF). To this we add one degree of freedom which is the intercept, degree that is cancelled out by the null model. This gives a total of 3DF as you have written. Is it correct ? $\endgroup$ – Marine Galantin May 20 '20 at 21:27
  • 1
    $\begingroup$ Correct. ......... $\endgroup$ – kjetil b halvorsen May 20 '20 at 23:48

Your Answer

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

Not the answer you're looking for? Browse other questions tagged or ask your own question.