Regressão logística BIC: qual é o N certo?

7

TL; DR: Qual é correto para o BIC na regressão logística, no binômio agregado ou no Bernoulli ?NN

ATUALIZAÇÕES NO INFERIOR

Suponha que eu tenha um conjunto de dados ao qual gostaria de aplicar a regressão logística. Por uma questão de exemplo, suponha que grupos com participantes cada, para um total de . O resultado é 0 ou 1. Por exemplo, o seguinte conjunto de dados (código R):j=5m=100n=500

library(dplyr)
library(tidyr)


set.seed(45)
d <- tibble(y = rbinom(500, 1, .5),
            x = factor(rep(LETTERS[1:5], each = 100))) 

Existem duas maneiras de representar isso: como é, acima, tratar todas as observações como uma variável aleatória de Bernoulli, ou agregar observações em grupos e tratar cada observação como binomial. O número de linhas no conjunto de dados será 500 na primeira instância e 5 na segunda.

Eu posso construir o conjunto de dados agregados:

d %>% 
  group_by(x, y) %>% 
  summarise(n = n()) %>%
  spread(y, n) %>%
  rename(f = `0`, s = `1`) %>%
  mutate(n = s + f) -> d_agg

Posso ajustar a regressão logística usando os dois conjuntos de dados em R:

g_bern  <- glm(y ~ x,          data=d,     family=binomial)
g_binom <- glm(cbind(s,f) ~ x, data=d_agg, family=binomial)

ATUALIZAÇÃO 2: Agora, cabemos apenas os modelos de interceptação:

g_bern0  <- glm(y ~ 1,          data=d,     family=binomial)
g_binom0 <- glm(cbind(s,f) ~ 1, data=d_agg, family=binomial)

e calcular a AIC:

> AIC(g_bern)  
# [1] 694.6011
> AIC(g_binom)  
# [1] 35.22172

que obviamente diferem por uma constante

2*sum(lchoose(d_agg$n, d_agg$s))  # [1] 659.3794

conforme o esperado (consulte: Regressão logística: Bernoulli vs. Variáveis ​​de resposta binomial ).

No entanto, os BICs diferem por essa constante AND um fator que depende do "número de observações" e o número de observações difere em cada:

> BIC(g_bern)    
# [1] 715.6742
> BIC(g_binom)  
# [1] 33.26891
> nobs(g_bern)   
# [1] 500
> nobs(g_binom)  
# [1] 5

Apenas para confirmar, podemos recalcular o BIC para ambos:

> -2*logLik(g_bern) + attr(logLik(g_bern),"df")*log(nobs(g_bern))
# 'log Lik.' 715.6742 (df=5)
> -2*logLik(g_binom) + attr(logLik(g_binom),"df")*log(nobs(g_binom))
# 'log Lik.' 33.26891 (df=5)

e, na verdade o único lugar estes dois números diferem é .N

ATUALIZAÇÃO 2: Quando tentamos avaliar o fator x, vemos uma discordância que SÓ é atribuível ao número de observações:

> BIC(g_bern0) - BIC(g_bern)
# [1] -17.66498
> BIC(g_binom0) - BIC(g_binom)
# [1] 0.7556999

ATUALIZAÇÃO 2: Como esperado, os AICs são consistentes:

> AIC(g_bern0) - AIC(g_bern)
# [1] -0.8065485
> AIC(g_binom0) - AIC(g_binom)
# [1] -0.8065485

Isso me surpreende, pois acho que R "saberia" qual dos dois usar para evitar ambiguidade. Tem a mesma informação nos dois casos.

Qual é o "certo"? Ou a BIC é realmente tão arbitrária?

ATUALIZAÇÃO: Não estou tentando comparar o Bernoulli ao modelo binomial. Este é apenas um exemplo de brinquedo. Eu tenho um conjunto de comparações onde importa qual configuração eu uso, porque as penalidades para são diferentes. Eu tenho dois conjuntos de comparações de modelos e o modelo vencedor muda com base na penalidade de , mesmo que esses pareçam ser os mesmos conjuntos de modelos .NN

ATUALIZAÇÕES 2 e 3: Adicionadas as comparações ao modelo somente de interceptação e alteradas as sementes aleatórias para obter uma diferença de sinal no BIC.

Molho de salada
fonte

Respostas:

6

