Regressão Linear e Regressão Logística

lm e glm

Jose Storopoli https://scholar.google.com/citations?user=xGU7H1QAAAAJ&hl=en (UNINOVE)https://www.uninove.br
April 19, 2021

Stats

Figure 1: Stats

Regressão Linear

Regressão linear permite com que você use uma ou mais variáveis discretas ou contínuas como variáveis independentes e mensurar o poder de associação com a variável dependente, que deve ser contínua.

Interpretações

Para compreender regressão linear podemos usar de três interpretações distintas mas complementares:

Interpretação Geométrica

Imagine que seus dados são pontos que vivem em um espaço multidimensional. A regressão é uma técnica para encontrar a melhor reta1 entre o conjunto de dados levando em conta todas as observações.

Isto é valido para qualquer espaço multidimensional, até para além de 3-D. Vamos mostrar um exemplo em 2-D da relação entre x e y, mas isto poder ser estendido para a relação x1, x2, … e y.

Uma relação entre variáveis representada por uma reta de tendência

Figure 2: Uma relação entre variáveis representada por uma reta de tendência

Vejam que regressão linear usando apenas uma variável dependente e uma variável independente é a mesma coisa que que correlação.

Interpretação Matemática

A interpretação matemática é vista como uma otimização: encontrar a melhor reta entre os pontos que minimiza o erro quadrático médio (mean squared error – MSE). Ao escolhermos a melhor reta, devemos escolher a melhor reta que minimiza as distâncias entre os pontos, sendo que podemos errar para mais ou para menos. Para evitarmos que os erros se cancelem, precisamos eliminar o sinal negativo de alguns erros e convertê-los para valores positivos. Para isso, pegamos todas os erros (diferenças entre o valor previsto pela reta e o valor verdadeiro) e elevamos ao quadrado (assim todo número negativo se tornará positivo e todo positivo se manterá positivo). Portanto, a regressão se torna a busca do menor valor de uma função erro (MSE).

A melhor reta que minimiza a distância dos erros

Figure 3: A melhor reta que minimiza a distância dos erros

Interpretação Estatística

A regressão linear usando uma única variável independente contínua se torna exatamente uma correlação. Agora quando empregamos mais de uma variável independente, a interpretação da regressão se torna: “O efeito de X em Y mantendo Z fixo”. Isto quer dizer que a regressão linear controla os efeitos das diferentes variáveis independentes ao calcular o efeito de uma certa variável independente. Esta é o que chamamos de interpretação estatística da regressão linear.

Por exemplo, digamos que você esteja em busca dos fatores que acarretam ataque cardíaco. Você coleta dados de pessoas que quantifiquem as seguintes variáveis: sono, stress, tabagismo, sedentarismo, entre outros… A regressão te permite mensurar o efeito de qualquer uma dessas variáveis na prevalência de ataque cardíaco controlando para outros efeitos. Em outras palavras, é possível mensurar o efeito de stress em ataque cardíaco, mantendo fixo os efeitos de sono, tabagismo, sedentarismo, etc… Isso permite você isolar o efeito de uma variável sem deixar que outras variáveis a influenciem na mensuração da sua relação com a variável dependente (no nosso caso: ataque cardíaco).

Exemplo - Score de QI de crianças

Para o nosso exemplo, usarei um dataset famoso chamado kidiq que está incluído no diretório datasets/. São dados de uma survey de mulheres adultas norte-americanas e seus respectivos filhos. Datado de 2007 possui 434 observações e 4 variáveis:

library(readr)
kidiq <- read_csv("datasets/kidiq.csv", col_types = "_didi")

Como especificar um modelo em R usando a sintaxe de “formula”

Podemos espeficiar modelos usando a sintaxe de formula:

y ~ x1 + x2 + ...
kidiq_1 <- lm(kid_score ~ mom_iq, data = kidiq)
summary(kidiq_1)

Call:
lm(formula = kid_score ~ mom_iq, data = kidiq)

Residuals:
   Min     1Q Median     3Q    Max 
-56.75 -12.07   2.22  11.71  47.69 

Coefficients:
            Estimate Std. Error t value             Pr(>|t|)    
