Exercícios — Regressão Linear e Máxima Verossimilhança

Aula 5 — Por Que Mínimos Quadrados É Máxima Verossimilhança Sob Ruído Gaussiano

Autor

Marcos M. Raimundo — Instituto de Computação, UNICAMP

Questões discursivas

Considere um ajuste de regressão linear simples, prevendo o preço mediano de imóveis (MedHouseVal) a partir da renda mediana da região (MedInc), no dataset California Housing — um modelo com apenas dois “vetores base”: $\varphi_1$ (a coluna de intercepto, um vetor de uns) e $\varphi_2$ (a coluna do atributo MedInc). Duas fontes clássicas descrevem a mesma ideia geométrica — a solução de mínimos quadrados como projeção ortogonal — usando notações ligeiramente diferentes: PRML escreve o subespaço gerado como $\mathcal S=\text{span}(\varphi_1,\ldots,\varphi_M)$; ESL descreve o mesmo objeto como o hiperplano gerado por $x_0,\ldots,x_p$. Usando o ajuste MedInc$\to$MedHouseVal acima como exemplo concreto, explique o que os vetores $\varphi_1,\varphi_2$ representam nesse ajuste específico, e por que $\hat y$ precisa estar dentro do subespaço gerado por eles.
O estimador de máxima verossimilhança $\hat\sigma^2_{\text{MLE}}= \text{RSS}/N$ é viesado para baixo — o mesmo padrão de viés que aparece em qualquer estimador de variância obtido por máxima verossimilhança que primeiro estima uma média (ou, num modelo de regressão, um hiperplano) a partir dos mesmos dados usados depois para medir a dispersão em torno dela. Considere um ajuste de regressão linear ao California Housing (MedInc$\to$MedHouseVal, com $2$ parâmetros: intercepto e inclinação) com $N=16\,640$ observações e $\text{RSS}\approx11\,710{,}64$ no ponto de ajuste ótimo. Calcule $\hat\sigma^2_{\text{MLE}}$ e o estimador não-viesado $\text{RSS}/(N-2)$, compare os dois valores numericamente, e explique por que a diferença relativa é tão pequena neste caso específico — e em que situação (em termos de $N$ e do número de parâmetros) essa diferença deixaria de ser desprezível.
Num ajuste de regressão linear do preço de imóveis (MedHouseVal) a partir da renda mediana (MedInc) no dataset California Housing, a variância dos resíduos cresce sistematicamente com MedInc — de $0{,}486$ no 1º quartil de renda a $0{,}849$ no 3º quartil —, uma violação real da suposição de homocedasticidade (variância do ruído constante, igual para todas as observações). Considerando a construção usual da verossimilhança sob ruído gaussiano homocedástico (uma única variância $\sigma^2$ comum a todas as observações — ver o bloco “Construção da Verossimilhança Conjunta” abaixo para a fórmula geral), explique exatamente qual passo da derivação deixaria de ser válido se $\sigma^2$ dependesse de $x$, e proponha (sem precisar resolver completamente) como a log-verossimilhança mudaria de forma.

Questões de Verdadeiro/Falso

Cada bloco de 4 itens trata do mesmo tema. A questão só é considerada correta se todos os 4 itens forem julgados corretamente (deixar em branco tem penalidade de 20% da nota da questão).

Nas questões desta seção, considere um ajuste de regressão linear $Y\mid X=x\sim\mathcal N(\beta^Tx,\sigma^2)$ (ruído gaussiano homocedástico, observações independentes), com solução de mínimos quadrados $\hat\beta=(X^TX)^{-1}X^Ty$ — que é também a projeção ortogonal de $y$ sobre o subespaço $\mathcal S$ gerado pelas colunas de $X$.

Teste 1 — Geometria da Projeção e Resíduos