O BIC (e o AIC) são medidas relativas para comparar modelos. No entanto, não faz sentido comparar o que, de outra forma, é o mesmo modelo entre usar uma resposta agregada versus uma resposta desagregada. Também não faria sentido comparar modelos que seriam diferentes (por exemplo, diferentes regressores), mas onde um modelo usa uma resposta agregada e o outro modelo usa uma versão desagregada da resposta. Desde que os dois modelos comparados representem a variável de resposta no mesmo formato, tudo ficará bem. Note-se que os dois formatos são em última análise, equivalente a eles contêm as mesmas informações e, principalmente, basta olhar diferente do lado de fora, ver: Formato de entrada para a resposta em glm binomial em R .

- Reinstate Monica
fonte
Apenas para esclarecer: não estou tentando comparar g_bern com g_binom. Este é um exemplo de brinquedo. Eu tenho um aplicativo em que essa penalidade é importante para qual modelo vence. Portanto, nem tudo está bem no sentido de que eles não são equivalentes. Esse é o ponto. (ver minha atualização)
Molho para salada
11
Então você precisa descrever a situação real que está enfrentando, não um exemplo de brinquedo diferente. As perguntas XY geralmente levam a respostas que não abordam o que o solicitante realmente precisa saber.
gung - Restabelece Monica
Eu adicionei uma atualização mostrando que as comparações com o modelo somente de interceptação são inconsistentes entre as abordagens devido ao número de observações ser diferente. Minha pergunta é qual das duas abordagens aparentemente equivalentes é "correta".
Molho para salada
11
@Saladdressing, não sigo sua atualização. Onde g_bern0, por exemplo, é definido? Se a ideia é que deveria ser um modelo somente de interceptação, não vejo problema. Os números reais do BIC não importam. É apenas uma medida relativa. Nos dois casos (ie, g_bern0& g_binom0), a diferença é negativa, implicando a mesma preferência (para o modelo somente de interceptação). Não há inconsistência.
gung - Restabelece Monica
A diferença no BIC entre as duas comparações foi de uma ordem de magnitude, que é enorme em um valor logarítmico, como o BIC. E esses devem ser os mesmos modelos. Mas de qualquer maneira, eu incrementei a semente. Agora, há uma inconsistência de sinal entre as duas comparações da BIC. Não consigo ver como "os dois formatos são equivalentes - eles contêm as mesmas informações e, na maioria das vezes, parecem diferentes do lado de fora" poderiam ser verdade, dado que um BIC é -18 e o outro é 1. Essa é uma diferença impressionante .
Molho para salada
3

Pergunta interessante! Chegando a isso de uma configuração aplicada, acho que você precisa se lembrar de que tanto o BIC quanto o AIC são medidas de ajuste relativo do modelo.

Em outras palavras, essas medidas não dizem muito quando você as examina para um único modelo, mas podem ajudá-lo a selecionar um modelo apropriado entre um conjunto de modelos concorrentes. Em particular:

  1. Se seu objetivo é encontrar o melhor entre os modelos concorrentes para a previsão da variável de resultado, selecione o modelo com o menor valor de AIC;
  2. Se seu objetivo é encontrar o 'melhor' entre os modelos concorrentes para entender e descrever os efeitos das variáveis ​​preditivas incluídas no modelo na variável de resultado, selecione o modelo com o menor valor de BIC.

Ao definir seu conjunto de modelos concorrentes, você deve garantir que os modelos sigam a mesma estrutura conceitual. Assim, você pode comparar vários modelos de regressão logística binomial ou vários modelos de logística binária, mas não uma mistura de ambos. (É importante comparar gostos com gostos, caso contrário você não saberá se um modelo venceu a competição com base em seus próprios méritos ou simplesmente porque alterou o procedimento de especificação / ajuste do modelo.)

Nesta perspectiva, a única coisa que importa é que R seja consistente ao calcular o AIC e o BIC em modelos do mesmo tipo (por exemplo, modelos de regressão logística binomial).

Apenas para esclarecer: g_bern é um modelo de regressão logística binária, enquanto g_binom é um modelo de regressão logística binomial. Embora ambos modelem a probabilidade de sucesso em uma tentativa, você não misturaria variações desses modelos ao definir seu conjunto de modelos concorrentes (pelos motivos explicados acima e também cobertos pelo @gung).

Isabella Ghement
fonte
Por favor, veja minha atualização. Apenas para esclarecer: não estou tentando comparar g_bern com g_binom. Este é um exemplo de brinquedo.
Molho para salada