Zaloguj prawdopodobieństwo dla GLM

W poniższym kodzie wykonuję regresję logistyczną zgrupowanych danych za pomocą glm i „ręcznie” za pomocą mle2. Dlaczego funkcja logLik w R daje mi logarytm logLik (fit.glm) = - 2,336, który jest inny niż logLik (fit.ml) = - 5,514, który otrzymuję ręcznie?

library(bbmle)

#successes in first column, failures in second
Y <- matrix(c(1,2,4,3,2,0),3,2)

#predictor
X <- c(0,1,2)

#use glm
fit.glm <- glm(Y ~ X,family=binomial (link=logit))
summary(fit.glm)

#use mle2
invlogit <- function(x) { exp(x) / (1+exp(x))}
nloglike <- function(a,b) {
  L <- 0
  for (i in 1:n){
     L <- L + sum(y[i,1]*log(invlogit(a+b*x[i])) + 
               y[i,2]*log(1-invlogit(a+b*x[i])))
  }
 return(-L) 
}  

fit.ml <- mle2(nloglike,
           start=list(
             a=-1.5,
             b=2),
           data=list(
             x=X,
             y=Y,
             n=length(X)),
           method="Nelder-Mead",
           skip.hessian=FALSE)
summary(fit.ml)

#log likelihoods
logLik(fit.glm)
logLik(fit.ml)


y <- Y
x <- X
n <- length(x)
nloglike(coef(fit.glm)[1],coef(fit.glm)[2])
nloglike(coef(fit.ml)[1],coef(fit.ml)[2])

r self-study generalized-linear-model

— Tomek
źródło

Częstym powodem takich różnic jest fakt, że prawdopodobieństwo definiuje się tylko do stałej multiplikatywnej : „ Dokładniej zatem funkcja prawdopodobieństwa jest dowolnym przedstawicielem klasy równoważności funkcji gdzie stała proporcjonalności nie może zależeć od i musi być taka sama dla wszystkich funkcji wiarygodności używanych w dowolnym porównanie. $\mathcal{L} \in \left\lbrace \alpha \; P_\theta: \alpha > 0 \right\rbrace, \,$ $α > 0$ $θ$ ”Prawdopodobieństwo logarytmiczne może z kolei zostać przesunięte o dowolną stałą. ... (ctd)

— Glen_b

(ctd) ... To nie znaczy, że to wyjaśnienie tej konkretnej różnicy, ale jest to częsty powód różnic między tym, jak różne funkcje dają różne prawdopodobieństwa.

— Glen_b

Niepoprawnie założyłem, że prawdopodobieństwo dziennika zostało zdefiniowane w jądrze pliku pdf i dlatego jest unikalne dla tego problemu.

— Tom

Warto to jednak zbadać, ponieważ czasem wyjaśnienie to coś innego.

— Glen_b

Wydaje się, że funkcja logLik w R oblicza to, co określa się w SAS jako „funkcję pełnego prawdopodobieństwa”, która w tym przypadku obejmuje współczynnik dwumianowy. Nie uwzględniłem współczynnika dwumianowego w obliczeniach mle2, ponieważ nie ma to wpływu na oszacowania parametrów. Po dodaniu tej stałej do prawdopodobieństwa dziennika w obliczeniach mle2, glm i mle2 zgadzają się.

— Tomek
źródło

(+1) Dziękujemy za sprawdzenie i opublikowanie rezolucji po jej ustaleniu. Twoje zdrowie.

— kardynał