Em estatística, um modelo aditivo generalizado (MAG; em inglês: generalized additive model, GAM) é um modelo linear generalizado no qual a variável resposta linear depende linearmente de funções suaves desconhecidas de algumas variáveis preditoras, e cujo interesse se concentra na inferência sobre essas funções suaves.
Os GAMs foram desenvolvidos originalmente por Trevor Hastie e Robert Tibshirani para combinar propriedades dos modelos lineares generalizados com as dos modelos aditivos. Eles podem ser interpretados como a generalização discriminativa do modelo gerativo naive Bayes.
O modelo relaciona uma variável resposta univariada, Y, a algumas variáveis preditoras, xi. Especifica-se para Y uma distribuição da família exponencial (por exemplo, as distribuições normal, binomial ou de Poisson), juntamente com uma função de ligação g (por exemplo, a função identidade ou a função logarítmica) que relaciona o valor esperado de Y às variáveis preditoras por meio de uma estrutura como
As funções fi podem ter uma forma paramétrica especificada (por exemplo, um polinômio ou um spline de regressão não penalizado de uma variável) ou podem ser especificadas de forma não paramétrica ou semiparamétrica, simplesmente como "funções suaves", a serem estimadas por métodos não paramétricos. Assim, um GAM típico pode usar uma função de suavização de diagrama de dispersão, como uma média ponderada local, para f1(x1), e um modelo de fatores para f2(x2). Essa flexibilidade, que permite ajustes não paramétricos com pressupostos menos rígidos sobre a relação real entre resposta e preditor, oferece a possibilidade de ajustes melhores aos dados do que os modelos puramente paramétricos, mas, possivelmente, com alguma perda de interpretabilidade.
Sabe-se desde a década de 1950 (por meio do teorema da representação de Kolmogorov–Arnold) que qualquer função contínua multivariada pode ser representada como somas e composições de funções univariadas,
Embora o teorema da representação de Kolmogorov–Arnold afirme a existência de uma função dessa forma, ele não fornece nenhum mecanismo para construí-la. Existem algumas demonstrações construtivas, mas elas tendem a exigir funções altamente complicadas (isto é, fractais) e, portanto, não são adequadas para abordagens de modelagem. Por isso, o modelo aditivo generalizado descarta a soma externa e exige, em vez disso, que a função pertença a uma classe mais simples,
em que é uma função monótona suave. Denotando por a inversa de , isso é tradicionalmente escrito como
Quando essa função aproxima a esperança de alguma quantidade observada, ela pode ser escrita como
Essa é a formulação padrão de um modelo aditivo generalizado. Demonstrou-se então[necessário esclarecer] que o algoritmo de backfitting sempre converge para essas funções.
A classe de modelos GAM é bastante ampla, dado que função suave é uma categoria bastante abrangente. Por exemplo, uma covariável pode ser multivariada e a correspondente pode ser uma função suave de várias variáveis, ou pode ser a função que associa o nível de um fator ao valor de um efeito aleatório. Outro exemplo é um termo de coeficiente variável (regressão geográfica), como , em que e são ambas covariáveis. Ou, se for ela própria a observação de uma função, pode-se incluir um termo como (por vezes chamado de termo de regressão de sinal). também pode ser uma função paramétrica simples, como as usadas em qualquer modelo linear generalizado. A classe de modelos foi generalizada em várias direções, notadamente para além das distribuições de resposta da família exponencial, para além da modelagem apenas da média e para além de dados univariados.
O método original de ajuste de GAMs estimava os componentes suaves do modelo usando suavizadores não paramétricos (por exemplo, splines de suavização ou suavizadores de regressão linear local) por meio do algoritmo de backfitting. O backfitting funciona por suavização iterativa de resíduos parciais e fornece um método de estimação modular muito geral, capaz de usar uma ampla variedade de métodos de suavização para estimar os termos . Uma desvantagem do backfitting é a dificuldade de integrá-lo à estimação do grau de suavidade dos termos do modelo, de modo que, na prática, o usuário precisa defini-los ou escolher entre um conjunto modesto de níveis de suavização predefinidos.
Se as forem representadas por splines de suavização, o grau de suavidade pode ser estimado como parte do ajuste do modelo usando validação cruzada generalizada, ou por máxima verossimilhança restrita (REML, por vezes chamada de "GML"), que explora a dualidade entre suavizadores spline e efeitos aleatórios gaussianos. Essa abordagem de spline completo tem custo computacional , em que é o número de observações da variável resposta, o que a torna pouco prática para conjuntos de dados moderadamente grandes. Métodos mais recentes trataram esse custo computacional seja pela redução prévia do tamanho da base usada na suavização (redução de posto), seja pela obtenção de representações esparsas dos suavizadores por meio de campos aleatórios de Markov, que se prestam ao uso de métodos de matrizes esparsas na computação. Esses métodos computacionalmente mais eficientes usam GCV (ou AIC ou critério semelhante) ou REML, ou adotam uma abordagem totalmente bayesiana para a inferência sobre o grau de suavidade dos componentes do modelo. A estimação do grau de suavidade via REML pode ser vista como um método de Bayes empírico.
Uma abordagem alternativa, com vantagens particulares em contextos de alta dimensionalidade, é o uso de boosting, embora isso normalmente exija bootstrapping para a quantificação da incerteza. Verificou-se que GAMs ajustados por bagging e boosting geralmente superam GAMs ajustados por métodos de spline.
Muitas implementações modernas de GAMs e de suas extensões são construídas em torno da abordagem de suavização de posto reduzido, pois ela permite uma estimação bem fundamentada da suavidade dos componentes suaves a um custo computacional comparativamente modesto, além de facilitar a implementação de diversas extensões do modelo de uma forma que é mais difícil com outros métodos. Em sua forma mais simples, a ideia consiste em substituir as funções suaves desconhecidas do modelo por expansões em bases
em que as são funções de base conhecidas, geralmente escolhidas por suas boas propriedades do ponto de vista da teoria da aproximação (por exemplo, B-splines ou splines de placa fina de posto reduzido), e os são coeficientes a serem estimados como parte do ajuste do modelo. A dimensão da base é escolhida grande o suficiente para que se espere um sobreajuste aos dados disponíveis (evitando, assim, o viés decorrente de uma simplificação excessiva do modelo), mas pequena o suficiente para preservar a eficiência computacional. Se , então o custo computacional da estimação do modelo dessa maneira será .
As só são identificáveis a menos de um termo de intercepto (poder-se-ia somar qualquer constante a e subtraí-la de sem alterar em nada as previsões do modelo); por isso, é necessário impor restrições de identificabilidade aos termos suaves para eliminar essa ambiguidade. A inferência mais precisa sobre as é geralmente obtida usando as restrições de soma zero
isto é, exigindo que a soma de cada avaliada nos valores observados de sua covariável seja igual a zero. Essas restrições lineares podem ser impostas mais facilmente por reparametrização na etapa de construção da base, de modo que, a seguir, se supõe que isso já tenha sido feito.
Substituídas todas as do modelo por essas expansões em bases, o GAM transforma-se em um modelo linear generalizado (MLG), com uma matriz do modelo que simplesmente contém as funções de base avaliadas nos valores observados de . No entanto, como as dimensões das bases, , foram escolhidas um pouco maiores do que se acredita ser necessário para os dados, o modelo está superparametrizado e sofrerá sobreajuste se for estimado como um MLG comum. A solução para esse problema é penalizar o afastamento da suavidade no processo de ajuste do modelo, controlando o peso atribuído às penalidades de suavização por meio de parâmetros de suavização. Por exemplo, considere-se a situação em que todos os suavizadores são funções univariadas. Escrevendo todos os parâmetros em um único vetor, , suponha-se que seja a desviância (o dobro da diferença entre a log-verossimilhança saturada e a log-verossimilhança do modelo). Minimizar a desviância pelo método usual de mínimos quadrados iterativamente reponderados resultaria em sobreajuste; por isso, busca-se o que minimiza
em que as penalidades da integral do quadrado da segunda derivada servem para penalizar a ondulação (falta de suavidade) das durante o ajuste, e os parâmetros de suavização controlam o equilíbrio entre a qualidade do ajuste e a suavidade do modelo. No exemplo, garantiria que a estimativa de fosse uma reta em .
Dada a expansão em base de cada , as penalidades de ondulação podem ser expressas como formas quadráticas nos coeficientes do modelo. Ou seja, pode-se escrever
em que é uma matriz de coeficientes conhecidos, calculável a partir da penalidade e da base, é o vetor de coeficientes de , e é simplesmente completada com zeros, de modo que a segunda igualdade se verifique e a penalidade possa ser escrita em termos do vetor completo de coeficientes . Muitas outras penalidades de suavização podem ser escritas da mesma forma e, dados os parâmetros de suavização, o problema de ajuste do modelo passa a ser
cuja solução pode ser obtida por uma versão penalizada do algoritmo usual de mínimos quadrados iterativamente reponderados (IRLS) para MLGs: o algoritmo permanece inalterado, exceto pelo fato de que a soma das penalidades quadráticas é adicionada ao objetivo de mínimos quadrados de trabalho em cada iteração.
A penalização tem vários efeitos sobre a inferência em relação a um MLG comum. Em primeiro lugar, as estimativas estão sujeitas a certo viés de suavização, que é o preço a pagar por limitar a variância do estimador por meio da penalização. No entanto, se os parâmetros de suavização forem escolhidos adequadamente, o viés de suavização (ao quadrado) introduzido pela penalização deve ser menor do que a redução de variância que ela produz, de modo que o efeito líquido seja uma redução do erro quadrático médio de estimação em relação à ausência de penalização. Um efeito relacionado da penalização é que a noção de graus de liberdade de um modelo precisa ser modificada para levar em conta a ação das penalidades na redução da liberdade de variação dos coeficientes. Por exemplo, se é a matriz diagonal dos pesos do IRLS na convergência e é a matriz do modelo GAM, então os graus de liberdade efetivos do modelo são dados por , em que
é a matriz de graus de liberdade efetivos. De fato, somando apenas os elementos diagonais de correspondentes aos coeficientes de , obtêm-se os graus de liberdade efetivos da estimativa de .
O viés de suavização complica a estimação intervalar nesses modelos, e a abordagem mais simples acaba por envolver uma perspectiva bayesiana. Compreender essa visão bayesiana da suavização também ajuda a entender as abordagens REML e totalmente bayesiana para a estimação dos parâmetros de suavização. Em certo nível, as penalidades de suavização são impostas porque se acredita que funções suaves são mais prováveis do que funções onduladas e, se isso for verdade, pode-se igualmente formalizar essa noção atribuindo uma distribuição a priori à ondulação do modelo. Uma distribuição a priori muito simples poderia ser
(em que é o parâmetro de escala do MLG, introduzido apenas por conveniência posterior), mas é possível reconhecer imediatamente nela uma distribuição a priori normal multivariada com média e matriz de precisão . Como a penalidade deixa passar algumas funções sem penalização (retas, dadas as penalidades do exemplo), tem posto deficiente, e a distribuição a priori é, na verdade, imprópria, com matriz de covariância dada pela pseudoinversa de Moore–Penrose de (a impropriedade corresponde a atribuir variância infinita aos componentes não penalizados de um suavizador).
Se essa distribuição a priori for combinada com a verossimilhança do MLG, verifica-se que a moda a posteriori de é exatamente o obtido acima pelo IRLS penalizado. Além disso, tem-se o resultado para grandes amostras
que pode ser usado para produzir intervalos de confiança/credibilidade para os componentes suaves, . As distribuições a priori gaussianas de suavidade também são a base da inferência totalmente bayesiana com GAMs, bem como dos métodos que estimam GAMs como modelos mistos, que são essencialmente métodos de Bayes empírico.
Até aqui, tratou-se da estimação e da inferência dados os parâmetros de suavização, , mas estes também precisam ser estimados. Uma possibilidade é adotar uma abordagem totalmente bayesiana, definindo distribuições a priori para os parâmetros de suavização (em escala logarítmica) e usando simulação estocástica ou métodos de aproximação de ordem elevada para obter informações sobre a distribuição a posteriori dos coeficientes do modelo. Uma alternativa é escolher os parâmetros de suavização de modo a otimizar um critério de erro de previsão, como a validação cruzada generalizada (GCV) ou o critério de informação de Akaike (AIC). Por fim, pode-se optar por maximizar a verossimilhança marginal (REML), obtida integrando os coeficientes do modelo, , na densidade conjunta de :
Como é simplesmente a verossimilhança de , isso pode ser visto como a escolha de que maximiza a verossimilhança média de amostras aleatórias extraídas da distribuição a priori. A integral anterior geralmente não tem solução analítica, mas pode ser aproximada com precisão bastante elevada pelo método de Laplace.
A inferência sobre os parâmetros de suavização é a parte computacionalmente mais exigente da estimação/inferência do modelo. Por exemplo, otimizar o GCV ou a verossimilhança marginal normalmente requer otimização numérica por um método de Newton ou quase-Newton, sendo que cada valor experimental do vetor de parâmetros de suavização (em escala logarítmica) exige uma iteração de IRLS penalizado para avaliar o correspondente, juntamente com os demais ingredientes do escore GCV ou da verossimilhança marginal aproximada de Laplace (LAML). Além disso, a obtenção das derivadas do GCV ou da LAML, necessárias para a otimização, envolve diferenciação implícita para obter as derivadas de em relação aos logaritmos dos parâmetros de suavização, o que exige certo cuidado para preservar a eficiência e a estabilidade numérica.
GAMs ajustados por backfitting foram originalmente disponibilizados pela função gam da linguagem S, posteriormente portada para a linguagem R como o pacote gam. O procedimento GAM do SAS também oferece GAMs ajustados por backfitting. O pacote recomendado em R para GAMs é o mgcv, sigla de mixed GAM computational vehicle, baseado na abordagem de posto reduzido com seleção automática dos parâmetros de suavização. O procedimento GAMPL do SAS é uma implementação alternativa. Em Python, há o pacote PyGAM, com funcionalidades semelhantes às do mgcv do R. Há também o pacote InterpretML, que implementa uma abordagem de bagging e boosting. Existem muitos pacotes alternativos. Entre os exemplos estão os pacotes R mboost, que implementa uma abordagem de boosting; gss, que oferece os métodos de suavização por spline completo; VGAM, que oferece GAMs vetoriais; e gamlss, que oferece modelos aditivos generalizados para locação, escala e forma. O BayesX e sua interface em R oferecem GAMs e extensões por meio de MCMC e de métodos de verossimilhança penalizada. O software INLA implementa uma abordagem totalmente bayesiana baseada em representações por campos aleatórios de Markov que exploram métodos de matrizes esparsas.
Como exemplo de como os modelos podem ser estimados na prática com software, considere-se o pacote R mgcv. Suponha-se que o ambiente de trabalho do R contenha os vetores y, x e z e que se queira estimar o modelo
No R, poderiam ser executados os comandos
library(mgcv) # carrega o pacote
b = gam(y ~ s(x) + s(z))
Assim como a maioria das funções de modelagem do R, gam espera que seja fornecida uma fórmula do modelo, especificando a estrutura do modelo a ajustar. A variável resposta é indicada à esquerda do ~, enquanto a especificação do preditor linear é dada à direita. gam constrói as bases e as penalidades dos termos suaves, estima o modelo, incluindo seus parâmetros de suavização, e, no estilo padrão do R, retorna um objeto de modelo ajustado, que pode então ser examinado por meio de várias funções auxiliares, como summary, plot, predict e AIC.
Esse exemplo simples utilizou várias configurações padrão das quais é importante estar ciente. Por exemplo, foram assumidas uma distribuição gaussiana e a ligação identidade, e o critério de seleção dos parâmetros de suavização foi o GCV. Além disso, os termos suaves foram representados por "splines de regressão de placa fina penalizados", e a dimensão da base de cada um foi fixada em 10 (o que implica um máximo de 9 graus de liberdade após a imposição das restrições de identificabilidade). Um segundo exemplo ilustra como controlar esses aspectos. Suponha-se que se queira estimar o modelo
usando seleção dos parâmetros de suavização por REML, e que se espere que seja uma função relativamente complicada, que se deseja modelar com um spline cúbico de regressão penalizado. Para , também é preciso decidir se e estão naturalmente na mesma escala, de modo que um suavizador isotrópico, como o spline de placa fina, seja adequado (especificado por "s(v,w)"), ou se estão de fato em escalas diferentes, exigindo penalidades e parâmetros de suavização separados para e , como os fornecidos por um suavizador de produto tensorial. Supondo-se que se opte pela segunda alternativa, o seguinte código R estimaria o modelo:
b1 = gam(y ~ x + s(t,bs="cr",k=100) + te(v,w),family=poisson,method="REML")
que usa uma base de tamanho 100 para o suavizador de . A especificação da distribuição e da função de ligação usa os objetos "family", que são padrão no ajuste de MLGs em R ou S. Efeitos aleatórios gaussianos também podem ser adicionados ao preditor linear.
Esses exemplos pretendem apenas dar uma ideia muito básica de como o software para GAMs é usado; para mais detalhes, deve-se consultar a documentação dos diversos pacotes e as referências abaixo.
Como em qualquer modelo estatístico, é importante verificar os pressupostos de um GAM. Os gráficos de resíduos devem ser examinados da mesma forma que em qualquer MLG. Ou seja, os resíduos de desviância (ou outros resíduos padronizados) devem ser examinados em busca de padrões que possam sugerir uma violação substancial dos pressupostos de independência ou da relação média-variância do modelo. Isso geralmente envolve plotar os resíduos padronizados contra os valores ajustados e as covariáveis, a fim de identificar problemas na relação média-variância ou padrões não captados, e pode envolver também o exame de correlogramas (FACs) e/ou variogramas dos resíduos para verificar violações da independência. Se a relação média-variância do modelo estiver correta, os resíduos escalonados devem ter variância aproximadamente constante. Como MLGs e GAMs podem ser estimados por quase-verossimilhança, os detalhes da distribuição dos resíduos além da relação média-variância têm importância relativamente menor.
Um problema mais comum nos GAMs do que em outros MLGs é o risco de concluir erroneamente que os dados apresentam inflação de zeros. A dificuldade surge quando os dados contêm muitos zeros que podem ser modelados por uma distribuição de Poisson ou binomial com valor esperado muito baixo: a flexibilidade da estrutura do GAM frequentemente permite representar uma média muito baixa em alguma região do espaço das covariáveis, mas a distribuição dos resíduos padronizados não se parecerá em nada com a normalidade aproximada que os cursos introdutórios de MLG ensinam a esperar, mesmo que o modelo esteja perfeitamente correto.
A verificação adicional que os GAMs introduzem é a necessidade de confirmar se os graus de liberdade escolhidos são adequados. Isso é particularmente crítico quando se usam métodos que não estimam automaticamente a suavidade dos componentes do modelo. Ao usar métodos com seleção automática dos parâmetros de suavização, ainda é necessário verificar se a dimensão da base escolhida não foi excessivamente pequena, embora isso seja improvável se os graus de liberdade efetivos da estimativa de um termo estiverem confortavelmente abaixo da dimensão de sua base. Em todo caso, a verificação de baseia-se no exame de padrões nos resíduos em relação a . Isso pode ser feito com resíduos parciais sobrepostos ao gráfico de , ou por meio da permutação dos resíduos para construir testes de padrões residuais.
Quando os parâmetros de suavização são estimados como parte do ajuste do modelo, grande parte do que tradicionalmente seria considerado seleção de modelos já está incorporada ao processo de ajuste: a estimação dos parâmetros de suavização já selecionou entre uma rica família de modelos de diferentes complexidades funcionais. Contudo, a estimação dos parâmetros de suavização normalmente não remove por completo um termo suave do modelo, pois a maioria das penalidades deixa algumas funções sem penalização (por exemplo, as retas não são penalizadas pela penalidade da derivada do spline apresentada acima). Assim, permanece a questão de saber se um termo deve ou não constar do modelo. Uma abordagem simples para essa questão é adicionar uma penalidade extra a cada termo suave do GAM, que penaliza os componentes do suavizador que, de outro modo, ficariam sem penalização (e apenas esses). Cada penalidade extra tem seu próprio parâmetro de suavização, e a estimação prossegue como antes, mas agora com a possibilidade de que os termos sejam completamente penalizados até zero. Em contextos de alta dimensionalidade, pode fazer mais sentido realizar essa tarefa com o lasso ou com a regularização elastic net. O boosting também realiza automaticamente a seleção de termos como parte do ajuste.
Uma alternativa é usar métodos tradicionais de regressão stepwise para a seleção de modelos. Esse também é o método padrão quando os parâmetros de suavização não são estimados como parte do ajuste; nesse caso, geralmente se permite que cada termo suave assuma um dentre um pequeno conjunto de níveis de suavidade predefinidos no modelo, e a escolha entre eles é feita de forma stepwise. Os métodos stepwise operam comparando iterativamente modelos com ou sem determinados termos (ou, possivelmente, com diferentes níveis de complexidade dos termos) e requerem medidas de ajuste do modelo ou de significância dos termos para decidir qual modelo selecionar em cada etapa. Por exemplo, podem-se usar valores-p para testar se cada termo é igual a zero, a fim de decidir quais termos candidatos remover de um modelo, e comparar os valores do critério de informação de Akaike (AIC) de modelos alternativos.
O cálculo de valores-p para suavizadores não é simples, devido aos efeitos da penalização, mas existem aproximações disponíveis. O AIC pode ser calculado de duas maneiras para GAMs. O AIC marginal baseia-se na verossimilhança marginal (ver acima), com os coeficientes do modelo integrados. Nesse caso, a penalidade do AIC baseia-se no número de parâmetros de suavização (e de eventuais parâmetros de variância) do modelo. No entanto, devido ao fato bem conhecido de que a REML não é comparável entre modelos com diferentes estruturas de efeitos fixos, normalmente não é possível usar esse AIC para comparar modelos com termos suaves diferentes (já que seus componentes não penalizados atuam como efeitos fixos). É possível basear o AIC na verossimilhança marginal em que apenas os efeitos penalizados são integrados (o número de coeficientes não penalizados passa então a ser somado à contagem de parâmetros da penalidade do AIC), mas essa versão da verossimilhança marginal sofre da tendência à suavização excessiva que motivou originalmente o desenvolvimento da REML. Diante desses problemas, os GAMs são frequentemente comparados por meio do AIC condicional, no qual se usa a verossimilhança do modelo (e não a verossimilhança marginal) no AIC, e a contagem de parâmetros é tomada como os graus de liberdade efetivos do modelo.
Demonstrou-se que versões ingênuas do AIC condicional tendem excessivamente a selecionar modelos maiores em algumas circunstâncias, uma dificuldade atribuível à desconsideração da incerteza dos parâmetros de suavização no cálculo dos graus de liberdade efetivos; no entanto, a correção dos graus de liberdade efetivos para esse problema restabelece um desempenho razoável.
O sobreajuste pode ser um problema nos GAMs, especialmente se houver autocorrelação residual não modelada ou sobredispersão não modelada. A validação cruzada pode ser usada para detectar e/ou reduzir problemas de sobreajuste em GAMs (ou em outros métodos estatísticos), e os softwares costumam permitir aumentar o nível de penalização para forçar ajustes mais suaves. A estimação de um número muito grande de parâmetros de suavização também tende a ser estatisticamente desafiadora, e sabe-se que os critérios de erro de previsão (GCV, AIC etc.) ocasionalmente subsuavizam de forma substancial, em especial com tamanhos amostrais moderados, sendo a REML um pouco menos problemática nesse aspecto.
Quando apropriado, modelos mais simples, como os MLGs, podem ser preferíveis aos GAMs, a menos que estes melhorem substancialmente a capacidade preditiva (em conjuntos de validação) para a aplicação em questão.