(Intercept)  25.7998     5.9174    4.36             0.000016 ***
mom_iq        0.6100     0.0585   10.42 < 0.0000000000000002 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 18 on 432 degrees of freedom
Multiple R-squared:  0.201, Adjusted R-squared:  0.199 
F-statistic:  109 on 1 and 432 DF,  p-value: <0.0000000000000002
kidiq_2 <- lm(kid_score ~ mom_iq + mom_hs, data = kidiq)
summary(kidiq_2)

Call:
lm(formula = kid_score ~ mom_iq + mom_hs, data = kidiq)

Residuals:
   Min     1Q Median     3Q    Max 
 -52.9  -12.7    2.4   11.4   49.5 

Coefficients:
            Estimate Std. Error t value             Pr(>|t|)    
(Intercept)  25.7315     5.8752    4.38             0.000015 ***
mom_iq        0.5639     0.0606    9.31 < 0.0000000000000002 ***
mom_hs        5.9501     2.2118    2.69               0.0074 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 18 on 431 degrees of freedom
Multiple R-squared:  0.214, Adjusted R-squared:  0.21 
F-statistic: 58.7 on 2 and 431 DF,  p-value: <0.0000000000000002
kidiq_3 <- lm(kid_score ~ mom_iq * mom_hs, data = kidiq)
summary(kidiq_3)

Call:
lm(formula = kid_score ~ mom_iq * mom_hs, data = kidiq)

Residuals:
   Min     1Q Median     3Q    Max 
-52.09 -11.33   2.07  11.66  43.88 

Coefficients:
              Estimate Std. Error t value      Pr(>|t|)    
(Intercept)    -11.482     13.758   -0.83        0.4044    
mom_iq           0.969      0.148    6.53 0.00000000018 ***
mom_hs          51.268     15.338    3.34        0.0009 ***
mom_iq:mom_hs   -0.484      0.162   -2.99        0.0030 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 18 on 430 degrees of freedom
Multiple R-squared:  0.23,  Adjusted R-squared:  0.225 
F-statistic: 42.8 on 3 and 430 DF,  p-value: <0.0000000000000002

Se você quiser plotar modelos de regressão há o pacote {ggeffect}.

Regressão Logística

Uma regressão logística se comporta exatamente como um modelo linear: faz uma predição simplesmente computando uma soma ponderada das variáveis independentes, mais uma constante. Porém ao invés de retornar um valor contínuo, como a regressão linear, retorna a função logística desse valor.

\[\operatorname{Logística}(x) = \frac{1}{1 + e^{(-x)}}\]

A função logística é uma gambiarra transformação que pega qualquer valor entre menos infinito \(-\infty\) e mais infinito \(+\infty\) e transforma em um valor entre 0 e 1. Veja na figura 4 uma representação gráfica da função logística.

library(dplyr)
library(ggplot2)
tibble(
  x = seq(-10, 10, length.out = 100),
  logit = 1 / (1 + exp(-x))) %>%
  ggplot(aes(x, logit)) +
  geom_line()
Função Logística

Figure 4: Função Logística

Ou seja, a função logística é a candidata ideal para quando precisamos converter algo contínuo sem restrições para algo contínuo restrito entre 0 e 1. Por isso ela é usada quando precisamos que um modelo tenha como variável dependente uma probabilidade (lembrando que qualquer numero real entre 0 e 1 é uma probabilidade válida). No caso de uma variável dependente binária, podemos usar essa probabilidade como a chance da variável dependente tomar valor de 0 ou de 1.

Exemplo - Propensão a mudar de poço de água contaminado

Para exemplo, usaremos um dataset chamado wells que está incluído no diretório datasets/. É uma survey com 3.200 residentes de uma pequena área de Bangladesh na qual os lençóis freáticos estão contaminados por arsênico. Respondentes com altos níveis de arsênico nos seus poços foram encorajados para trocar a sua fonte de água para uma níveis seguros de arsênico.

Possui as seguintes variáveis:

wells <- read_csv("datasets/wells.csv", col_types = "iddii")
wells1 <- glm(switch ~ arsenic + dist + educ + assoc,
              data = wells,
              family = binomial)
summary(wells1)

Call:
glm(formula = switch ~ arsenic + dist + educ + assoc, family = binomial, 
    data = wells)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-2.594  -1.198   0.754   1.063   1.674  

Coefficients:
            Estimate Std. Error z value             Pr(>|z|)    
