Notei que em R, Poisson e regressões binomiais negativas (RN) sempre parecem se encaixar nos mesmos coeficientes para preditores categóricos, mas não contínuos.
Por exemplo, aqui está uma regressão com um preditor categórico:
data(warpbreaks)
library(MASS)
rs1 = glm(breaks ~ tension, data=warpbreaks, family="poisson")
rs2 = glm.nb(breaks ~ tension, data=warpbreaks)
#compare coefficients
cbind("Poisson"=coef(rs1), "NB"=coef(rs2))
Aqui está um exemplo com um preditor contínuo, em que Poisson e NB se ajustam a diferentes coeficientes:
data(cars)
rs1 = glm(dist ~ speed, data=cars, family="poisson")
rs2 = glm.nb(dist ~ speed, data=cars)
#compare coefficients
cbind("Poisson"=coef(rs1), "NB"=coef(rs2))
(É claro que esses dados não são contados e os modelos não são significativos ...)
Então, recodifico o preditor em um fator e os dois modelos novamente se encaixam nos mesmos coeficientes:
library(Hmisc)
speedCat = cut2(cars$speed, g=5)
#you can change g to get a different number of bins
rs1 = glm(cars$dist ~ speedCat, family="poisson")
rs2 = glm.nb(cars$dist ~ speedCat)
#compare coefficients
cbind("Poisson"=coef(rs1), "NB"=coef(rs2))
No entanto, a regressão binomial negativa de Joseph Hilbe dá um exemplo (6.3, pág. 118-119) em que um preditor categórico, sexo, se encaixa com coeficientes ligeiramente diferentes pelos Poisson ( ) e NB ( ). Ele diz: “As taxas de incidência entre os modelos de Poisson e RN são bastante semelhantes. Isso não é surpreendente, dada a proximidade de [correspondente a em R] a zero. ”b = 0,881 α 1 / θ
Entendo isso, mas nos exemplos acima, summary(rs2)
diz-nos que é estimado em 9,16 e 7,93, respectivamente.
Então, por que os coeficientes são exatamente os mesmos? E por que apenas para preditores categóricos?
Editar # 1
Aqui está um exemplo com dois preditores não ortogonais. De fato, os coeficientes não são mais os mesmos:
data(cars)
#make random categorical predictor
set.seed(1); randomCats1 = sample( c("A","B","C"), length(cars$dist), replace=T)
set.seed(2); randomCats2 = sample( c("C","D","E"), length(cars$dist), replace=T)
rs1 = glm(dist ~ randomCats1 + randomCats2, data=cars, family="poisson")
rs2 = glm.nb(dist ~ randomCats1 + randomCats2, data=cars)
#compare coefficients
cbind("Poisson"=coef(rs1), "NB"=coef(rs2))
Além disso, a inclusão de outro preditor faz com que os modelos se ajustem a diferentes coeficientes, mesmo quando o novo preditor é contínuo. Então, isso tem algo a ver com a ortogonalidade das variáveis fictícias que criei no meu exemplo original?
rs1 = glm(dist ~ randomCats1 + speed, data=cars, family="poisson")
rs2 = glm.nb(dist ~ randomCats1 + speed, data=cars)
#compare coefficients
cbind("Poisson"=coef(rs1), "NB"=coef(rs2))
fonte
Respostas:
Você descobriu uma propriedade íntima, mas genérica, dos GLMs, ajustada pela máxima probabilidade . O resultado desaparece quando se considera o caso mais simples de todos: ajustar um único parâmetro a uma única observação!
Uma resposta frase : Se tudo o que interessa é encaixe meios separados para subconjuntos disjuntos de amostra, em seguida, MLG irá sempre originar μ J = ˉ y j para cada subconjunto j , de modo que a estrutura de erro real e parametrização da densidade tanto tornar irrelevante para a estimativa (ponto)!μ^j= y¯j j
Um pouco mais : Ajustar fatores categóricos ortogonais à máxima probabilidade é equivalente a ajustar meios separados para separar subconjuntos de nossa amostra, portanto, isso explica por que Poisson e GLMs binomiais negativos produzem as mesmas estimativas de parâmetros. De fato, o mesmo acontece se usamos regressão de Poisson, negbin, Gaussian, Gaussian inversa ou Gamma (veja abaixo). No caso de Poisson e negbin, a função de link padrão é o link de , mas esse é um arenque vermelho; Embora isso produza as mesmas estimativas brutas de parâmetro, veremos abaixo que essa propriedade realmente não tem nada a ver com a função de link.registro
Quando estamos interessados em uma parametrização com mais estrutura, ou que depende de preditores contínuos, a estrutura de erro assumida se torna relevante devido à relação de variância média da distribuição no que se refere aos parâmetros e à função não linear usada para modelar a condicional significa.
GLMs e famílias de dispersão exponencial: curso intensivo
Uma família de dispersão exponencial na forma natural é aquela que a densidade logarítmica é da forma
Aqui é o parâmetro natural e ν é o parâmetro de dispersão . Se v fosse conhecido, seria apenas uma família exponencial de um parâmetro padrão. Todos os GLMs considerados abaixo assumem um modelo de erro dessa família.θ ν ν
Considere uma amostra de uma única observação dessa família. Se nos encaixamos por máxima verossimilhança, temos que y = b ' ( θ ) , independentemente do valor de ν . Este estende-se facilmente ao caso de uma amostra iid uma vez que as probabilidades de log adicionar, obtendo-se ˉ y = b ' ( θ ) .θ y= b′( θ^) ν y¯= b′( θ^)
Mas também sabemos, devido à boa regularidade da densidade do log em função de , que ∂θ
Então, de fato b ′ ( θ ) = E Y = μ .
Uma vez que as estimativas de probabilidade máxima são invariante sob transformações, isto significa que para esta família de densidades.y¯= μ^
O que há de tão diferente nos preditores contínuos?
Quando os preditores são contínuos ou são categóricos, mas não podem ser reduzidos a uma forma ortogonal, a probabilidade não é mais fatorada em termos individuais com uma média separada, dependendo de um parâmetro separado. Neste ponto, a estrutura de erro e função de ligação não entram em jogo.
Dessa forma, a função de link e o modelo de erro assumido tornam-se relevantes para a estimativa.
Exemplo: o modelo de erro (quase) não importa
A partir da tabela, podemos ver que as estimativas dos parâmetros são idênticas , embora alguns desses GLMs sejam para dados discretos e outros sejam contínuos, e alguns sejam para dados não negativos, enquanto outros não.
Exemplo: a função de link (quase) não importa
A ressalva no cabeçalho refere-se simplesmente ao fato de que as estimativas brutas variarão com a função de link, mas as estimativas implícitas dos parâmetros médios não.
Código R
fonte
family=negative.binomial(theta=2)
"?Para ver o que está acontecendo aqui, é útil primeiro fazer a regressão sem a interceptação, pois uma interceptação em uma regressão categórica com apenas um preditor não tem sentido:
Como Poisson e regressões binomiais negativas especificam o log do parâmetro da média, para a regressão categórica, exponenciando os coeficientes, você obterá o parâmetro da média real para cada categoria:
Estes parâmetros correspondem às médias reais sobre os diferentes valores de categoria:
fonte
y~X+0
e tente novamente. :-)