Rede Neural para Regressão

Nesta primeira aula prática vamos trabalhar com um banco de dados sobre o custo de segurados para uma seguradora. O objetivo é ajustar um modelo de redes neurais Perceptron de camada única a fim de prever o custo de cada segurado para a seguradora a partir da observação de alguns atributos.

A base de dados pode ser baixada pelo link https://www.kaggle.com/datasets/teertha/ushealthinsurancedataset. Após baixar a base para o seu computador você terá um arquivo insurance.csv. Salve ele no mesmo diretório do seu arquivo R.

O primeiro passo é usar o comando getwd() para verificar se o diretório corrente é o mesmo onde a base foi salva. Caso não seja use o comando setwd() para modificar o diretório corrente para o local onde o arquivo de dados se encontra. Por exemplo:

setwd("D:/Dropbox/Jessica/Trabalho/UFF/Disciplinas/Cursos/Estatistica/GET00113 - Redes Neurais/Aulas Praticas/Atividade 1")

Leitura

Para ler o arquivo insurance.csv no R use o comando a seguir.

base = read.csv2("insurance.csv",sep = ",",dec=".")

Para transformar a base em um tibble:

library(tidyverse)
base = tibble::as_tibble(base)
base
output
## # A tibble: 1,338 × 7
##      age sex      bmi children smoker region    charges
##    <int> <chr>  <dbl>    <int> <chr>  <chr>       <dbl>
##  1    19 female  27.9        0 yes    southwest  16885.
##  2    18 male    33.8        1 no     southeast   1726.
##  3    28 male    33          3 no     southeast   4449.
##  4    33 male    22.7        0 no     northwest  21984.
##  5    32 male    28.9        0 no     northwest   3867.
##  6    31 female  25.7        0 no     southeast   3757.
##  7    46 female  33.4        1 no     southeast   8241.
##  8    37 female  27.7        3 no     northwest   7282.
##  9    37 male    29.8        2 no     northeast   6406.
## 10    60 female  25.8        0 no     northwest  28923.
## # … with 1,328 more rows

A base apresenta as seguintes variáveis: age, sex, bmi, children, smoker, region, charges. São numéricas quantitativas as variáveis age,bmi, children e charges. As outras são variáveis qualitativas: sex, smoker, region.

Divisão entre treino e teste

O próximo passo é separar a base em treino (70%) e teste (30%). A base de treino será usada para fazer análises estatísticas e ajustar os modelos. A base de teste será usada mais adiante, para medir a qualidade do ajuste fora da amostra.

Para fazer essa partição vamos usar a função createDataPartition o pacote caret.

library(caret)

A função createDataPartition é usada para criar uma partição nos índices, depois dividimos a base em dois grupos de acordo com os índices selecionados. Por isso o primeiro argumento da função são os números de 1 até o índice da última linha da base que queremos partir.

O argumento p indica a porcentagem da base que será selecionada.

set.seed(123456789)
n = dim(base)[1]
indices_treino = createDataPartition(1:n,p=0.7)[[1]]
base_treino = base |> slice(indices_treino)
base_teste  = base[-indices_treino,]

Padronização das variáveis.

As variáveis quantitativas muitas vezes estão em escalas bem diferentes, o que atrapalha o ajuste do modelo. Esse problema é facilmente resolvido com uma transformação dessas variáveis antes do ajuste de forma que todas fiquem com média zero e variância 1.

A padronização é feita considerando os dados treino, porque apenas estes vamos considerar conhecidos para o treinamento. É importante guardar os parâmetros da transformação, porque eles serão usados novamente mais pra frente.

media = apply(base_treino |> select(age,bmi,charges),2,mean)
desvp = apply(base_treino |> select(age,bmi,charges),2,sd)
base_treino_ = base_treino |> mutate(age = (age - media["age"])/(desvp["age"]),
                                     bmi = (bmi - media["bmi"])/(desvp["bmi"]),
                                     charges = (charges - media["charges"])/(desvp["charges"])
                                      )