(Intercept) -0.15671    0.09960   -1.57                 0.12    
arsenic      0.46702    0.04160   11.23 < 0.0000000000000002 ***
dist        -0.00896    0.00105   -8.57 < 0.0000000000000002 ***
educ         0.04245    0.00959    4.43            0.0000095 ***
assoc       -0.12430    0.07697   -1.61                 0.11    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 4118.1  on 3019  degrees of freedom
Residual deviance: 3907.8  on 3015  degrees of freedom
AIC: 3918

Number of Fisher Scoring iterations: 4

Vamos pegar o exp() dos coeficientes:

exp(wells1$coefficients)
(Intercept)     arsenic        dist        educ       assoc 
       0.85        1.60        0.99        1.04        0.88 
library(ggplot2)
wells %>% 
  ggplot(aes(arsenic)) +
  geom_histogram()

Regressão de Poisson - Dados de Contagem

Uma regressão de Poisson se comporta exatamente como um modelo linear: faz uma predição simplesmente computando uma soma ponderada das variáveis independentes \(\mathbf{X}\) pelos coeficientes estimados \(\boldsymbol{\beta}\), mais uma constante \(\alpha\). Porém ao invés de retornar um valor contínuo \(y\), como a regressão linear, retorna o logarítmo natural desse valor \(\log(y)\).

\[ \log(y)= \theta_0 + \theta_1 x_1 + \theta_2 x_2 + \dots + \theta_n x_n \]

que é o mesmo que

\[ y = e^{(\theta_0 + \theta_1 x_1 + \theta_2 x_2 + \dots + \theta_n x_n)} \]

A função \(e^x\) é chamada de função exponencial. Veja a figura 5 para uma intuição gráfica da função exponencial:

ggplot(data = tibble(
  x = seq(0, 10, length.out = 100),
  y =  exp(x)
  ),
  aes(x, y)) +
  geom_line(size = 2, color = "steelblue") +
  ylab("Exponencial(x)")
Função Exponencial

Figure 5: Função Exponencial

Regressão de Poisson é usada quando a nossa variável dependente só pode tomar valores positivos, geralmente em contextos de dados de contagem.

Exemplo Poisson - Exterminação de baratas

Para exemplo, usaremos um dataset chamado roaches que está incluído no diretório datasets/. É uma base de dados com 262 observações sobre a eficácia de um sistema de controle de pragas em reduzir o número de baratas (roaches) em apartamentos urbanos.

Possui as seguintes variáveis:

roaches <- read_csv("datasets/roaches.csv", col_types = "idiid")
roaches1 <- glm(
  y ~ roach1 + treatment + senior + exposure2,
  data = roaches,
  family = poisson
)
summary(roaches1)

Call:
glm(formula = y ~ roach1 + treatment + senior + exposure2, family = poisson, 
    data = roaches)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-18.35   -5.29   -4.18    0.29   27.67  

Coefficients:
              Estimate Std. Error z value             Pr(>|z|)    
(Intercept)  2.9665993  0.0433859   68.38 < 0.0000000000000002 ***
roach1       0.0065252  0.0000902   72.38 < 0.0000000000000002 ***
treatment   -0.5148722  0.0246802  -20.86 < 0.0000000000000002 ***
senior      -0.3787607  0.0335552  -11.29 < 0.0000000000000002 ***
exposure2    0.1620458  0.0362256    4.47            0.0000077 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 16475  on 261  degrees of freedom
Residual deviance: 11431  on 257  degrees of freedom
AIC: 12195

Number of Fisher Scoring iterations: 6

Vamos pegar o exp() dos coeficientes:

exp(roaches1$coefficients)
(Intercept)      roach1   treatment      senior   exposure2 
      19.43        1.01        0.60        0.68        1.18 

Como plotar modelos com o {ggeffects}

