Regressão Logística e Modelos Lineares Generalizados (GLMs)

Aula 6 — Classificação Paramétrica, Entropia Cruzada e a Família Exponencial

Autor

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

Data de Publicação

21 de setembro de 2026

Slides Lista de aulas

1 Revisão e Introdução

Na Aula 5, estabelecemos a base paramétrica do aprendizado supervisionado contínuo: sob a suposição de ruído aditivo gaussiano homocedástico \(\epsilon \sim \mathcal{N}(0, \sigma^2)\), o modelo estatístico assume que o alvo \(Y\) condicionado às características \(\mathbf{x}\) segue uma distribuição normal com média linear \(Y \mid X=\mathbf{x} \sim \mathcal{N}(\mathbf{w}^T\mathbf{x}, \sigma^2)\). Provamos analiticamente que maximizar a log-verossimilhança conjunta sob independência equivale rigorosamente a minimizar a soma dos quadrados dos resíduos (RSS). O estimador de Mínimos Quadrados Ordinários (OLS) é, portanto, o estimador de Máxima Verossimilhança (MLE) do hiperplano ótimo: \(\hat{\mathbf{w}} = (X^T X)^{-1}X^T \mathbf{y}\).

Hoje avançamos da regressão contínua para a classificação paramétrica. O objetivo agora é prever um rótulo categórico binário \(y \in \{0, 1\}\). No entanto, em vez de aplicar um corte rígido determinístico ou utilizar modelos não-paramétricos como árvores de decisão (Aula 3), nosso propósito é estimar a probabilidade a posteriori de classe \(p(C_1 \mid \mathbf{x}) = \mathbb{E}[Y \mid \mathbf{x}]\) mantendo a simplicidade e a interpretabilidade de uma combinação linear dos atributos.

O obstáculo imediato é conceitual: um hiperplano linear \(\mathbf{w}^T\mathbf{x}\) projeta valores ao longo de toda a reta real \((-\infty, +\infty)\), enquanto probabilidades requerem estritamente o intervalo unitário \([0, 1]\). Tentar ajustar uma reta de regressão linear diretamente aos rótulos \(0\) e \(1\) gera probabilidades espúrias (negativas ou maiores que \(1\)) e sofre de sensibilidade extrema a observações distantes da fronteira. A solução estatisticamente principled consiste em trocar a distribuição de ruído gaussiana pela distribuição de Bernoulli, e mapear a combinação linear para o intervalo de probabilidades através de uma função de ativação estritamente monotônica: a função sigmoide logística.

Hoje responderemos a quatro perguntas centrais:

  1. Por que a regressão linear falha ao ser aplicada a problemas de classificação com rótulos binários?
  2. Como a distribuição de Bernoulli induz a perda de Entropia Cruzada a partir do princípio da Máxima Verossimilhança?
  3. Por que a regressão logística não possui solução fechada, mas pode ser resolvida de maneira exata e eficiente por Newton-Raphson via Mínimos Quadrados Ponderados Iterativos (IRLS)?
  4. Qual é o arcabouço unificador dos Modelos Lineares Generalizados (GLMs) que conecta regressão linear, regressão logística e redes neurais numa mesma família estatística?

1.0.1 O Problema Motivador: Tentando Usar Mínimos Quadrados em Classificação

Para constatar concretamente o fracasso de Mínimos Quadrados em classificação binária, examinamos o dataset real Adult Census Income. Nosso alvo \(y\) indica se uma pessoa ganha mais de 50 mil dólares ao ano (\(y=1\)) ou não (\(y=0\)). Se ajustarmos uma regressão linear comum \(y = w_0 + w_1 \cdot \text{idade}\), obtemos a reta da figura abaixo.

A inspeção visual revela duas falhas capitais na reta linear:

  1. Violação das leis de probabilidade: Para indivíduos muito jovens (por exemplo, abaixo de 10 anos na extrapolação matemática da reta), a reta prevê valores negativos. Em amostras com grande dispersão, valores superiores a \(1{,}0\) ocorrem rotineiramente.
  2. Penalização assimétrica irracional: Um indivíduo de 75 anos que ganha \(>50\text{K}\) tem rótulo \(y=1\). A reta OLS prevê cerca de \(0{,}51\). No entanto, se tivéssemos um ponto extremo com \(y=1\) e valor previsto \(1{,}30\), o resíduo seria \(1 - 1{,}30 = -0{,}30\), e o quadrado desse erro penalizaria a reta, forçando-a a mudar de inclinação simplesmente porque um ponto foi “correto demais”. Mínimos quadrados penaliza previsões com erro quadrático mesmo quando elas apontam com extrema convicção para o lado correto da fronteira.
DicaPausa Ativa 1 — Regressão Linear Aplicada a Rótulos Binários

A reta de mínimos quadrados ajustada sobre dados \(y \in \{0, 1\}\) é frequentemente chamada de “Modelo de Probabilidade Linear” na econometria. Pense sobre os pressupostos estatísticos do erro \(e_i = y_i - \hat{y}_i\).

  • □ Se a variável de resposta \(y_i\) assume exclusivamente os valores \(0\) ou \(1\), os resíduos da regressão linear \(e_i = y_i - \mathbf{w}^T\mathbf{x}_i\) não podem ter distribuição gaussiana, pois assumem apenas dois valores condicionais discretos para cada \(\mathbf{x}_i\).
  • □ Como a variância da distribuição de Bernoulli é dada por \(\mu(1 - \mu)\), a variância do termo de ruído em dados binários depende diretamente de \(\mathbf{x}_i\), violando a suposição de homocedasticidade exigida pelo Teorema de Gauss-Markov.
  • □ A inclusão de um atributo extremo perfeitamente classificado (um indivíduo de idade avançada que ganha muito acima de 50 mil) não afeta a inclinação da reta de regressão linear, pois seu resíduo é interpretado pelo modelo como acerto perfeito.
  • □ Mapear os valores contínuos da reta linear através de uma função que simplesmente trunca previsões menores que \(0\) em \(0\) e maiores que \(1\) em \(1\) recupera a equivalência com a Máxima Verossimilhança de Bernoulli.