Agora que os dados estão padronizados precisamos verificar se existe grande correlação entre as covariáveis. Caso exista algumas terão que ser excluídas da base. A ideia é que as covariáveis usadas em qualquer modelo não tenham correlação maior que \(0,75\). O comando abaixo apresenta uma matriz com todas as informações necessárias para uma análise de correlação.

library(ggplot2)
library(GGally)
ggpairs(base_treino)

Perceptron 1

Problema 1: Ajustar um modelo de redes neurais perceptron de camada única, para prever o custo (charge) de um segurado do sexo feminino que mora em northeast em função da sua idade (age) e índice de massa corporal (bmi).

Perguntas:

  • Este é um problema de regressão, classificação ou agrupamento?
output
## Regressão
  • Quem é a variável alvo nesta aplicação?
output
## O custo de cada segurado, `charge`
  • Quais os atributos que vamos usar?
output
## A idade (age) e o índice de massa corporal (bmi).
  • Todas as linhas da base original entram na base a ser analisada para esta aplicação?
output
## Não, serão considerados somente as instâncias que representam mulheres da região de northeast.

A base de dados

Segue o comando para cirar a base a ser usada, com as linhas e colunas da base bruta que nos interessa.

base_treino_1_ = base_treino_ |> filter(sex=="female", region=="northeast") |> select(charges,age,bmi)
base_treino_1_
output
## # A tibble: 106 × 3
##    charges    age     bmi
##      <dbl>  <dbl>   <dbl>
##  1 -0.219   0.895  0.0113
##  2 -0.0164  1.46   0.856 
##  3  2.02   -0.372  0.196 
##  4 -0.836  -1.50   1.29  
##  5 -0.934  -1.50   0.795 
##  6  0.0151  1.32   0.180 
##  7  0.660  -1.50  -0.0962
##  8 -0.836  -0.935 -0.311 
##  9  0.0158  1.53  -1.40  
## 10 -0.591  -0.513  1.04  
## # … with 96 more rows

O Modelo de Redes Neurais

O modelo perceptron de camada única terá 2 variáveis de entrada: \(x_1=\)age e \(x_2 =\)bmi .A variável de saida, aquela que queremos prever, será \(y=\)charges. Esse modelo tem 3 constantes para serem estimadas: \(w_1\), \(w_2\) e \(\Theta\).

Figura: Rede Perceptron camada única para Regressão

Este ainda é um exemplo muito inicial de redes neurais, por isso ainda é possível escrever a a função que define a previsão da variável alvo \(\hat{y}\) em termos das variáveis independentes \(x_1\) e \(x_2\) e também dos parâmetros estimados \(\hat{w}_1x_1\), \(\hat{w}_2x_2\) e \(\hat\Theta\).

\[ \hat{y} = g \left( \hat{w}_1x_1 + \hat{w}_2x_2 + \hat\Theta \right) \]

Ajuste do modelo

Quando dizemos que vamos ajustar ou treinar o modelo de redes neurais, na prática, o que vamos fazer é encontrar estimativas para os parâmetros desconhecidos. No caso desta aplicação, estimativas para \(w_1\), \(w_2\) e \(\Theta\).

Como essas estimativas são encontradas será assunto para aulas futuras, mas vale nesse momento dizer que a ideia principal é resolver um problema de otiumização que busca minimizar alguma medida de erro no ajuste do modelo.

Para encontrar as estimativas dos parâmetros desconhecidos da rede perceptron vamos usar a função neuralnet do pacote neuralnet.

library(neuralnet)

Como é a nossa primeira aula prática vamos usar os parâmetros padrões da função sempre que possível, nas próximas exploramos melhor as suas variações. A função neuralnet é usada da seguinte maneira:

perceptron1 = neuralnet(formula = charges ~ age + bmi,
                       data = base_treino_1_,
                       hidden=0)

Assim é criado um objeto perceptron1 da classe nn que guarda as seguintes informações sobre o ajuste da rede neural.