Podemos usar o pacote [{ggeffects}](https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmodel.html) para *plotar* objetoslmeglm`

$mom_iq


$mom_hs

plot(ggeffect(wells1))
$arsenic


$dist


$educ


$assoc

plot(ggeffect(roaches1))
$roach1


$treatment


$senior


$exposure2

Ambiente

R version 4.1.0 (2021-05-18)
Platform: aarch64-apple-darwin20 (64-bit)
Running under: macOS Big Sur 11.4

Matrix products: default
BLAS:   /Library/Frameworks/R.framework/Versions/4.1-arm64/Resources/lib/libRblas.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.1-arm64/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods  
[7] base     

other attached packages:
 [1] ggeffects_1.1.0   broom_0.7.8       furrr_0.2.3      
 [4] future_1.21.0     purrr_0.3.4       ggridges_0.5.3   
 [7] ggExtra_0.9       gghighlight_0.3.2 ggrepel_0.9.1    
[10] patchwork_1.1.1   forcats_0.5.1     plotly_4.9.4.1   
[13] repurrrsive_1.0.0 ggplot2_3.3.5     stringr_1.4.0    
[16] tidyr_1.1.3       janitor_2.1.0     dplyr_1.0.7      
[19] readr_1.4.0       magrittr_2.0.1    tibble_3.1.2     

loaded via a namespace (and not attached):
 [1] nlme_3.1-152       lubridate_1.7.10   bit64_4.0.5       
 [4] insight_0.14.2     RColorBrewer_1.1-2 httr_1.4.2        
 [7] rprojroot_2.0.2    backports_1.2.1    tools_4.1.0       
[10] bslib_0.2.5.1      sjlabelled_1.1.8   utf8_1.2.1        
[13] R6_2.5.0           DBI_1.1.1          lazyeval_0.2.2    
[16] mgcv_1.8-35        colorspace_2.0-2   withr_2.4.2       
[19] tidyselect_1.1.1   downlit_0.2.1      bit_4.0.4         
[22] compiler_4.1.0     textshaping_0.3.5  cli_3.0.0         
[25] labeling_0.4.2     bookdown_0.22      sass_0.4.0        
[28] scales_1.1.1       systemfonts_1.0.2  digest_0.6.27     
[31] rmarkdown_2.9      jpeg_0.1-8.1       pkgconfig_2.0.3   
[34] htmltools_0.5.1.1  parallelly_1.26.1  dbplyr_2.1.1      
[37] fastmap_1.1.0      highr_0.9          htmlwidgets_1.5.3 
[40] rlang_0.4.11       rstudioapi_0.13    RSQLite_2.2.7     
[43] shiny_1.6.0        jquerylib_0.1.4    farver_2.1.0      
[46] generics_0.1.0     jsonlite_1.7.2     crosstalk_1.1.1   
[49] distill_1.2        Matrix_1.3-4       Rcpp_1.0.6        
[52] munsell_0.5.0      fansi_0.5.0        lifecycle_1.0.0   
[55] stringi_1.6.2      yaml_2.2.1         snakecase_0.11.0  
[58] plyr_1.8.6         grid_4.1.0         blob_1.2.1        
[61] parallel_4.1.0     listenv_0.8.0      promises_1.2.0.1  
[64] crayon_1.4.1       miniUI_0.1.1.1     lattice_0.20-44   
[67] splines_4.1.0      hms_1.1.0          knitr_1.33        
[70] pillar_1.6.1       codetools_0.2-18   glue_1.4.2        
[73] evaluate_0.14      data.table_1.14.0  httpuv_1.6.1      
[76] png_0.1-7          vctrs_0.3.8        gtable_0.3.0      
[79] assertthat_0.2.1   cachem_1.0.5       xfun_0.24         
[82] mime_0.11          xtable_1.8-4       later_1.2.0       
[85] ragg_1.1.3         viridisLite_0.4.0  memoise_2.0.0     
[88] globals_0.14.0     ellipsis_0.3.2    

  1. tecnicamente reta aqui se refere um hiperplano que é subespaço de dimensão \(n-1\) de um espaço de dimensão \(n\). Por exemplo, uma reta é um hiperplano 1-D de uma plano 2-D; um plano 2-D é um hiperplano de um plano 3-D; e assim por diante…↩︎

Corrections

If you see mistakes or want to suggest changes, please create an issue on the source repository.

Reuse

Text and figures are licensed under Creative Commons Attribution CC BY-SA 4.0. Source code is available at https://github.com/storopoli/Linguagem-R, unless otherwise noted. The figures that have been reused from other sources don't fall under this license and can be recognized by a note in their caption: "Figure from ...".

Citation

For attribution, please cite this work as

Storopoli (2021, April 19). Linguagem R: Regressão Linear e Regressão Logística. Retrieved from https://storopoli.io/Linguagem-R/5-Regressao.html

BibTeX citation

@misc{storopoli2021regressaoR,
  author = {Storopoli, Jose},
  title = {Linguagem R: Regressão Linear e Regressão Logística},
  url = {https://storopoli.io/Linguagem-R/5-Regressao.html},
  year = {2021}
}