2 Intuição — A Transformação Sigmoide e a Razão de Chances

Para construir um modelo estatístico consistente, precisamos manter duas propriedades desejáveis:

  1. A simplicidade de representação: Queremos continuar operando com uma combinação linear de atributos \(a(\mathbf{x}) = \mathbf{w}^T\mathbf{x} = w_0 + w_1 x_1 + \dots + w_D x_D\), chamada de ativação linear;
  2. A validade probabilística: Queremos que o modelo produza uma probabilidade \(\mu(\mathbf{x}) = p(C_1 \mid \mathbf{x}) \in (0, 1)\).

A conexão entre a reta real e o intervalo \((0, 1)\) emerge naturalmente da ideia de chances (odds). Em probabilidade, as chances a favor de um evento são definidas pela razão entre a probabilidade de ocorrência e a probabilidade de não-ocorrência: \[\text{odds} = \frac{p}{1 - p}\]

Se \(p = 0{,}8\), as chances são de \(4\) para \(1\) (\(\text{odds} = 4\)). Se \(p = 0{,}2\), as chances são de \(0{,}25\). Enquanto \(p\) varia no intervalo \((0, 1)\), as chances variam em \((0, +\infty)\). Se tomarmos o logaritmo natural das chances — quantidade conhecida em estatística como logit: \[a = \text{logit}(p) = \ln \left( \frac{p}{1 - p} \right)\] constatamos que o logit varia livremente em todo o domínio real \((-\infty, +\infty)\). Se \(p \to 0\), \(\ln(p/(1-p)) \to -\infty\); se \(p \to 1\), \(\ln(p/(1-p)) \to +\infty\); e se \(p = 0{,}5\), as chances são unitárias e o logit é exatamente \(\ln(1) = 0\).

A hipótese fundacional da Regressão Logística é postular que o logaritmo das chances é uma função linear das características: \[\ln \left( \frac{\mu(\mathbf{x})}{1 - \mu(\mathbf{x})} \right) = \mathbf{w}^T\mathbf{x}\]

Para obter a probabilidade \(\mu(\mathbf{x})\) a partir da ativação linear \(a = \mathbf{w}^T\mathbf{x}\), basta inverter a função logit: \[\frac{\mu}{1 - \mu} = e^a \implies \mu = e^a(1 - \mu) \implies \mu(1 + e^a) = e^a \implies \mu = \frac{e^a}{1 + e^a} = \frac{1}{1 + e^{-a}}\]

Essa função inversa é a célebre função sigmoide logística (ou sigmoide padrão): \[\sigma(a) \equiv \frac{1}{1 + e^{-a}}\]

2.0.1 Propriedades Analíticas da Sigmoide

A função sigmoide possui duas propriedades analíticas cruciais que simplificam enormemente a matemática da regressão logística e de redes neurais:

  1. Simetria: \[\sigma(-a) = \frac{1}{1 + e^a} = \frac{e^{-a}}{e^{-a}(1 + e^a)} = \frac{e^{-a}}{1 + e^{-a}} = 1 - \frac{1}{1 + e^{-a}} = 1 - \sigma(a)\] Isso garante que \(p(C_0 \mid \mathbf{x}) = 1 - \sigma(\mathbf{w}^T\mathbf{x}) = \sigma(-\mathbf{w}^T\mathbf{x})\), conferindo tratamento estritamente simétrico a ambas as classes.

  2. Derivada Notável: A derivada de \(\sigma(a)\) em relação a \(a\) expressa-se elegantemente em termos da própria função: \[\frac{d\sigma}{da} = \frac{d}{da}(1 + e^{-a})^{-1} = -(1 + e^{-a})^{-2}(-e^{-a}) = \frac{1}{1 + e^{-a}} \cdot \frac{e^{-a}}{1 + e^{-a}} = \sigma(a)(1 - \sigma(a))\] Essa identidade é a responsável direta pelo cancelamento algébrico que recupera a elegância do gradiente no ajuste por máxima verossimilhança.

3 A Verossimilhança de Bernoulli e a Entropia Cruzada

Tendo estabelecido que a saída do modelo \(\mu_n = \sigma(\mathbf{w}^T\mathbf{x}_n)\) representa a probabilidade condicional \(p(y_n = 1 \mid \mathbf{x}_n, \mathbf{w})\), procedemos agora ao desenvolvimento formal do estimador pelo princípio da Máxima Verossimilhança, seguindo a formulação do PRML (§4.3.2).

3.0.1 Premissas Formais do Modelo

  1. Rótulos Binários: As variáveis de resposta observadas são binárias: \(y_n \in \{0, 1\}\) para \(n = 1, \dots, N\).
  2. Independência Amostral: Condicionado aos atributos \(\mathbf{x}_n\) e aos pesos \(\mathbf{w}\), as observações são estatisticamente independentes.
  3. Distribuição Condicional de Bernoulli: Para cada observação individual, a probabilidade condicional segue estritamente a lei de Bernoulli: \[Y_n \mid \mathbf{x}_n, \mathbf{w} \sim \text{Bernoulli}(\mu_n), \quad \text{onde } \mu_n = \sigma(\mathbf{w}^T\mathbf{x}_n)\]

