| This post was kindly contributed by Ridículas - go there to comment and to read the full post. |
Quero falar de um erro comum de interpretação de resultados em análise de regressão através de um exemplo. Considere que você tem pessoas com diferentes pesos e diferentes alturas. Facilmente você aceita que uma pessoa mais alta tem maior peso, certo? Ou seja, existe uma correlação positiva entre peso e altura. Pois bem, vamos simular observações desse experimento e em seguida ajustar uma regressão linear simples para relação peso-altura.
require(MASS)
set.seed(12345)
da <- as.data.frame(mvrnorm(80, c(altura=30, peso=60),
matrix(c(2,1.9,1.9,3), 2,2)))
plot(peso~altura, da)
m0 <- lm(peso~altura, data=da) # ajusta a reta
summary(m0) # estimativa dos parâmetros
abline(m0) # adiciona uma reta ao gráfico
vcov(m0) # covariância das estimativas
Aqui está o ponto que eu quero comentar: a interpretação da matriz de covariância. Perceba que a covariância foi negativa. Já vi gente interpretando isso da seguinte forma: espera-se que pessoas mais altas tenham menor peso. Duplamente errado. Primeiro nos sabemos por experiência que a correlação de peso e altura é positiva. Segundo, a matriz de covariância se refere as estimativas dos parâmetros e não as variáveis envolvidas. Nunca a interprete dessa forma.
Então como interpretar? Bem, a matriz de covariância das estimativas, é um reflexo da função objetivo ao redor da solução. A função objetivo nesse caso é minimizar a soma de quadrados (SQ). Então se eu aumento o valor do parâmetro b0, o parâmetro b1 tem que diminuir para compensar, ou seja, para diminuir a SQ. Ocorre um efeito compensatório aqui. Além do mais, esse efeito pode ser eliminado facilmente aplicando uma reparametrização, como por exemplo, centrar a covariável na média.
da$alturac <- with(da, altura-mean(altura)) # centrando a covariável par(mfrow=c(1,2)) plot(peso~altura, da) abline(m0) plot(peso~alturac, da) m1 <- lm(peso~alturac, data=da) # ajusta a reta summary(m1) # estimativa dos parâmetros abline(m1) # adiciona uma reta ao gráfico vcov(m1) # covariância das estimativas require(ellipse) par(mfrow=c(1,2)) # regiões de confiança conjunta plot(ellipse(vcov(m0)), type="l") plot(ellipse(vcov(m1)), type="l") cov2cor(vcov(m0)) cov2cor(vcov(m1))
Perceba que os resultados são os mesmos em termos de estatísticas de teste e medidas de ajuste. Isso é esperado por eu só fiz uma traslação da altura. Mas o importante é que agora a matriz de covariância tem covariâncias praticamente nulas, que é resultado da translação. O intercepto é então o valor esperado para alturac igual a zero (centro dos dados). Veja porque essa covariância é zero: se eu aumentar b0, não adiante eu alterar b1 que não vai diminuir a SQ, e vice-versa, porque agora eles são ortogonais. Ortogonalidade entre parâmetros é uma característica desejável pois permite você inferir sobre um deles sem considerar o outro. Além disso, tem vantagens do ponto de vista de estimação por métodos numéricos de otimização. Em outras palavras, pegando conceitos de probabilidade, se a distribuição amostral de duas variáveis aleatórias é normal bivariada com covariância nula, a distribuição condicional de A|B é igual a distribuição marginal de A pela independência. Reforçando, eu não preciso conhecer valores de B para descrever A. Considerando tudo que foi argumentado, é sempre preferível que você adote o modelo que apresente menor covariância entre os parâmetros. Até a próxima ridícula.