Funções de perda e funções de escore

Como muitos devem saber, a regressão linear é um método estística para a estimação da média condicional de uma distribuição. Desta forma, para sua estimação, minimiza-se o MSE (mean square error) ou erro médio quadrático. Caso se deseje estimar a mediana, deve-se minimizar o MAE (mean absolute error) ou erro médio absoluto.

Genericamente falando, no entanto, pode-se optar por minimizar qualquer função de perda que se deseje. No entanto, a depender da distribuição real dos dados, uma função será mais apropriada do que as outras. Por exemplo, se os dados seguem a distribuição normal, a eficiência na estimação da média é muito superior à eficiência para a estimação da mediana. Por outro lado, se a distribuição dos dados seguir a distribuição de Laplace, a estimação da mediana é muito mais eficiente.

Há que se diferenciar, ainda, as funções de perda das funções de escore(https://en.wikipedia.org/wiki/Scoring_rule). Em outras palavras, há de se diferenciar estimação de predição. Na comunidade de Machine Learning geralmente parte-se direto dos dados para a previsão de valores, sem passar por qualquer estimação de parâmetros. Ou seja, escolhe-se uma função de escore apropriada e aplica-se um método (redes neurais, por exemplo) e busca-se a minimização do erro de predição.

Na aprendizagem de máquina, no entanto, perde-se a capacidade de explicar o fenômeno. Desta forma, em muitas ramos da ciência, como na Engenharia de Avaliações, opta-se primeiramente pelo ajuste de um modelo que explique o funcionamento do mercado para, subsequentemente, com este modelo, efetuar as previsões de valores necessárias.

Pois bem. Feita esta introdução, se os dados seguem a distribuição normal, como já dito acima, a estimação da média é mais eficiente. Um analista deveria preferir, portanto, um modelo de regressão linear, que minimiza o MSE, a um modelo de regressão quantílica, que minimiza o MAE, já que assim estará utilizando um método mais eficiente de estimação.

Se os dados são normais, no entanto, moda, média e mediana coincidem. Se a homocedasticidade for verificada, pode-se estimar qualquer que seja o percentil através da regressão linear (inclusive a mediana), bastando o ajuste apropriado do intercepto.

Os problemas começam a aparecer quando os dados não seguem a distribuição normal.

(Este artigo argumenta que boa parte dos dados que pensamos apresentar distribuição normal, na verdade, podem se ajustar melhor à distribuição de Laplace. Mas isto é assunto para um outro dia.)

Hoje pretendo discutir um caso que é muito comum na Engenharia de Avaliações, em que os dados apresentam uma distribuição lognormal.

Existem outras alternativas para a modelagem de dados lognormais, mais imaginemos que devamos escolher apenas entre regressão linear e regressão quantílica. Se os dados apresentam distribuição lognormal, então basta transformar os dados com a função logarítimca para a obtenção da normalidade. Uma vez normalizados os dados, pode-se utilizar tranquilamente a regressão linear e toda a inferência clássica. Isto é mais eficiente do que a regressão quantílica, mesmo para a estimação da mediana da distribuição.

Mas é neste ponto que começa a confusão entre função de escore e função de perda.

Se os dados tem distribuição lognormal, a utilização da função de perda quadrática (MSE) é mais eficiente na estimação, mas isto não significa que a previsão de valores deva ser feita, necessariamente, com a média da distribuição. Pode ser que se deseje, para previsão, minimizar o MAE. Esta, em última análise, é uma decisão do cliente. Ele é quem deve definir se, no seu caso, é mais importante que se erre menos em termos absolutos, o que penaliza da mesma forma os valores mais extremos e menos extremos, ou se é mais importante que se erre menos em termos quadráticos, o que acarreta em dar um peso maior aos erros mais extremos.

Se o cliente opta pela minimização do MSE, o analista deve adotar a média da distribuição lognormal para efetuação das previsões. Já se a opção do cliente for pelo menor MAE, o analista deve optar pela mediana da distribuição.

E quanto às outras funções de escore? Existe uma infinidade delas.

Com algumas pode ser que não haja solução teórica, apenas numérica. O cliente pode optar por qualquer uma delas, cabendo ao analista escolher o método estatístico mais adequado.

Dentre as funções de escore que mais estão me intrigando no momento está a MAPE, (mean absolute percentage error), ou erro médio percentual absoluto. Esta consiste em uma minimização do erro percentual absoluto cometido na previsão das observações.

Ainda não está claro pra mim se e quando poderia ser melhor utilizar esta função de escore em detrimento das outras, mais tradicionais, como a MSE e a MAE.

(Aliás, o problema de delegar ao cliente a escolha da função de escore mais apropriada é este: o cliente sabe escolher qual a melhor função para ele?)

Esta resposta, no entanto, sugere que, ao menos para dados lognormais, a moda da distribuição minimiza o MAPE.

Apenas para ilustrar, fiz umas simulações com dados lognormais: imagine-se que sejam dados de lotes com área variando uniformemente de 360 a 600 \(m^2\), com valores unitários seguindo uma distribuição lognormal.

set.seed(1)
dados <- data.frame(
  Area = runif(100, min = 360, max = 600)
)
dados$LVU <- 12 - .0075*dados$Area + rnorm(100, mean = 0, sd = .25)
dados$VU <- exp(dados$LVU)

Figure 1: Valores Unitários de Lotes com distribuição lognormal.

Se ajustamos um modelo de regressão linear, sem a tranformação dos dados, obtemos:

linFit <- lm(VU~Area, dados)
summary(linFit)

## 
## Call:
## lm(formula = VU ~ Area, data = dados)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2792.9  -770.8  -196.8   495.8  6287.5 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 21575.896   1126.245   19.16   <2e-16 ***
## Area          -34.442      2.306  -14.94   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1473 on 98 degrees of freedom
## Multiple R-squared:  0.6948, Adjusted R-squared:  0.6917 
## F-statistic: 223.2 on 1 and 98 DF,  p-value: < 2.2e-16

Percebam que trata-se de um modelo razoavelmente bom, porém, nem são necessárias maiores análises para saber que as hipóteses da inferência clássica não se verificam: basta a análise da Figura 2, em que diversos pontos se localizam bem distantes da reta de regressão, o que sugere uma distribuição com cauda longa.

Figure 2: Modelo linear.

Para previsão de valores centrais, no entanto, o modelo linear não deixa muito a desejar, como pode-se ver na Figura 3:

powerPlot(linFit, axis = "inverted")

Figure 3: Poder de Predição para o modelo linear.

No entanto, se ajustamos um modelo com a variável dependente transformada, temos um modelo melhor ajustado (na escala logarítmica), onde se verificam os pressupostos básicos da inferência clássica.

loglinFit <- lm(log(VU)~Area, dados)
s <- summary(loglinFit)
s

## 
## Call:
## lm(formula = log(VU) ~ Area, data = dados)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.46244 -0.14056 -0.02177  0.13107  0.62915 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 11.8380400  0.1798516   65.82   <2e-16 ***
## Area        -0.0071746  0.0003682  -19.49   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2353 on 98 degrees of freedom
## Multiple R-squared:  0.7949, Adjusted R-squared:  0.7928 
## F-statistic: 379.7 on 1 and 98 DF,  p-value: < 2.2e-16

Mas, na Engenharia de Avaliações, geralmente o que importa é a previsão de valores na escala original. Os gráficos das Figuras 4 a 6 abaixo mostram como a previsão de valores com a média, moda e mediana minimizam o MSE, o MAPE e o MAE, respectivamente. Nestes gráficos, a linha vermelha é a bissetriz e a linha azul é a regressão dos valores previstos aos valores observados.

yhatMedia <- exp(fitted(loglinFit) + s$sigma^2/2)
powerPlot(dados$VU, yhatMedia, axis = "inverted")

Figure 4: Poder de predição com a média.

yhatMediana <- exp(fitted(loglinFit))
powerPlot(dados$VU, yhatMediana, axis = "inverted")

Figure 5: Poder de predição com a mediana.

yhatModa <- exp(fitted(loglinFit) - s$sigma^2)
powerPlot(dados$VU, yhatModa, axis = "inverted")

Figure 6: Poder de predição com a moda.

Vamos ver como nos sairíamos se optássemos pela mediana, porém estimando-a através da regressão quantílica.

library(quantreg)
medFit <- rq(VU~Area, data = dados, tau = .5)
powerPlot(medFit, axis = "inverted")

Figure 7: Poder de predição com a mediana, obtida através da regressão quantílica.

Conforme esperado, percebe-se que a estimação com a regressão linear é mais eficiente para este tipo de dados. A estimação da mediana à partir da regressão linear minimizou o erro médio absoluto mais eficientemente do que a própria regressão quantílica, que é um método que busca minimizar o MAE diretamente.