Sob a distribuição de Bernoulli, a função de massa de probabilidade de uma única observação \(y_n\) pode ser escrita de forma algébrica compacta e unificada: \[p(y_n \mid \mathbf{w}) = \mu_n^{y_n} (1 - \mu_n)^{1 - y_n}\] Observe a elegância dessa expressão: quando \(y_n = 1\), o segundo termo se anula \((1-\mu_n)^0 = 1\) e obtemos \(\mu_n\); quando \(y_n = 0\), o primeiro termo se anula e resta \(1 - \mu_n\).

3.0.2 Construção da Verossimilhança Conjunta

Como as observações são assumidas condicionalmente independentes, a função de verossimilhança conjunta de todo o conjunto de dados \(\mathcal{D} = \{(\mathbf{x}_n, y_n)\}_{n=1}^N\) é o produtório das probabilidades individuais: \[p(\mathbf{y} \mid X, \mathbf{w}) = \prod_{n=1}^N p(y_n \mid \mathbf{w}) = \prod_{n=1}^N \mu_n^{y_n} (1 - \mu_n)^{1 - y_n}\]

Aplicando o logaritmo natural para obter a log-verossimilhança \(\ln L(\mathbf{w})\): \[\ln p(\mathbf{y} \mid X, \mathbf{w}) = \sum_{n=1}^N \left[ y_n \ln \mu_n + (1 - y_n) \ln(1 - \mu_n) \right]\]

3.0.3 A Perda de Entropia Cruzada Negativa

Em aprendizado de máquina e otimização, convenciona-se formular problemas computacionais como a minimização de uma função de perda (ou erro) \(E(\mathbf{w})\). Tomando o negativo da log-verossimilhança (NLL — Negative Log-Likelihood), definimos a função de erro de Entropia Cruzada Binária (Binary Cross-Entropy): \[E(\mathbf{w}) = -\ln p(\mathbf{y} \mid X, \mathbf{w}) = -\sum_{n=1}^N \left[ y_n \ln \mu_n + (1 - y_n) \ln(1 - \mu_n) \right]\]

Esta função possui uma interpretação imediata à luz da Teoria da Informação (estudada na Aula 3): ela mede a divergência de Kullback-Leibler entre a distribuição empírica dos dados \(p(y_n)\) e a distribuição prevista pelo modelo \(q(y_n) = \mu_n\). Se o modelo prevê \(\mu_n = 0{,}99\) para um ponto com \(y_n = 1\), a perda é ínfima (\(-\ln(0{,}99) \approx 0{,}01\)). Se o modelo comete um erro catastrófico prevendo \(\mu_n \to 0\) para um ponto cujo rótulo é \(y_n = 1\), a perda explode para \(-\ln(\mu_n) \to +\infty\).

DicaPausa Ativa 2 — A Verossimilhança de Bernoulli e a Entropia Cruzada

Considere a função de perda \(E(\mathbf{w}) = -\sum_{n=1}^N [y_n \ln \mu_n + (1 - y_n)\ln(1 - \mu_n)]\), onde \(\mu_n = \sigma(\mathbf{w}^T\mathbf{x}_n)\).

  • □ Se multiplicarmos todos os rótulos \(y_n\) por uma constante arbitrária \(c > 0\), a perda de entropia cruzada continua sendo uma log-verossimilhança estritamente válida para uma distribuição de probabilidade normalizada.
  • □ No caso de uma amostra com desbalanceamento severo (exemplo: 99% dos dados pertencem à classe \(0\)), o preditor ingênuo \(\mu_n = 0{,}01\) constante para todos os pontos atinge perda finita e baixa, mas reflete apenas a entropia da priori de classe.
  • □ Em um sistema de detecção de fraudes financeiras onde falsos negativos custam 50 vezes mais que falsos positivos, a minimização direta da perda de entropia cruzada clássica ajusta automaticamente os pesos para refletir essa assimetria de custo de decisão.
  • □ A perda de entropia cruzada binária é estritamente equivalente à soma dos erros quadráticos quando \(\mu_n\) se aproxima de \(0{,}5\).

4 Otimização: Gradiente, Convexidade e IRLS

Na Aula 5, encontramos o ponto de mínimo igualando o gradiente de \(\text{RSS}(\mathbf{w})\) a zero, obtendo uma equação linear fechada. Veremos agora por que isso não é possível em Regressão Logística, como calcular o gradiente de forma analítica e por que a geometria da função de erro garante convergência robusta.

4.0.1 Derivação do Gradiente Passo a Passo

Queremos calcular o vetor gradiente \(\nabla E(\mathbf{w}) = \frac{\partial E}{\partial \mathbf{w}}\). A função de erro para um único ponto \(n\) é: \[E_n(\mathbf{w}) = -\left[ y_n \ln \mu_n + (1 - y_n) \ln(1 - \mu_n) \right]\] onde \(\mu_n = \sigma(a_n)\) e a ativação linear é \(a_n = \mathbf{w}^T\mathbf{x}_n\). Pela regra da cadeia multivariada: \[\frac{\partial E_n}{\partial \mathbf{w}} = \frac{\partial E_n}{\partial \mu_n} \cdot \frac{\partial \mu_n}{\partial a_n} \cdot \frac{\partial a_n}{\partial \mathbf{w}}\]

