Generalized linear model with binomial distribution

Question

just a basic question, but is there a difference between generalized linear model with binomial distribution and binary logistic regression? Or it´s the same?

Answer 1

Binary logistic regression is a generalized linear model with the Bernoulli distribution. The Bernoulli distribution is just a special case of the binomial distribution. So they're not "the same" necessarily, but one is a special case of the other.

This is why in R there is no Bernoulli option in the glm() function. You use binomial for both binary data and binomial count data.

Answer 2

It comes down to what you mean by the same.

One way to understand two statistical procedures to be the same is that they always produce the same results. In this case, "results" ought to be understood as the same multivariate distribution of parameter estimates: that implies all predictions will be the same, as well as all standard errors, confidence intervals, and so forth.

Because a sum of independent Bernoulli variables is Binomial, and any Binomial variable can be expressed as the sum of independent Bernoulli variables (with the same parameter), the likelihoods will be identical and therefore Maximum Likelihood estimates will be the same in this sense.

Here is a demonstration in R. It begins by generating Bernoulli data for three values $x=0,1,2$:

logit.inv <- function(x) 1 / (1 + exp(-x))
x <- rep(0:2, 20)
beta <- c(-2,2)
set.seed(17)
y <- rbinom(length(x), 1, logit.inv(beta[1] + beta[2]*x))
library(data.table)
X <- data.table(x=x, y=y)

Each value of $x$ is sampled $20$ times, producing a dataset of 60 Bernoulli observations stored in X. Here they are. Superimposed on the plot as a black curve is the underlying graph of the Bernoulli parameter:

We may summarize these 60 observations by recording how many "successes" (values of $1$) and how many "failures" (values of $0$) were observed for each $x$.

X.s <- X[, .(successes=sum(y), failures=sum(1-y)), by=x]

This produces a dataset X.s with just three rows. Here is a visualization as a stacked bar chart:

From this we could recover a dataset equivalent to the first simply by generating $20$ new rows for each value of $x$. Some of those would have values of $0$ for $y$ and the remainder would have values of $1$. Their numbers are given by the failure and success counts, respectively.

The procedures to fit these datasets--the first with a Bernoulli GLM, the second with a Binomial GLM (and both with logistic links)--yield identical estimates and variance-covariance matrices. Here is partial output from the commands summary(fit) and summary(fit.s):

Call:
glm(formula = y ~ x, family = "binomial", data = X)

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -1.8071     0.5779  -3.127  0.00177 ** 
x             2.2882     0.5485   4.172 3.02e-05 ***


Call:
glm(formula = cbind(successes, failures) ~ x, family = "binomial", data = X.s)

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -1.8071     0.5779  -3.127  0.00177 ** 
x             2.2882     0.5485   4.172 3.02e-05 ***

The estimates and standard errors are identical. Both fits also produce the same variance-covariance matrix, as you can check using vcov(fit) and vcov(fit.s) (it's not worth reporting the details here).

However, the software computes different statistics intended for comparing similar models. For instance, the deviance summaries and associated AIC statistics differ. Therefore, although you will obtain identical results using either method, make sure to use a common approach when comparing models.

You probably noticed there was an enormous efficiency in data storage when the original sixty-row dataset X was collapsed into a three-row summary X.s. One way to understand the general Binomial GLM is that it's a convenient way to perform analysis of Bernoulli data -- that is, 0/1 responses -- without having to list the data explicitly.