□ Se o vetor observado $y$ já estivesse exatamente dentro do subespaço gerado pelas colunas de $X$ (por exemplo, se os dados não tivessem ruído algum), a projeção $\hat y$ coincidiria com o próprio $y$, e o resíduo seria o vetor nulo.
□ Se o número de colunas de $X$ (incluindo o intercepto) fosse igual ao número de observações $N$, e essas colunas fossem linearmente independentes, o subespaço $\mathcal S$ ocuparia todo o espaço $\mathbb R^N$, e o ajuste por mínimos quadrados reproduziria $y$ exatamente, com resíduo nulo em todos os pontos.
□ Num algoritmo de recomendação que aproxima as avaliações de um usuário como combinação linear de “perfis latentes” pré-definidos, a mesma lógica de projeção ortogonal sobre o subespaço gerado pelos perfis garante que a aproximação encontrada terá erro zero sempre que houver pelo menos dois perfis latentes disponíveis.
□ Como o resíduo $y-\hat y$ é ortogonal a cada coluna de $X$ individualmente, ele é necessariamente ortogonal a qualquer vetor fora do subespaço $\mathcal S$ também.

Teste 2 — RSS, Equações Normais e Posto de X

□ Se todas as colunas de atributos de $X$ (exceto o intercepto) fossem multiplicadas por uma mesma constante $c\ne 0$ antes do ajuste, o vetor $\hat\beta$ mudaria (os coeficientes correspondentes seriam divididos por $c$), mas a previsão $\hat y=X\hat\beta$ e o RSS mínimo permaneceriam exatamente os mesmos.
□ No limite em que $X$ tem uma única coluna (só o intercepto, sem nenhum atributo), a equação normal $X^T(y-X\hat\beta)=0$ se reduz a dizer que $\hat\beta_0$ é exatamente a média amostral $\bar y$.
□ Num ajuste de regressão linear que busca prever o preço mediano de um imóvel (dataset California Housing) a partir da renda mediana da região, a mesma equação normal determinaria os coeficientes independentemente de qual atributo do dataset for usado como preditor — e a extrapolação da reta ajustada para valores de renda muito fora do intervalo observado no treino teria a mesma confiabilidade estatística que dentro do intervalo observado.
□ Como $\text{RSS}(\beta)$ é uma função quadrática (e portanto convexa) de $\beta$, ela tem sempre um único mínimo global estritamente, mesmo quando as colunas de $X$ são linearmente dependentes ($X^TX$ singular).

Teste 3 — Premissas do Modelo Gaussiano Homocedástico

□ Se assumíssemos que os erros $\epsilon_n$ são gaussianos mas NÃO independentes entre observações (por exemplo, erros correlacionados entre imóveis vizinhos geograficamente), a verossimilhança conjunta deixaria de ser simplesmente o produto das densidades marginais de cada observação.
□ Se a variância do ruído $\sigma^2$ fosse igual para todas as observações mas extremamente grande (tendendo a infinito), o modelo linear perderia praticamente toda sua capacidade preditiva, mesmo que $\beta$ estivesse corretamente especificado.
□ Num modelo que prevê o tempo de resposta de um servidor a partir da carga de requisições, supor ruído gaussiano homocedástico seria tão razoável quanto no caso de preços de imóveis, desde que a variabilidade do tempo de resposta realmente não dependa sistematicamente do nível de carga — e, uma vez feita essa suposição, ela nunca precisaria ser verificada empiricamente depois do ajuste, pois o próprio ajuste por mínimos quadrados a confirma automaticamente.
□ Como um modelo de regressão linear assume $Y\mid X\sim\mathcal N(\beta^TX,\sigma^2)$, isso implica que a distribuição marginal (não condicional) de $Y$ também deve ser gaussiana, para o modelo fazer sentido.

Teste 4 — Construção da Verossimilhança Conjunta

Sob o modelo $Y\mid X=x\sim\mathcal N(\beta^Tx,\sigma^2)$ com observações independentes, a log-verossimilhança dos $N$ dados é $\ell(\beta,\sigma^2)=-\frac N2\ln(2\pi)-N\ln\sigma-\frac{\text{RSS}(\beta)}{2\sigma^2}$, com $\text{RSS}(\beta)=\sum_n(y_n-\beta^Tx_n)^2$.