Calculamos cada termo separadamente:

  1. Derivada da perda em relação à probabilidade \(\mu_n\): \[\frac{\partial E_n}{\partial \mu_n} = -\left[ \frac{y_n}{\mu_n} - \frac{1 - y_n}{1 - \mu_n} \right] = -\frac{y_n(1 - \mu_n) - (1 - y_n)\mu_n}{\mu_n(1 - \mu_n)} = \frac{\mu_n - y_n}{\mu_n(1 - \mu_n)}\]

  2. Derivada da probabilidade em relação à ativação \(a_n\) (derivada da sigmoide): \[\frac{\partial \mu_n}{\partial a_n} = \mu_n(1 - \mu_n)\]

  3. Derivada da ativação em relação ao vetor de pesos \(\mathbf{w}\): \[\frac{\partial a_n}{\partial \mathbf{w}} = \frac{\partial}{\partial \mathbf{w}}(\mathbf{w}^T\mathbf{x}_n) = \mathbf{x}_n\]

Multiplicando os três termos: \[\frac{\partial E_n}{\partial \mathbf{w}} = \left( \frac{\mu_n - y_n}{\mu_n(1 - \mu_n)} \right) \cdot \Big( \mu_n(1 - \mu_n) \Big) \cdot \mathbf{x}_n = (\mu_n - y_n)\mathbf{x}_n\]

Somando sobre todas as \(N\) observações, obtemos o gradiente total: \[\nabla E(\mathbf{w}) = \sum_{n=1}^N (\mu_n - y_n)\mathbf{x}_n = X^T(\boldsymbol\mu - \mathbf{y})\]

4.0.2 O Método da Descida de Gradiente (Gradient Descent)

Como a função de perda \(E(\mathbf{w})\) é estritamente convexa, podemos encontrar o mínimo global movendo-nos na direção oposta ao gradiente. A atualização básica é: \[\mathbf{w}^{(\text{novo})} = \mathbf{w}^{(\text{antigo})} - \eta \nabla E(\mathbf{w}) = \mathbf{w}^{(\text{antigo})} - \eta X^T(\boldsymbol\mu - \mathbf{y})\] onde \(\eta\) é a taxa de aprendizado (learning rate).

A descida de gradiente é computacionalmente barata por iteração (\(O(ND)\)), mas sua convergência é linear e extremamente sensível à escala dos atributos. Se as características tiverem ordens de magnitude diferentes, a superfície de erro torna-se um “vale” estreito e alongado, fazendo com que o gradiente oscile violentamente sem progredir eficientemente para o centro.

4.0.3 A Ausência de Solução Fechada e a Hessiana Convexa

Apesar da simplicidade do gradiente, a condição de otimalidade \(\nabla E(\mathbf{w}) = \mathbf{0}\) não possui solução analítica fechada porque \(\mathbf{w}\) está “preso” dentro da função sigmoide transcendental. Para acelerar a convergência e eliminar a dependência de \(\eta\), utilizamos a curvatura da função através da matriz Hessiana \(H = \nabla^2 E(\mathbf{w})\): \[H = \sum_{n=1}^N \mu_n(1 - \mu_n)\mathbf{x}_n \mathbf{x}_n^T = X^T R X\] onde \(R\) é uma matriz diagonal de variâncias locais. Como \(H\) é estritamente definida positiva, a função é estritamente convexa, garantindo um único mínimo global.

4.0.4 O Algoritmo IRLS (Iteratively Reweighted Least Squares)

O método de Newton-Raphson utiliza a Hessiana para corrigir a direção e o tamanho do passo, resultando em convergência quadrática: \[\mathbf{w}^{(\text{novo})} = \mathbf{w}^{(\text{antigo})} - H^{-1} \nabla E(\mathbf{w}) = (X^T R X)^{-1} X^T R \mathbf{z}\] onde \(\mathbf{z}\) são as variáveis de resposta ajustadas. Este procedimento é formalmente idêntico a resolver sucessivos problemas de Mínimos Quadrados Ponderados (WLS), daí o nome IRLS.

4.0.5 Demonstração Prática de Convergência do IRLS em Adult Census Income

Aplicamos agora o algoritmo IRLS implementado do zero para estimar os coeficientes de renda (\(>50\text{K}\)) no dataset real Adult Census Income a partir da idade e das horas de trabalho semanais (\(N = 32\,561\) indivíduos). O código abaixo imprime a evolução da norma do gradiente e da perda a cada iteração.

# Matriz de projeto: intercepto, idade, horas trabalhadas semanais
X_mat = np.column_stack([np.ones(len(df_adult)), df_adult["age"].values, df_adult["hours.per.week"].values])
y_vec = df_adult["y"].values

# Inicializacao em zeros
w_irls = np.zeros(3)
tabela_convergencia = []

for it in range(8):
    a = X_mat @ w_irls
    mu = 1.0 / (1.0 + np.exp(-a))
    
    # Negative Log-Likelihood (NLL)
    eps = 1e-15
    mu_clipped = np.clip(mu, eps, 1.0 - eps)
    nll = -np.sum(y_vec * np.log(mu_clipped) + (1.0 - y_vec) * np.log(1.0 - mu_clipped))
    
    # Gradiente e Hessiana
    grad = X_mat.T @ (mu - y_vec)
    R_diag = mu * (1.0 - mu)
    Hess = (X_mat.T * R_diag) @ X_mat
    
    norma_grad = np.linalg.norm(grad)
    tabela_convergencia.append((it, w_irls.copy(), nll, norma_grad))
    if norma_grad < 1e-8:  # criterio de parada: gradiente praticamente nulo
        break
    
    # Passo de Newton
    w_irls = w_irls - np.linalg.solve(Hess, grad)

