A regressão logística é um algoritmo de modelagem preditiva usado quando a variável \(Y\) é categórica binária. Ou seja, pode ter apenas dois valores como 1 ou 0. O objetivo é determinar uma equação matemática que pode ser usada para prever a probabilidade do evento 1. Uma vez que a equação é estabelecida, ela pode ser usada para prever o \(Y\) quando apenas os \(X\)s são conhecidos.
Anteriormente, vimos o que é regressão linear e como usá-la para prever variáveis \(Y\) contínuas. Na regressão linear, a variável \(Y\) é sempre uma variável contínua. Se supor que a variável \(Y\) era categórica, você não pode usar o modelo de regressão linear dela.
Um ponto importante a se notar aqui é que \(Y\) pode ter apenas 2 classes e não mais do que isso. Se \(Y\) tiver mais de 2 classes, ela se tornará uma classificação multiclasse e você não pode mais usar a regressão logística para isso.
No entanto, a regressão logística é uma técnica de modelagem preditiva clássica e ainda permanece uma escolha popular para modelar variáveis categóricas binárias. Outra vantagem da regressão logística é que ela calcula pontualmente a probabilidade de previsão de um evento.
Agora vamos ver como implementar a regressão logística usando o conjunto de dados BreastCancer no pacote mlbench. Você terá que instalar o pacote mlbench para isso.
O objetivo aqui é modelar e prever se um determinado espécime (linha no conjunto de dados) é benigno ou maligno, com base em 9 outras características celulares. Então, vamos carregar os dados e manter apenas os casos completos.
library(mlbench)
data(BreastCancer, package="mlbench")
bc <- BreastCancer[complete.cases(BreastCancer), ] # criar cópia
str(bc)
## 'data.frame': 683 obs. of 11 variables:
## $ Id : chr "1000025" "1002945" "1015425" "1016277" ...
## $ Cl.thickness : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 5 5 3 6 4 8 1 2 2 4 ...
## $ Cell.size : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 1 4 1 8 1 10 1 1 1 2 ...
## $ Cell.shape : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 1 4 1 8 1 10 1 2 1 1 ...
## $ Marg.adhesion : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 1 5 1 1 3 8 1 1 1 1 ...
## $ Epith.c.size : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 2 7 2 3 2 7 2 2 2 2 ...
## $ Bare.nuclei : Factor w/ 10 levels "1","2","3","4",..: 1 10 2 4 1 10 10 1 1 1 ...
## $ Bl.cromatin : Factor w/ 10 levels "1","2","3","4",..: 3 3 3 3 3 9 3 3 1 2 ...
## $ Normal.nucleoli: Factor w/ 10 levels "1","2","3","4",..: 1 2 1 7 1 7 1 1 1 1 ...
## $ Mitoses : Factor w/ 9 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 5 1 ...
## $ Class : Factor w/ 2 levels "benign","malignant": 1 1 1 1 1 2 1 1 1 1 ...
Breast Cancer Dataset: Conjunto de dados de câncer de mama.
O conjunto de dados tem 699 observações e 11 colunas. A coluna Class é a variável de resposta (dependente) e informa se um determinado tecido é maligno (malignant) ou benigno (benign).
Exceto Id, todas as outras colunas são fatores. Isso é um problema quando você modela esse tipo de dados.
Porque, quando você constrói um modelo logístico com variáveis de fator como recursos, ele converte cada nível do fator em uma variável binária fictícia de 1 e 0.
Por exemplo, a forma da célula Cell shape é um fator com 10 níveis. Quando você usa glm para modelar Class como uma função do formato da célula, o formato da célula será dividido em 9 variáveis categóricas binárias diferentes antes de construir o modelo.
Se você pretende construir um modelo logístico sem realizar nenhuma etapa preparatória, o seguinte é o que você pode fazer. Mas não vamos seguir isso, pois há certas coisas a serem feitas antes de construir o modelo logit.
glm(Class ~ Cell.shape, family="binomial", data = bc)
##
## Call: glm(formula = Class ~ Cell.shape, family = "binomial", data = bc)
##
## Coefficients:
## (Intercept) Cell.shape.L Cell.shape.Q Cell.shape.C Cell.shape^4
## 4.189 20.911 6.848 5.763 -1.267
## Cell.shape^5 Cell.shape^6 Cell.shape^7 Cell.shape^8 Cell.shape^9
## -4.439 -5.183 -3.013 -1.289 -0.860
##
## Degrees of Freedom: 682 Total (i.e. Null); 673 Residual
## Null Deviance: 884.4
## Residual Deviance: 243.6 AIC: 263.6
Vamos ver como pode ser o código para construir um modelo logístico. Voltarei a esta etapa mais tarde, pois há algumas etapas de pré-processamento a serem realizadas antes de construir o modelo.
No modelo acima, Class é modelado apenas como uma função de Cell.shape.
Mas observe a saída, o Cell.Shape foi dividido em 9 variáveis diferentes. Isso ocorre porque, como Cell.Shape é armazenado como uma variável de fator, glm cria 1 variável binária (também conhecida como variável fictícia ou dummy) para cada um dos 10 níveis categóricos de Cell.Shape.
Claramente, a partir do significado de Cell.Shape, parece haver algum tipo de ordenação dentro dos níveis categóricos de Cell.Shape. Ou seja, um valor de formato de célula de 2 é maior do que o formato de célula 1 e assim por diante.
Este é o caso com outras variáveis no conjunto de dados. Portanto, é preferível convertê-los em variáveis numéricas e remover a coluna id.
Se fosse uma variável categórica pura sem ordenação interna, como, digamos, o sexo do paciente, você pode deixar essa variável como um fator em si.
# removendo a coluna id
bc <- bc[,-1]
# Pode-se conver fatores em numéricos usando os comandos abaixo
#for(i in 1:9) {
# bc[, i] <- as.numeric(as.character(bc[, i]))
#}
Outro ponto importante a ser observado. Ao converter um fator em uma variável numérica, você deve sempre convertê-lo em caractere e depois em numérico, caso contrário, os valores podem ser confundidos.
Agora todas as colunas são numéricas.
Também gostaria de codificar a variável de resposta em uma variável de fator de 1 e 0. Porém, esta é apenas uma etapa opcional.
Portanto, sempre que a Class for maligna, será 1, senão, será 0. Então, estamos convertendo-a em um fator.
bc$Class <- factor(ifelse(bc$Class == "malignant", 1, 0))
str(bc)
## 'data.frame': 683 obs. of 10 variables:
## $ Cl.thickness : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 5 5 3 6 4 8 1 2 2 4 ...
## $ Cell.size : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 1 4 1 8 1 10 1 1 1 2 ...
## $ Cell.shape : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 1 4 1 8 1 10 1 2 1 1 ...
## $ Marg.adhesion : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 1 5 1 1 3 8 1 1 1 1 ...
## $ Epith.c.size : Ord.factor w/ 10 levels "1"<"2"<"3"<"4"<..: 2 7 2 3 2 7 2 2 2 2 ...
## $ Bare.nuclei : Factor w/ 10 levels "1","2","3","4",..: 1 10 2 4 1 10 10 1 1 1 ...
## $ Bl.cromatin : Factor w/ 10 levels "1","2","3","4",..: 3 3 3 3 3 9 3 3 1 2 ...
## $ Normal.nucleoli: Factor w/ 10 levels "1","2","3","4",..: 1 2 1 7 1 7 1 1 1 1 ...
## $ Mitoses : Factor w/ 9 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 5 1 ...
## $ Class : Factor w/ 2 levels "0","1": 1 1 1 1 1 2 1 1 1 1 ...
A variável de resposta Class é agora uma variável de fator e todas as outras colunas são numéricas.
Tudo bem, as classes de todas as colunas estão definidas. Vamos prosseguir para a próxima etapa.
Antes de construir o regressor logístico, você precisa dividir aleatoriamente os dados em amostras de treinamento e teste.
Como a variável de resposta é uma variável categórica binária, você precisa se certificar de que os dados de treinamento tenham uma proporção aproximadamente igual de classes.
table(bc$Class)
##
## 0 1
## 444 239
As classes ‘benigno’ e ‘maligno’ são divididas aproximadamente na proporção de 1:2.
É evidente que existe um desequilíbrio de classes. Portanto, antes de construir o modelo logístico, precisamos construir as amostras de modo que os 1s e 0s estejam em proporções aproximadamente iguais.
Essa preocupação normalmente é tratada com algumas técnicas chamadas: Down Sampling, Up Sampling, outras.
Então, o que é Down Sampling e Up Sampling?
Na amostragem reduzida, a classe majoritária é amostrada aleatoriamente para ter o mesmo tamanho que a classe menor. Isso significa que, ao criar o conjunto de dados de treinamento, as linhas com a classe benigna serão selecionadas menos vezes durante a amostragem aleatória.
Da mesma forma, em UpSampling, as linhas da classe minoritária, ou seja, maligno, são amostradas repetidamente até atingir o mesmo tamanho da classe majoritária (benigno).
Mas no caso de amostragem híbrida, pontos de dados artificiais são gerados e são sistematicamente adicionados em torno da classe minoritária. Isso pode ser implementado usando os pacotes SMOTE e ROSE.
No entanto, para este exemplo, vou mostrar como fazer up e down sampling.library(caret)
'%ni%' <- Negate('%in%') # define 'not in' func
options(scipen=999) # prevents printing scientific notations.
# Prep Training and Test data.
set.seed(100)
trainDataIndex <- createDataPartition(bc$Class, p=0.7, list = F) # 70% training data
trainData <- bc[trainDataIndex, ]
testData <- bc[-trainDataIndex, ]
No fragmento acima, carreguamos o pacote de caret e usamos a função createDataPartition para gerar os números de linha para o conjunto de dados de treinamento. Ao definir p = 0.7, escolhemos 70% das linhas para ir para trainData e os 30% restantes para ir para testData.
table(trainData$Class)
##
## 0 1
## 311 168
Existem aproximadamente 2 vezes mais amostras benignas. Então, vamos reduzir a resolução usando a função downSample do mesmo pacote.
Para fazer isso, somente precisamos fornecer as variáveis \(X\) e \(Y\) como argumentos.
# Down Sample
set.seed(100)
down_train <- downSample(x = trainData[, colnames(trainData) %ni% "Class"],
y = trainData$Class)
table(down_train$Class)
##
## 0 1
## 168 168
Benigno e maligno estão agora na mesma proporção. O %ni% é a negação da função %in% e usamos aqui para selecionar todas as colunas, exceto a coluna Class.
A função downSample requer o \(y\) como uma variável de fator, é por isso que convertimos Class em um fator nos dados originais.
Excelente! Agora, vamos fazer o upsampling usando a função upSample. Ele segue uma sintaxe semelhante a downSample.
# Up Sample.
set.seed(100)
up_train <- upSample(x = trainData[, colnames(trainData) %ni% "Class"],
y = trainData$Class)
table(up_train$Class)
##
## 0 1
## 311 311
Como esperado, benigno e maligno estão agora na mesma proporção. Usaremos a versão downSampled do conjunto de dados para construir o modelo logit na próxima etapa.
# Construindo o Modelo Logístico
logitmod <- glm(Class ~ Cl.thickness + Cell.size + Cell.shape, family = "binomial", data=down_train)
summary(logitmod)
##
## Call:
## glm(formula = Class ~ Cl.thickness + Cell.size + Cell.shape,
## family = "binomial", data = down_train)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.9152 -0.1237 0.0000 0.0000 3.2886
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 21.9416 4019.6894 0.005 0.996
## Cl.thickness.L 33.9297 6086.7759 0.006 0.996
## Cl.thickness.Q 0.8495 4875.2763 0.000 1.000
## Cl.thickness.C 6.6813 4869.4753 0.001 0.999
## Cl.thickness^4 -7.0850 5601.3115 -0.001 0.999
## Cl.thickness^5 6.0546 7466.2277 0.001 0.999
## Cl.thickness^6 -4.3105 5819.2129 -0.001 0.999
## Cl.thickness^7 -9.0908 4183.1343 -0.002 0.998
## Cl.thickness^8 -13.4431 5429.1383 -0.002 0.998
## Cl.thickness^9 -9.0762 3805.9905 -0.002 0.998
## Cell.size.L 14.6655 9175.0629 0.002 0.999
## Cell.size.Q -21.3572 7128.8892 -0.003 0.998
## Cell.size.C 3.0727 5359.4600 0.001 1.000
## Cell.size^4 32.7610 10439.1068 0.003 0.997
## Cell.size^5 13.2450 12068.3401 0.001 0.999
## Cell.size^6 -15.6867 10704.1590 -0.001 0.999
## Cell.size^7 -2.8623 7425.0501 0.000 1.000
## Cell.size^8 10.4855 6654.1545 0.002 0.999
## Cell.size^9 -7.4842 7029.6069 -0.001 0.999
## Cell.shape.L 21.4241 10161.6068 0.002 0.998
## Cell.shape.Q 4.3431 5535.7296 0.001 0.999
## Cell.shape.C -0.6347 5174.5686 0.000 1.000
## Cell.shape^4 -0.1887 10567.5848 0.000 1.000
## Cell.shape^5 0.6414 13106.3506 0.000 1.000
## Cell.shape^6 -2.5925 11011.4241 0.000 1.000
## Cell.shape^7 -9.8297 7232.4121 -0.001 0.999
## Cell.shape^8 -11.8856 5255.5429 -0.002 0.998
## Cell.shape^9 -8.1257 3141.5453 -0.003 0.998
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 465.795 on 335 degrees of freedom
## Residual deviance: 61.453 on 308 degrees of freedom
## AIC: 117.45
##
## Number of Fisher Scoring iterations: 21
O logitmod agora está construído. Agora pode ser usado para prever a resposta em testData.
pred <- predict(logitmod, newdata = testData, type = "response")
Agora, pred contém a probabilidade de que a observação seja maligna para cada observação.
A prática comum é considerar o corte de probabilidade como 0.5. Se a probabilidade de \(Y\) for maior do que 0.5, ele pode ser classificado como um evento (maligno).
Portanto, se pred for maior que 0.5, é maligno, do contrário, é benigno.
y_pred_num <- ifelse(pred > 0.5, 1, 0)
y_pred <- factor(y_pred_num, levels=c(0, 1))
y_act <- testData$Class
Vamos calcular a precisão, que nada mais é do que a proporção de y_pred que corresponde a y_act.
mean(y_pred == y_act)
## [1] 0.9460784
Significa que temos uma taxa de acerto de 94%.
Tudo bem, porque precisamos cuidar do desequilíbrio de classe. Para entender isso, vamos supor que você tenha um conjunto de dados onde 95% dos valores de \(Y\) pertencem à classe benigna e 5% pertencem à classe maligna.
Se tivessemos previsto cegamente todos os pontos de dados como benignos, alcançariamos uma porcentagem de precisão de 95%. O que parece muito alto. Mas obviamente isso é errado. O que importa é o quão bem você prevê as classes malignas e isso requer que as classes benignas e malignas estejam equilibradas.