□ Se, em vez de tomar o logaritmo da verossimilhança, trabalhássemos diretamente com o produto de densidades $\prod_n\mathcal N(y_n\mid\beta^Tx_n,\sigma^2)$, maximizar esse produto em relação a $\beta$ ainda levaria exatamente ao mesmo $\hat\beta$ obtido via log-verossimilhança, porque o logaritmo é uma função estritamente crescente.
□ No caso extremo de $N=1$ (uma única observação), a “log-verossimilhança” ainda pode ser escrita pela mesma fórmula geral, mas maximizá-la em relação a $\beta$ (com mais de um parâmetro livre) não determina um único $\hat\beta$, porque há infinitas retas passando exatamente por um ponto.
□ Num modelo que assume que o tempo entre chegadas de clientes numa fila segue uma distribuição exponencial (não gaussiana), a mesma lógica de multiplicar densidades individuais (assumindo observações independentes) e tomar o logaritmo se aplicaria — ainda que a forma final da função a maximizar fosse diferente da gaussiana acima, o valor numérico do máximo obtido seria sempre diretamente comparável ao valor máximo de uma log-verossimilhança gaussiana ajustada a um problema de regressão, permitindo escolher qual dos dois modelos descreve melhor seus respectivos dados só comparando os números brutos.
□ Como a log-verossimilhança acima é $\ell(\beta,\sigma^2)=-\frac N2\ln(2\pi)-N\ln\sigma-\frac{\text{RSS}(\beta)}{2\sigma^2}$, maximizar $\ell$ em relação a $\beta$ e maximizar $\ell$ em relação a $\sigma^2$ são operações que precisam ser feitas simultaneamente, nunca uma depois da outra, porque os dois parâmetros aparecem na mesma expressão.

Teste 5 — O Teorema Central: OLS = MLE

□ Se a variância $\sigma^2$ fosse conhecida e fixa de antemão (não estimada), o valor de $\hat\beta$ que maximiza a log-verossimilhança ainda seria exatamente o mesmo $\hat\beta_{\text{OLS}}$, porque $\sigma^2$ entra na log-verossimilhança apenas como um fator multiplicativo positivo do termo que depende de $\beta$.
□ Se, por hipótese, $\text{RSS}(\hat\beta)=0$ para o $\hat\beta$ que minimiza a soma de quadrados (ajuste perfeito, sem nenhum resíduo), o estimador $\hat\sigma^2_{\text{MLE}}=\text{RSS}(\hat\beta)/N$ seria exatamente zero, e a log-verossimilhança avaliada nesse ponto tenderia a $+\infty$ conforme $\sigma\to0$ nessa mesma expressão, revelando uma degenerescência do modelo gaussiano nesse caso extremo.
□ Se, em outro domínio qualquer, alguém provar que maximizar uma verossimilhança assumida é algebricamente idêntico a minimizar uma soma de quadrados, essa prova por si só já implica que o ruído daquele modelo específico foi assumido gaussiano, porque é a densidade gaussiana (e nenhuma outra densidade comum) que produz esse termo quadrático especificamente na log-verossimilhança.
□ Como $\hat\beta_{\text{MLE}}=\hat\beta_{\text{OLS}}$ sob um modelo de regressão linear com ruído gaussiano homocedástico, conclui-se que máxima verossimilhança e mínimos quadrados são sempre a mesma coisa, para qualquer modelo estatístico, não apenas para regressão linear sob ruído gaussiano.

Teste 6 — Viés do MLE de $\sigma^2$

□ Se, em vez de dividir por $N$, dividíssemos $\text{RSS}(\hat\beta)$ por $N-2$ (número de observações menos o número de parâmetros do modelo — por exemplo, $2$ num ajuste de regressão linear simples, com intercepto e uma única inclinação), o estimador resultante teria valor esperado exatamente igual a $\sigma^2$ verdadeiro — mas deixaria de ser o estimador de máxima verossimilhança.
□ Conforme $N\to\infty$ (mantendo o número de parâmetros de $\beta$ fixo), a diferença relativa entre $\hat\sigma^2_{\text{MLE}}=\text{RSS}/N$ e o estimador não-viesado $\text{RSS}/(N-2)$ tende a zero, porque $N/(N-2)\to1$.
□ O mesmo padrão de viés — dividir por $N$ em vez de $N$ menos o número de parâmetros estimados — reaparece em qualquer estimador de variância obtido por máxima verossimilhança que primeiro estime uma média (ou, aqui, um hiperplano) a partir dos mesmos dados usados para medir a dispersão em torno dela; e por isso, na prática, é sempre preferível reportar o estimador não-viesado em vez do MLE, em qualquer aplicação, sem exceção.
□ Como $\hat\sigma^2_{\text{MLE}}$ é viesado para baixo, ele é necessariamente um estimador pior do que o estimador não-viesado $\text{RSS}/(N-2)$ em qualquer critério razoável de qualidade, incluindo o erro quadrático médio do próprio estimador de $\sigma^2$.