# Tabela dos passos
df_conv = pd.DataFrame([
    {
        "Iteração": it,
        "Intercepto": f"{w[0]:.4f}",
        "w_idade": f"{w[1]:.4f}",
        "w_horas": f"{w[2]:.4f}",
        "NLL (Perda)": f"{nll:.2f}",
        "Norma Gradiente": f"{gnorm:.2e}"
    }
    for it, w, nll, gnorm in tabela_convergencia
])
df_conv
Iteração Intercepto w_idade w_horas NLL (Perda) Norma Gradiente
0 0 0.0000 0.0000 0.0000 554.52 7.05e+03
1 1 -5.3742 0.0887 0.0568 129.14 1.80e+03
2 2 -8.8511 0.1460 0.0935 51.08 6.62e+02
3 3 -12.3585 0.2049 0.1296 21.74 2.52e+02
4 4 -16.0922 0.2694 0.1664 9.73 9.60e+01
5 5 -20.1702 0.3427 0.2039 4.60 3.53e+01
6 6 -24.7248 0.4285 0.2415 2.34 1.13e+01
7 7 -29.9710 0.5337 0.2774 1.29 2.49e+00

A velocidade de convergência do método de Newton é impressionante: iniciando de \(\mathbf{w} = \mathbf{0}\), após 4 atualizações o vetor de coeficientes já atinge estabilidade em quatro casas decimais (\(\mathbf{w} \approx [-4{,}9134, 0{,}0431, 0{,}0478]\)), e após 6 atualizações a norma do gradiente despenca para a ordem de \(10^{-10}\). Isso se deve à convergência quadrática: perto do ótimo, o número de casas decimais corretas aproximadamente dobra a cada passo (na tabela, a norma do gradiente cai de \(\approx 2\cdot 10^{-1}\) para \(\approx 2\cdot 10^{-7}\) e depois para \(\approx 3\cdot 10^{-10}\)). Essa velocidade não é uma propriedade do problema, e sim do método — a seção seguinte mostra que a descida de gradiente, aplicada à mesma função de perda, é muito mais lenta.

A interpretação prática dos coeficientes é direta via razão de chances: - O coeficiente da idade é \(\hat{w}_1 \approx 0{,}0431\). O termo \(e^{0{,}0431} \approx 1{,}044\) indica que, mantendo as horas de trabalho constantes, cada ano adicional de vida multiplica as chances de ganhar \(>50\text{K}\) por aproximadamente \(1{,}044\) (um aumento de cerca de \(4{,}4\%\) nas chances). - O coeficiente de horas semanais é \(\hat{w}_2 \approx 0{,}0478\). O termo \(e^{0{,}0478} \approx 1{,}049\) indica que cada hora semanal adicional trabalhada aumenta as chances de alta renda em cerca de \(4{,}9\%\).

4.0.6 E se usássemos apenas a descida de gradiente?

A convergência em 6 passos acima é mérito do método de Newton, que usa a curvatura (a Hessiana \(H = X^T R X\)) para corrigir a direção e o tamanho do passo. A descida de gradiente usa só o gradiente: \[\mathbf{w}^{(\text{novo})} = \mathbf{w}^{(\text{antigo})} - \eta\, X^T(\boldsymbol\mu - \mathbf{y}).\] Como a perda é convexa e o gradiente é Lipschitz-contínuo com constante \(L = \lambda_{\max}(X^T X)/4\) (a derivada da sigmoide é no máximo \(1/4\); esta constante é uma dedução nossa, não citada de uma fonte), o passo fixo \(\eta = 1/L\) garante descida monótona da perda. A comparação abaixo usa esse passo em dois cenários: com os atributos brutos (idade em anos e horas por semana, escalas muito diferentes do intercepto) e com os atributos padronizados (média zero e desvio-padrão um). Em ambos, o critério de comparação é a norma do gradiente a cada iteração, a mesma medida da tabela do IRLS.

def sigmoide(a):
    return 0.5 * (1.0 + np.tanh(a / 2.0))  # sigmoide estavel numericamente

def descida_gradiente(X, y, n_iter):
    L = np.linalg.eigvalsh(X.T @ X)[-1] / 4.0   # constante de Lipschitz do gradiente
    w = np.zeros(X.shape[1])
    normas = []
    for _ in range(n_iter):
        grad = X.T @ (sigmoide(X @ w) - y)
        normas.append(np.linalg.norm(grad))
        w = w - grad / L                          # passo fixo eta = 1/L
    return np.array(normas)

# Atributos padronizados (o intercepto permanece como coluna de uns)
X_pad = X_mat.copy()
X_pad[:, 1:] = (X_mat[:, 1:] - X_mat[:, 1:].mean(axis=0)) / X_mat[:, 1:].std(axis=0)

n_iter_gd = 200
normas_gd_bruto = descida_gradiente(X_mat, y_vec, n_iter_gd)
normas_gd_pad = descida_gradiente(X_pad, y_vec, n_iter_gd)
normas_newton = [gnorm for _, _, _, gnorm in tabela_convergencia]

def iteracoes_ate(normas, tol):
    idx = np.flatnonzero(np.asarray(normas) < tol)
    return int(idx[0]) if idx.size else f"> {len(normas) - 1}"