names(perceptron1)
output
##  [1] "call"                "response"            "covariate"          
##  [4] "model.list"          "err.fct"             "act.fct"            
##  [7] "linear.output"       "data"                "exclude"            
## [10] "net.result"          "weights"             "generalized.weights"
## [13] "startweights"        "result.matrix"

Os parâmetros estimados estão guardados em perceptron1$weights e é possível desenhar a rede criada com o comando plot(perceptron1, rep="best").

perceptron1$weights
output
## [[1]]
## [[1]][[1]]
##            [,1]
## [1,] 0.03992513
## [2,] 0.28898887
## [3,] 0.16109294
plot(perceptron1, rep="best")
plot

Previsão na base de treino

O próprio objeto perceptron1 criado pela função neuralnet guarda os valores previstos pelo modelo para os dados usados no ajuste. Estes valores estão em perceptron1$net.result, que é uma lista de tamanho 1. Esta lista guarda uma matriz com uma única coluna, que são ps valores previstos para os dados de treinamento. Lembre-se que para comparar com os dados reais é preciso retornar para a escala inicial.

prev_treino_1_ = perceptron1$net.result[[1]][,1]
prev_treino_1 = prev_treino_1_*desvp["charges"] + media["charges"]

Uma vez tendo os valores previstos na base de treino é possível calcular o erro médio quadrático (MSE) para a base de treino. Esta será uma medida de qualidade de ajuste da rede perceptron.

n = dim(base_treino_1_)[1]
charges_treino_1 = 
  (base_treino |> filter(sex=="female", region=="northeast"))$charges
MSE_treino_1 = (1/n)*sum((charges_treino_1 - prev_treino_1)^2)
MSE_treino_1
output
## [1] 98396281

Previsão para novos dados - base de teste

Um dos objetivos do modelo pode ser prever o custo de novos segurados, por exemplo. Para medir a sua capacidade de previsão para novas observações vamos usar a base de teste, que contém dados não usados no ajuste do modelo.

E como fazer essa previsão? Temos que reproduzir tudo que foi feito na base de treino no novo dado.

O primeiro passo desta reprodução é a padronização realizada nos dados. Neste momento vamos precisar das constantes que usamos no dado de treino e usar as mesmas nos dados de teste.

base_teste_ = base_teste |> mutate(age = (age - media["age"])/(desvp["age"]),
                                   bmi = (bmi - media["bmi"])/(desvp["bmi"]),
                                   charges = (charges - media["charges"])/(desvp["charges"])
)

Agora que já temos a base de teste padronizada, vamos selecionar as linhas e colunas de interesse.

base_teste_1_ = base_teste_ |> 
  filter(sex=="female", region=="northeast") |> select(charges,age,bmi)

A previsão para dados que não foram usados no ajuste é feita pela função compute.

prev_teste_1_ = (perceptron1 |> compute(base_teste_1_))$net.result[,1]
prev_teste_1 = prev_teste_1_*desvp["charges"] + media["charges"]

Agora que temos previsões para dados fora da amostra podemos medir a qualidade do ajuste fora da amostra usando o erro médio quadrático na base de teste.

n = dim(base_teste_1_)[1]
charges_teste_1 = 
  (base_teste |> filter(sex=="female", region=="northeast"))$charges
MSE_teste_1 = (1/n)*sum((charges_teste_1 - prev_teste_1)^2)
MSE_teste_1
output
## [1] 93295861
barplot(
  cbind(c(MSE_treino_1,MSE_teste_1)),
  beside = T,
  names.arg = c("perceptron1"),
  main = "Erro Médio Quadrático (MSE)",
  col=c("tomato","skyblue"),
  legend.text = c("treino","teste"))
text(1.5,MSE_treino_1*.9,format(MSE_treino_1,digits = 3,scientific = T))
text(2.5,MSE_teste_1*.9,format(MSE_teste_1,digits = 3,scientific = T))
plot

Perceptron 2

A segunda aplicação será bem pareceda com a primeira, a diferença é que vamos incluir nas variáveis preditivas a variável que indica se um segurado é ou não fumante.