Teste 7 — Verificação Numérica e Otimização

□ Se uma otimização numérica (por exemplo, via algoritmo Nelder-Mead) da log-verossimilhança de um modelo de regressão linear gaussiano tivesse sido inicializada num ponto de partida muito distante de $\hat\beta_{\text{OLS}}$ (por exemplo, $\beta_0=100,\beta_1=-50$, quando o ótimo verdadeiro está perto da origem), o algoritmo ainda convergiria ao mesmo ponto ótimo, porque a log-verossimilhança gaussiana em função de $\beta$ é uma função côncava sem outros máximos locais.
□ Se uma grade (grid search) usada para visualizar a superfície de log-verossimilhança de um modelo de regressão linear gaussiano tivesse uma resolução extremamente grosseira (poucos pontos, muito espaçados), o pico da grade poderia ficar visivelmente deslocado do verdadeiro $\hat\beta_{\text{OLS}}$, mesmo que a forma fechada continuasse exata — e, por essa razão, sempre que uma forma fechada exata existir, o grid search se torna uma ferramenta cientificamente inútil, sem nenhum valor pedagógico ou de verificação.
□ Em problemas de otimização com muitos parâmetros (por exemplo, uma rede neural com milhões de pesos), a mesma ideia de verificar uma solução fechada contra uma busca numérica deixaria de ser prática simplesmente por causa do volume de contas — não porque o princípio de comparação deixe de fazer sentido.
□ Como uma otimização numérica confirmou o resultado da forma fechada num exemplo específico de ajuste, isso prova matematicamente (não apenas ilustra empiricamente) que $\hat\beta_{\text{MLE}}=\hat\beta_{\text{OLS}}$ para qualquer conjunto de dados sob o modelo gaussiano.

Teste 8 — Limitações: Heterocedasticidade, Censura e Não-Linearidade

O alvo MedHouseVal do dataset California Housing é artificialmente censurado: todo imóvel cujo preço mediano real seria maior que $\$500\,001$ (a unidade do dataset é centenas de milhares de dólares, então o valor aparece como $5{,}00001$) é registrado com exatamente esse valor de teto — cerca de $4{,}49\%$ das observações são afetadas. Além disso, num ajuste de regressão linear de MedHouseVal a partir de MedInc nesse dataset, a variância dos resíduos cresce sistematicamente com MedInc (heterocedasticidade real, não apenas hipotética).

□ Se MedHouseVal não fosse artificialmente censurado em $5{,}00001$ (ou seja, se os $4{,}49\%$ de observações no teto tivessem seus valores reais, possivelmente maiores), o $\hat\beta_1$ estimado para MedInc provavelmente mudaria, porque essas observações censuradas concentram-se desproporcionalmente entre os imóveis de renda mais alta.
□ No extremo em que a censura afetasse 100% das observações (todos os valores exibidos como exatamente $5{,}00001$), a variância amostral de MedHouseVal seria zero, e o ajuste de mínimos quadrados encontraria $\hat\beta_1=0$.
□ Num estudo que mede o tempo de sobrevivência de pacientes após um tratamento, mas o estudo termina antes de alguns pacientes falecerem (censura à direita, análoga à censura de MedHouseVal no teto), a mesma distorção de ajustar mínimos quadrados diretamente sobre os tempos observados (ignorando a censura) se aplicaria, subestimando o verdadeiro tempo médio de sobrevivência — mas esse problema desapareceria automaticamente com uma amostra suficientemente grande, já que o viés de censura, ao contrário do viés de estimadores de variância, diminui conforme $N$ aumenta.
□ Como a variância dos resíduos de um ajuste de regressão linear cresce com MedInc (heterocedasticidade real, já confirmada empiricamente nesse ajuste), o modelo linear ajustado por mínimos quadrados se torna inútil para prever MedHouseVal, e nenhuma previsão feita a partir dele deveria ser usada.