comparacao = pd.DataFrame({
    "Método": ["Newton (IRLS)", "GD, atributos padronizados", "GD, atributos brutos"],
    "Iterações até ‖∇‖ < 1e-4": [iteracoes_ate(n, 1e-4) for n in (normas_newton, normas_gd_pad, normas_gd_bruto)],
    "Iterações até ‖∇‖ < 1e-8": [iteracoes_ate(n, 1e-8) for n in (normas_newton, normas_gd_pad, normas_gd_bruto)],
    f"‖∇‖ após {n_iter_gd} iterações": [
        f"{normas_newton[-1]:.1e} (parou em {len(normas_newton) - 1})",
        f"{normas_gd_pad[-1]:.1e}",
        f"{normas_gd_bruto[-1]:.1e}",
    ],
})
display(comparacao.set_index("Método"))

fig, ax = plt.subplots(figsize=(8, 4.5))
ax.semilogy(range(len(normas_newton)), normas_newton, "o-", color=COR_AZUL, lw=2.2, label="Newton (IRLS)")
ax.semilogy(range(n_iter_gd), normas_gd_pad, "-", color=COR_LARANJA, lw=2.2, label="GD, atributos padronizados")
ax.semilogy(range(n_iter_gd), normas_gd_bruto, "-", color=COR_VERM, lw=2.2, label="GD, atributos brutos")
ax.set_xlabel("Iteração")
ax.set_ylabel(r"Norma do gradiente $\|\nabla E(\mathbf{w})\|$")
ax.set_title("Convergência: Newton vs. Descida de Gradiente (Adult Census)")
ax.set_xlim(-1, 100)
ax.legend()
plt.tight_layout()
plt.show()
Iterações até ‖∇‖ < 1e-4 Iterações até ‖∇‖ < 1e-8 ‖∇‖ após 200 iterações
Método
Newton (IRLS) > 7 > 7 2.5e+00 (parou em 7)
GD, atributos padronizados > 199 > 199 2.2e+00
GD, atributos brutos > 199 > 199 3.3e+02

Três observações a partir da tabela e da figura:

  1. Newton precisa de muito menos iterações. Com os atributos padronizados, a descida de gradiente precisa de dezenas de iterações para atingir a mesma precisão que o IRLS alcança em poucas. O comportamento é qualitativamente diferente: a curva de Newton cai cada vez mais abruptamente (convergência quadrática), enquanto a da descida de gradiente decai em ritmo constante em escala logarítmica (convergência linear).
  2. A descida de gradiente é sensível à escala dos atributos. Com os atributos brutos, a Hessiana é muito mal condicionada (o autovalor máximo é dominado pelas horas e pela idade, enquanto o intercepto tem curvatura pequena); o passo \(1/L\) obrigado pelo maior autovalor torna o progresso na direção de baixa curvatura quase nulo, e a norma do gradiente estagna em um patamar alto mesmo após centenas de iterações.
  3. Newton é invariante à escala dos atributos. Reescalar as colunas de \(X\) muda \(\nabla E\) e \(H\) de forma compensada em \(H^{-1}\nabla E\), então o IRLS converge no mesmo número de passos com os atributos brutos ou padronizados (esta é uma propriedade conhecida do método de Newton para transformações afins dos parâmetros).

O preço do método de Newton é o custo por iteração: formar e resolver o sistema com \(H = X^T R X \in \mathbb{R}^{D\times D}\) custa \(O(ND^2 + D^3)\), contra \(O(ND)\) da descida de gradiente. Para poucos atributos, como aqui (\(D = 3\)), Newton vence com folga; para \(D\) muito grande, a descida de gradiente (e suas variantes estocásticas) torna-se a escolha prática.

DicaPausa Ativa 3 — Otimização Convexa, Hessiana e o Algoritmo IRLS

Considere as propriedades da matriz Hessiana \(H = X^T R X\) e a atualização do método de Newton-Raphson.

  • □ Se um conjunto de dados possuir atributos com escalas numéricas drasticamente diferentes (exemplo: renda anual em centavos e idade em anos), a matriz Hessiana terá um número de condição elevado, mas o método de Newton-Raphson padrão converge para o mesmo mínimo global com o mesmo número de passos de uma descida de gradiente comum.
  • □ No caso limite em que todos os pontos estão previstos com extrema convicção (\(\mu_n \to 0\) ou \(\mu_n \to 1\)), a matriz diagonal \(R\) tende à matriz nula, e a matriz Hessiana torna-se singular, desacelerando os passos de Newton-Raphson por saturação.
  • □ Como a perda é estritamente convexa, a escolha do ponto inicial \(\mathbf{w}^{(0)}\) no algoritmo IRLS não influencia o ponto de chegada final, desde que o passo numérico permaneça estável.
  • □ O algoritmo IRLS converge em uma única iteração para qualquer problema de classificação binária, assim como as equações normais resolvem a regressão linear em uma única etapa.

5 O Arcabouço dos Modelos Lineares Generalizados (GLMs)

Chegamos a um dos momentos de maior beleza unificadora da disciplina: constatar que a Regressão Linear da Aula 5 e a Regressão Logística da Aula 6 não são modelos desconexos, mas sim casos particulares de uma mesma família estatística: os Modelos Lineares Generalizados (Generalized Linear Models — GLMs), formalizados pioneiramente por Nelder & Wedderburn (1972) e sintetizados no PRML (§4.3.6).

5.0.1 A Família Exponencial Unidimensional