Problema 2: Ajustar um modelo de redes neurais perceptron de camada única para prever o custo (charge) de um segurado do sexo feminino que mora em northeast em função da sua idade (age), do índice de massa corporal (bmi) e do hábito de fumar (smoker).

A base de dados

A única diferença da base de dados desta aplicação para a anterior é que nesta também será incluída a variável smoker.

base_treino_2_ = base_treino_ |> filter(sex=="female", region=="northeast") |> select(charges,age,bmi,smoker)
base_treino_2_
output
## # A tibble: 106 × 4
##    charges    age     bmi smoker
##      <dbl>  <dbl>   <dbl> <chr> 
##  1 -0.219   0.895  0.0113 no    
##  2 -0.0164  1.46   0.856  no    
##  3  2.02   -0.372  0.196  yes   
##  4 -0.836  -1.50   1.29   no    
##  5 -0.934  -1.50   0.795  no    
##  6  0.0151  1.32   0.180  no    
##  7  0.660  -1.50  -0.0962 no    
##  8 -0.836  -0.935 -0.311  no    
##  9  0.0158  1.53  -1.40   no    
## 10 -0.591  -0.513  1.04   no    
## # … with 96 more rows

O Modelo de Redes Neurais

Este modelo perceptron de camada única terá 3 variáveis de entrada: \(x_1=\)age, \(x_2 =\)bmi e \(x_3=\)smoker .A variável de saída, aquela que queremos prever, continua sendo \(y=\)charges. Esse modelo tem 4 constantes para serem estimadas: \(w_1\), \(w_2\), \(w_3\) e \(\Theta\).

Neste momento temos também que escolher a função de ativação do único nó da rede neural. Como a aplicação é de regressão, vamos continuar com a função identidade.

Figura: Rede Perceptron camada única para Regressão

Tratamento para covariáveis categóricas

Veja que neste exemplo \(x_3=\)smoker é uma variável categórica. Para ela entrar no modelo é preciso convertê-la em variável indicadora.