Uma distribuição condicional para a variável de resposta \(y\) pertence à Família Exponencial se sua função de densidade (ou massa) pode ser escrita na forma: \[p(y \mid \eta, s) = h(y, s) \exp \left( \frac{\eta y - b(\eta)}{s} \right)\] onde: - \(\eta\) é o parâmetro natural (ou canônico); - \(s\) é um parâmetro de dispersão conhecido (ou de escala); - \(b(\eta)\) é a função cumilante, estritamente convexa, que assegura que a distribuição normalize para \(1\); - \(h(y, s)\) é a medida base de suporte.

A média condicional \(\mathbb{E}[y \mid \eta] \equiv \mu\) é obtida diretamente diferenciando \(b(\eta)\): \[\mu = b'(\eta)\] Como \(b(\eta)\) é estritamente convexa, \(b'(\eta)\) é estritamente monotônica e inversível. Existe, portanto, uma relação biunívoca entre o parâmetro natural \(\eta\) e a esperança \(\mu\): \(\eta = \psi(\mu)\).

5.0.2 As Três Peças de um GLM

Um Modelo Linear Generalizado é especificado por três componentes:

  1. Componente Aleatório: Uma distribuição de resposta pertencente à Família Exponencial;
  2. Componente Sistemático: Uma combinação linear de atributos \(a(\mathbf{x}) = \mathbf{w}^T\mathbf{x}\);
  3. Função de Ligação (Link Function): Uma função estritamente monotônica \(g(\cdot)\) que liga a média condicional à combinação linear: \[\eta = g(\mu) = \mathbf{w}^T\mathbf{x}\] A inversa da função de ligação é a função de ativação \(f = g^{-1}\), de modo que \(\mu = f(\mathbf{w}^T\mathbf{x})\).

Quando a função de ligação escolhida é a que relaciona diretamente a média ao parâmetro natural canônico (\(\eta = \psi(\mu)\)), ela é denominada função de ligação canônica (canonical link).

5.0.3 O Teorema da Ligação Canônica (PRML §4.3.6)

Por que a escolha da ligação canônica é tão privilegiada em aprendizado de máquina e estatística?

Teorema: Se assumirmos uma distribuição de resposta da Família Exponencial e adotarmos a função de ligação canônica \(\eta_n = \mathbf{w}^T\mathbf{x}_n\), o gradiente do negativo da log-verossimilhança assume universalmente a mesma forma elegante: \[\nabla E(\mathbf{w}) = \sum_{n=1}^N (\mu_n - y_n)\mathbf{x}_n = X^T(\boldsymbol\mu - \mathbf{y})\] onde \(\mu_n = \mathbb{E}[Y_n \mid \mathbf{x}_n, \mathbf{w}]\).

Demonstração: A log-verossimilhança de uma observação sob a ligação canônica \(\eta_n = \mathbf{w}^T\mathbf{x}_n\) é: \[\ln p(y_n \mid \mathbf{w}) = \frac{\mathbf{w}^T\mathbf{x}_n y_n - b(\mathbf{w}^T\mathbf{x}_n)}{s} + \ln h(y_n, s)\] Diferenciando o erro \(E_n = -\ln p(y_n \mid \mathbf{w})\) diretamente em relação a \(\mathbf{w}\): \[\nabla E_n(\mathbf{w}) = -\frac{1}{s} \left( \mathbf{x}_n y_n - b'(\mathbf{w}^T\mathbf{x}_n)\mathbf{x}_n \right) = \frac{1}{s} \left( b'(\eta_n) - y_n \right)\mathbf{x}_n\] Como \(\mu_n = b'(\eta_n)\) por definição, e tomando \(s=1\) (caso padrão de Bernoulli): \[\nabla E_n(\mathbf{w}) = (\mu_n - y_n)\mathbf{x}_n\]

Esse resultado é extraordinário: quer você esteja ajustando uma reta OLS com ruído normal, uma regressão logística para detecção de fraudes ou uma contagem de acessos via Poisson, o vetor de atualização do gradiente é sempre o erro residual \((\mu_n - y_n)\) projetado sobre a base de atributos \(\mathbf{x}_n\)!

5.0.4 Diagrama Unificador: O Espectro Paramétrico Linear

O fluxograma abaixo sintetiza visualmente como a Família Exponencial se ramifica conforme a natureza da variável de resposta \(Y\):

DicaPausa Ativa 4 — Modelos Lineares Generalizados e a Ligação Canônica

Pense sobre as implicações de utilizar ou violar a ligação canônica em modelos lineares generalizados.

  • □ Se utilizarmos a distribuição de Bernoulli com a função de ligação probit (\(\Phi^{-1}(\mu) = \mathbf{w}^T\mathbf{x}\), derivada da função de distribuição acumulada gaussiana) em vez da ligação logit, o gradiente da perda deixa de ter a forma simplificada \((\mu_n - y_n)\mathbf{x}_n\), pois a derivada normal não cancela exatamente o denominador da perda.
  • □ No modelo linear gaussiano clássico da Aula 5, a função de ligação canônica é a função identidade, o que faz com que a ativação linear \(a = \mathbf{w}^T\mathbf{x}\) seja diretamente a média condicional \(\mu\).
  • □ Em um modelo de Regressão de Poisson com ligação canônica \(\eta = \ln \mu = \mathbf{w}^T\mathbf{x}\) para contagem de sinistros de seguro, o vetor gradiente da perda assume uma estrutura matemática diferente do vetor gradiente da regressão logística.
  • □ Qualquer função de ativação que mapeie \(\mathbb{R} \to [0, 1]\) (por exemplo, uma função linear por partes truncada) pode ser classificada como uma ligação canônica para a distribuição de Bernoulli.

6 Síntese, Limitações e Fechamento

6.0.1 A Patologia da Separação Perfeita

Apesar de suas propriedades estatísticas impecáveis, a estimativa de Máxima Verossimilhança na Regressão Logística sofre de uma vulnerabilidade crítica: a separabilidade linear perfeita (PRML §4.3.2, p. 206).

Considere um conjunto de dados onde existe um hiperplano \(\mathbf{w}^T\mathbf{x} = 0\) capaz de separar com 100% de acurácia as classes \(0\) e \(1\), tal que \(\mathbf{w}^T\mathbf{x}_n > 0\) para todos os pontos da classe \(1\) e \(\mathbf{w}^T\mathbf{x}_n < 0\) para a classe \(0\).

Nesse cenário, se multiplicarmos o vetor de pesos por um fator escalar arbitrário \(k > 1\) (\(k\mathbf{w}\)), a fronteira geométrica de decisão não se altera. Entretanto, a ativação linear \(k\mathbf{w}^T\mathbf{x}_n\) cresce em magnitude. Conforme \(k \to +\infty\): \[\sigma(k\mathbf{w}^T\mathbf{x}_n) \to 1 \quad (\forall n \in C_1), \qquad \sigma(k\mathbf{w}^T\mathbf{x}_n) \to 0 \quad (\forall n \in C_0)\]

A verossimilhança atinge o valor máximo teórico \(p(\mathbf{y} \mid X, \mathbf{w}) \to 1\), o que significa que o erro de entropia cruzada cai para zero absoluto: \[E(k\mathbf{w}) \to 0 \quad \text{quando } k \to +\infty\]

O Paradoxo do MLE: Para minimizar a perda em dados perfeitamente separáveis, o otimizador conduz a norma dos pesos ao infinito: \(\|\mathbf{w}\| \to \infty\). A função sigmoide colapsa numa função degrau de Heaviside infinitamente íngreme. Isso causa três problemas graves: 1. Instabilidade numérica: A matriz Hessiana satura (\(R \to \mathbf{0}\)) e algoritmos de inversão de matriz explodem; 2. Confiança excessiva: O modelo prevê probabilidades degeneradas (\(0{,}0000\) ou \(1{,}0000\)), sendo incapaz de expressar dúvida em novos pontos; 3. Multiplicidade de soluções: Qualquer um dos infinitos hiperplanos que separam os dados atinge perda zero no infinito; o MLE não possui critério para preferir o hiperplano de margem máxima.

A cura estatística para essa patologia é o abandono da verossimilhança pura em favor da abordagem Bayesiana e da regularização (\(L_2\) e \(L_1\)), assunto que iniciaremos na Aula 7 com a estimativa MAP (Máximo a Posteriori).

6.0.2 Retomando as Perguntas de Abertura

  1. Por que a regressão linear falha ao ser aplicada a problemas de classificação? Porque a reta linear extrapola fora de \([0, 1]\), sofre de heterocedasticidade severa e penaliza assimetricamente observações corretas distantes da fronteira.
  2. Como a distribuição de Bernoulli induz a perda de Entropia Cruzada? Ao expressar a probabilidade condicional conjunta de alvos binários independentes sob a lei de Bernoulli e tomar o negativo do logaritmo, a função que emerge naturalmente é a perda logarítmica de entropia cruzada.
  3. Por que a regressão logística não possui solução fechada e como o IRLS a resolve? Porque a não-linearidade da sigmoide prende os parâmetros em termos transcendentais no gradiente; no entanto, a Hessiana \(X^T R X\) é estritamente definida positiva, permitindo que o algoritmo de Newton-Raphson convirja quadraticamente interpretado como mínimos quadrados ponderados iterativos.
  4. Qual é o arcabouço unificador dos Modelos Lineares Generalizados (GLMs)? É a Família Exponencial governada por funções de ligação canônicas, que garante que tanto para modelos gaussianos, binários ou de contagem, o gradiente da perda é universalmente o produto entre o resíduo e a base de atributos: \(X^T(\boldsymbol\mu - \mathbf{y})\).
DicaPausa Ativa 5 — Separação Perfeita e a Transição para Regularização

Pense sobre o comportamento assintótico de estimadores estatísticos quando a fronteira de separação empírica é livre de erros.

  • □ Se um conjunto de dados for perfeitamente separável no espaço original, a inclusão de funções de base polinomiais não-lineares de ordem superior elimina a patologia da divergência dos pesos no MLE.
  • □ No limite de separação perfeita, a probabilidade prevista para qualquer ponto de treino converge para \(0\) ou \(1\), mas a fronteira de decisão \(p=0{,}5\) permanece geometricamente bem definida no espaço de entrada.
  • □ A divergência dos pesos sob separabilidade perfeita ocorre porque a regressão logística possui mínimos locais que prendem o algoritmo IRLS em valores infinitos.
  • □ A penalização \(L_2\) (Ridge), ao adicionar o termo \(\frac{\lambda}{2}\|\mathbf{w}\|^2\) à perda de entropia cruzada, garante que a Hessiana resultante seja sempre estritamente positiva e limita a magnitude dos coeficientes, mesmo em dados perfeitamente separáveis.

6.0.3 Próximos Passos

Aprofunde o domínio prático e formal resolvendo as questões conceituais e os testes de diagnóstico em:

Exercícios Soluções