\[ \tilde{x}_3 = \left\{ \begin{array}{ll} 1 & , \hbox{se } x_3 {= yes}\\ 0 & , \hbox{se } x_3 {= no}\\ \end{array} \right. \]

Depois de criada essa nova variável \(\tilde{x}_3\) ela substitui a variável categórica \(x_3\) e a função que define a previsão da variável resposta (\(\hat{y}\)) em termos das covariáveis fica assim:

\[ \hat{y} = g\left(w_1x_1 + w_2x_2 + w_3\tilde{x}_3 + \Theta \right) \]

Esse tratamento pode ser feito manualmente ou, o que é mais simples, podemos usar a função model.matrix que já transforma todas as variáveis categóricas em binárias “dummy”. A função neuralnet é um exemplo raro onde o R não faz essa transformação automaticamente.

matriz_treino_2_ <- model.matrix( ~ charges + age + bmi + smoker, data = base_treino_2_)
head(matriz_treino_2_)
output
##   (Intercept)     charges        age        bmi smokeryes
## 1           1 -0.21899526  0.8949865 0.01133022         0
## 2           1 -0.01643970  1.4581713 0.85595675         0
## 3           1  2.02227148 -0.3721793 0.19561237         1
## 4           1 -0.83577939 -1.4985489 1.28594844         0
## 5           1 -0.93426399 -1.4985489 0.79452937         0
## 6           1  0.01509283  1.3173751 0.18025553         0

Veja que foi criada uma nova variável smokeryes e a variável smoker foi eliminada. A base matriz_treino_2 será aquela que vamos usar para treinar o modelo.

Ajuste do modelo

Nesta etapa vamos encontrar estimativas para os parâmetros desconhecidos.

perceptron2 = neuralnet(formula = charges ~ age + bmi + smokeryes,
                       data = matriz_treino_2_,
                       hidden=0)
names(perceptron2)
output
##  [1] "call"                "response"            "covariate"          
##  [4] "model.list"          "err.fct"             "act.fct"            
##  [7] "linear.output"       "data"                "exclude"            
## [10] "net.result"          "weights"             "generalized.weights"
## [13] "startweights"        "result.matrix"

Previsão na base de treino

Vejamos os valores previstos para os custos dos indivíduos na base de treino.

prev_treino_2_ = perceptron2$net.result[[1]][,1]
prev_treino_2 = prev_treino_2_*desvp["charges"] + media["charges"]

Agora o cálculo do MSE.

n = dim(matriz_treino_2_)[1]
charges_treino_2 = 
  (base_treino |> filter(sex=="female", region=="northeast"))$charges
MSE_treino_2 = (1/n)*sum((charges_treino_2 - prev_treino_2)^2)
MSE_treino_2
output
## [1] 38303263

Previsão para novos dados - base de teste

Para realizar as previsões na base de teste temos que reproduzir tudo que foi feito na base de treino. A padronização dos dados já foi feita e o objeto base_teste_ guarda a base de teste padronizada. Seguimos então para a seleção das linhas e colunas de interesse.

base_teste_2_ = base_teste_ |> 
  filter(sex=="female", region=="northeast") |> select(charges,age,bmi,smoker)

No caso da base conter variáveis categóricas, antes de realizar a previsão é necessário fazer a transformação das variáveis categóricas em indicadoras também na base de teste.

matriz_teste_2_ <- model.matrix( ~ charges + age + bmi + smoker, 
                                 data = base_teste_2_)

A previsão para dados que não foram usados no ajuste é feita pela função compute.

prev_teste_2_ = (perceptron2 |> compute(matriz_teste_2_))$net.result[,1]
prev_teste_2 = prev_teste_2_*desvp["charges"] + media["charges"]

Agora o erro médio quadrático na base de teste.

n = dim(matriz_teste_2_)[1]
charges_teste_2 = 
  (base_teste |> filter(sex=="female", region=="northeast"))$charges
MSE_teste_2 = (1/n)*sum((charges_teste_2 - prev_teste_2)^2)
MSE_teste_2
output
## [1] 42349147
barplot(
  cbind(c(MSE_treino_1,MSE_teste_1),c(MSE_treino_2,MSE_teste_2)),
  beside = T,
  names.arg = c("perceptron1","perceptron2"),
  main = "Erro Médio Quadrático (MSE)",
  col=c("tomato","skyblue"),
  legend.text = c("treino","teste"))
text(1.5,MSE_treino_1*.9,format(MSE_treino_1,digits = 3,scientific = T))
text(2.5,MSE_teste_1*.9,format(MSE_teste_1,digits = 3,scientific = T))
text(4.5,MSE_treino_2*.9,format(MSE_treino_2,digits = 3,scientific = T))
text(5.5,MSE_teste_2*.9,format(MSE_teste_2,digits = 3,scientific = T))
plot

Perceptron 3

A terceira aplicação da aula de hoje será mais uma generalização da primeira. Vamos agora não nos restingir aos segurado do sexo feminino que moram em northeast.

Problema 3: Criar um modelo de redes neurais para prever o custo (charge) de um segurado qualquer em função da sua idade, gênero, índice de massa corporal, número de filhos, hábito de fumar e região onde mora.

A base de dados

Nessa aplicação usaremos a base completa, que já está guardada nos objetos base_treino e base_teste.

O Modelo de Redes Neurais

Este modelo perceptron de camada única terá 6 variáveis de entrada: \(x_1=\)age, \(x_3=\)sex, \(x_3 =\)bmi, \(x_4=\)children, \(x_5=\)smoker e \(x_6=\)region .A variável de saida, aquela que queremos prever, continua sendo \(y=\)charges. Esse modelo tem 7 constantes para serem estimadas: \(w_1\), \(w_2\), \(w_3\), \(w_4\), \(w_5\), \(w_6\) e \(\Theta\).

Figura: Rede Perceptron camada única para Regressão

Tratamento para covariáveis categóricas

Veja que agora temos mais variáveis categóricas, são elas: sex, smoker e region. Todas serão trasformadas em variáveis indicadoras.

Para cada variável categórica com 2 classes é criada 1 indicadora. Mas se ela tiver mais classes, serão criadas mais variáveis para conseguir contemplar todas as possibilidades. A variável region, por exemplo, tem 4 classes e serão criadas 3 variáveis indicadoras para substituí-la.

matriz_treino_3_ <- model.matrix( ~  age + sex + bmi + children + smoker + region + charges, 
                                  data = base_treino_)
head(matriz_treino_3_)
output
##   (Intercept)        age sexmale        bmi children smokeryes regionnorthwest
## 1           1 -1.4281508       0 -0.4542247        0         1               0
## 2           1 -1.4985489       1  0.4946667        1         0               0
## 3           1 -0.7945679       1  0.3701955        3         0               0
## 4           1 -0.5129755       1 -0.2958067        0         0               1
## 5           1 -0.5833736       0 -0.8033909        0         0               0
## 6           1 -0.1609850       0 -0.4800889        3         0               1
##   regionsoutheast regionsouthwest    charges
## 1               0               1  0.2881276
## 2               1               0 -0.9747148
## 3               1               0 -0.7478011
## 4               0               0 -0.7963349
## 5               1               0 -0.8055178
## 6               0               0 -0.5118794

Ajuste do modelo

perceptron3 = neuralnet(formula = charges ~ age + sexmale + bmi + children + smokeryes  + regionnorthwest +regionsoutheast + regionsouthwest,
                       data = matriz_treino_3_,
                       hidden=0)
plot(perceptron3, rep="best")
plot

Previsão na base de treino

Vejamos os valores previstos para os custos dos indivíduos na base de treino.

prev_treino_3_ = perceptron3$net.result[[1]][,1]
prev_treino_3 = prev_treino_3_*desvp["charges"] + media["charges"]

Agora o cálculo do MSE.

n = dim(matriz_treino_3_)[1]
charges_treino_3 = base_treino$charges
MSE_treino_3 = (1/n)*sum((charges_treino_3 - prev_treino_3)^2)
MSE_treino_3
output
## [1] 37921827

Previsão para novos dados - base de teste

Primeiro o tratamento das variáveis categóricas.

matriz_teste_3_ <- model.matrix( ~  age + sex + bmi + children + smoker + region + charges, 
                                 data = base_teste_)

A previsão para dados que não foram usados no ajuste é feita pela função compute.

prev_teste_3_ = (perceptron3 |> compute(matriz_teste_3_))$net.result[,1]
prev_teste_3 = prev_teste_3_*desvp["charges"] + media["charges"]

Agora o erro médio quadrático na base de teste.

n = dim(matriz_teste_3_)[1]
charges_teste_3 = base_teste$charges
MSE_teste_3 = (1/n)*sum((charges_teste_3 - prev_teste_3)^2)
MSE_teste_3
output
## [1] 33551342
barplot(
  cbind(c(MSE_treino_1,MSE_teste_1),c(MSE_treino_2,MSE_teste_2),c(MSE_treino_3,MSE_teste_3)),
  beside = T,
  names.arg = c("perceptron1","perceptron2","perceptron3"),
  main = "Erro Médio Quadrático (MSE)",
  col=c("tomato","skyblue"),
  legend.text = c("treino","teste"))
text(1.5,MSE_treino_1*.9,format(MSE_treino_1,digits = 3,scientific = T))
text(2.5,MSE_teste_1*.9,format(MSE_teste_1,digits = 3,scientific = T))
text(4.5,MSE_treino_2*.9,format(MSE_treino_2,digits = 3,scientific = T))
text(5.5,MSE_teste_2*.9,format(MSE_teste_2,digits = 3,scientific = T))
text(7.5,MSE_treino_3*.9,format(MSE_treino_3,digits = 3,scientific = T))
text(8.5,MSE_teste_3*.9,format(MSE_teste_3,digits = 3,scientific = T))
plot