Análise exploratória do conjunto de dados iris com R
Publicado em: Sep 30, 2020 | Última modificação: Oct 14, 2020

O conjunto de dados Iris é um dos mais explorados e conhecidos na análise de dados. Foi introduzido por Ronald Fisher em seu artigo "The use of multiple measurements in taxonomic problems" de 1936. O artigo contém três espécies de flores (setosa, virginica e versicolor), com quatro características medidas em cada uma: o comprimento e largura da pétala e da sépala. Essas medidas quantificam a diversidade morfológica dentro das espécies. Cada uma das medidas é dada em centímetros.

Importando os dados

Vamos importar as bibliotecas necessárias:

library(tidyverse)
library(gridExtra)

O conjunto de dados está disponível no R através do dataframe iris. Vamos definir uma variável d para armazenar o conjunto de dados como um tibble do R. Tibbles são uma "versão moderna" dos dataframes. A ideia é que eles mantenham boa parte das características do dataframe com poucas modificações (por exemplo, eles não convertem strings para factors como os dataframes).

d <- iris %>% as_tibble
d %>% glimpse
Rows: 150
Columns: 5
$ Sepal.Length <dbl> 5.1, 4.9, 4.7, 4.6, 5.0, 5.4, 4.6, 5.0, 4.4, 4.9, 5.4, 4…
$ Sepal.Width  <dbl> 3.5, 3.0, 3.2, 3.1, 3.6, 3.9, 3.4, 3.4, 2.9, 3.1, 3.7, 3…
$ Petal.Length <dbl> 1.4, 1.4, 1.3, 1.5, 1.4, 1.7, 1.4, 1.5, 1.4, 1.5, 1.5, 1…
$ Petal.Width  <dbl> 0.2, 0.2, 0.2, 0.2, 0.2, 0.4, 0.3, 0.2, 0.2, 0.1, 0.2, 0…
$ Species      <fct> setosa, setosa, setosa, setosa, setosa, setosa, setosa, …

Daí podemos ver que o conjunto de dados possui 150 observações, analisadas através de cinco atributos. Desses cinco atributos, quatro são numéricos (Sepal.Length, Sepal.Width, Petal.Length e Petal.Width) e uma é categórica (Species). Vamos checar também a falta de algum dos atributos:

any(is.na(d))
[1] FALSE

Assim, este conjunto de dados não carece de nenhum. Normalmente esse não é o caso!

Sumário estatístico

O sumário estatístico nos possibilita ter um panorama rápido dos dados. Podemos checar os outliers olhando para os valores máximo e mínimo de cada coluna em relação à média.

d %>% summary
  Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
 Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
 1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
 Median :5.800   Median :3.000   Median :4.350   Median :1.300  
 Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
 3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
 Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500  
       Species  
 setosa    :50  
 versicolor:50  
 virginica :50

Daqui podemos ver também que as observações foram feitas em 50 indivíduos (flores) de cada uma das três espécies (setosa, versicolor e virginica).

A informação de amplitude das variáveis numéricas pode ser visualizada através de um gráfico de densidade e histograma:

theme_set(
    theme_gray(base_size = 22)
)

plot.num.var.barplot <- function (attr) {
  col <- enquo(attr)
  text.height <- 0.6
  d %>% ggplot(aes(x = !!col)) +
    geom_histogram(
      aes(y = ..density..),
      color = "darkgray",
      fill  = "darkblue",
      alpha = 0.6,
      binwidth = 0.3,
      position = "identity"
    ) +
    geom_density(
      alpha = 0.7,
      size = 1,
      color = "dodgerblue4"
    ) +
    geom_vline(
      aes(xintercept = mean(!!col)),
      color = "grey28",
      linetype = "dashed",
      size = 1
    ) + 
    geom_label(
      aes(
        x = mean(!!col),
        y = text.height - 0.1,
        label = "Média", 
      ), 
      colour = "black", 
      hjust = 0
    ) +
    geom_vline(
      aes(xintercept = median(!!col)),
      color = "green",
      linetype = "dashed",
      size = 1
    ) +
    geom_label(
      aes(
        x = median(!!col), 
        y = text.height,
        label = "Mediana", 
      ), 
      colour = "green", 
      hjust = 0,
      vjust = 0
    )
}

g1 <- plot.num.var.barplot(Sepal.Length)
g2 <- plot.num.var.barplot(Sepal.Width)
g3 <- plot.num.var.barplot(Petal.Length)
g4 <- plot.num.var.barplot(Petal.Width)

grid.arrange(g1, g2, g3, g4, nrow = 2)

20200930-111134-sumario-estatistico-1.png

Do gráfico podemos ver por exemplo que os atributos Petal.Length e Petal.Width são bimodais, isto é, apresentam dois grupos ou populações características. Nestes dois casos, note como a média e a mediana não são próximas.

Esta análise preliminar não leva em consideração a distinção entre cada espécie. Olhando para cada tipo de flor separadamente podemos ganhar maior compreensão dos dados:

plot.num.boxplot <- function(attr) { 
col <- enquo(attr)

d %>% ggplot(aes(x = Species, y = !!col, fill = Species)) + 
  geom_boxplot(
    alpha = 0.3,
    outlier.colour = "red",
    outlier.fill = "red",
    outlier.size = 2
  ) + 
  theme(legend.position = "none")
}

g1 <- plot.num.boxplot(Sepal.Length)
g2 <- plot.num.boxplot(Sepal.Width)
g3 <- plot.num.boxplot(Petal.Length)
g4 <- plot.num.boxplot(Petal.Width)

grid.arrange(g1, g2, g3, g4, nrow = 2)

20200930-111134-sumario-estatistico-2.png

Daí podemos observar que a espécie setosa é a que possui um número maior de outliers (pontos "fora da curva").

Análise de correlação

Uma vez que temos uma maior "intimidade" com o conjunto de dados, podemos prosseguir e explorar as relações entre os diferentes atributos. Para isso vamos calcular a matriz de correlação para os atributos do nosso conjunto:

cor(d[,1:4])
             Sepal.Length Sepal.Width Petal.Length Petal.Width
Sepal.Length    1.0000000  -0.1175698    0.8717538   0.8179411
Sepal.Width    -0.1175698   1.0000000   -0.4284401  -0.3661259
Petal.Length    0.8717538  -0.4284401    1.0000000   0.9628654
Petal.Width     0.8179411  -0.3661259    0.9628654   1.0000000

Daí podemos ver que o comprimento da sépala é fortemente correlacionado com o comprimento e a largura da pétala: isso indica que as flores com sépalas mais compridas tendem a ter pétalas maiores e mais largas.

Também podemos observar uma correlação negativa entre a largura da sépala com o comprimento e largura das pétalas, indicando que sépalas mais largas tendem a ter pétalas menores e mais afiladas.

Podemos reforçar essas hipóteses olhando para um gráfico de dispersão (scatterplot):

theme_set(
    theme_gray(base_size = 22)
)

makePairs <- function(data) 
{
  grid <- expand.grid(x = 1:ncol(data), y = 1:ncol(data))
  grid <- subset(grid, x != y)
  all <- do.call("rbind", lapply(1:nrow(grid), function(i) {
    xcol <- grid[i, "x"]
    ycol <- grid[i, "y"]
    data.frame(xvar = names(data)[ycol], yvar = names(data)[xcol], 
               x = data[, xcol], y = data[, ycol], data)
  }))
  all$xvar <- factor(all$xvar, levels = names(data))
  all$yvar <- factor(all$yvar, levels = names(data))
  densities <- do.call("rbind", lapply(1:ncol(data), function(i) {
    data.frame(xvar = names(data)[i], yvar = names(data)[i], x = data[, i])
  }))
  list(all=all, densities=densities)
}
 
# expand iris data frame for pairs plot
gg1 = makePairs(iris[,-5])
 
mega_iris = data.frame(gg1$all, Species=rep(iris$Species, length=nrow(gg1$all)))
 
# pairs plot
ggplot(mega_iris, aes_string(x = "x", y = "y")) + 
  facet_grid(xvar ~ yvar, scales = "free") + 
  geom_point(aes(colour=Species), na.rm = TRUE, alpha=0.8) + 
  stat_density(aes(x = x, y = ..scaled.. * diff(range(x)) + min(x)), 
               data = gg1$densities, position = "identity", 
               colour = "grey20", geom = "line") + geom_smooth(method="lm")

20200930-111134-analise-correlacao.png

Modelagem dos dados

Para este conjunto de dados, sabemos de antemão as espécies de flores para cada conjunto de observações. Mesmo assim é interessante pensar como poderíamos prever a espécie de flor dado um conjunto de medidas para a pétala e a sépala. Para isso vamos utilizar o algoritmo da árvore de decisão.

Como o propósito aqui é meramente ilustrativo, vamos desconsiderar por ora os problemas relacionados à amostragem do conjunto de dados, vamos tomar uma pequena amostra aleatória como conjunto de treino:

set.seed(1234)
iris_train <- d %>% slice_sample(n = 50)
iris_test <- anti_join(d, iris_train)
iris_train %>% glimpse
Rows: 50
Columns: 5
$ Sepal.Length <dbl> 5.2, 5.7, 6.3, 6.5, 6.3, 6.4, 6.8, 7.9, 6.2, 7.1, 5.5, 5…
$ Sepal.Width  <dbl> 3.5, 2.6, 3.3, 3.2, 3.4, 2.8, 3.2, 3.8, 2.9, 3.0, 2.5, 2…
$ Petal.Length <dbl> 1.5, 3.5, 6.0, 5.1, 5.6, 5.6, 5.9, 6.4, 4.3, 5.9, 4.0, 3…
$ Petal.Width  <dbl> 0.2, 1.0, 2.5, 2.0, 2.4, 2.2, 2.3, 2.0, 1.3, 2.1, 1.3, 1…
$ Species      <fct> setosa, versicolor, virginica, virginica, virginica, vir…

Construindo o modelo

Vamos construir a árvore de decisão usando o algoritmo rpart (Recursive Partitioning And Regression Trees) algorithm:

library(rpart)

# building the classification tree with rpart
tree <- rpart(
  Species ~ Sepal.Length + Sepal.Width + Petal.Length + Petal.Width,
  data = iris_train,
  method = "class"
)
# Visualize the decision tree with rpart.plot
library(rpart.plot)

rpart.plot(tree, nn=TRUE, box.palette="Blues")

20200930-111134-arvore-decisao.png

Testando o modelo e calculando a eficiência

pred <- predict(object=tree, iris_test %>% select(-Species), type="class")
t <- table(iris_test$Species, pred) 
caret::confusionMatrix(t)
Confusion Matrix and Statistics

            pred
             setosa versicolor virginica
  setosa         35          0         0
  versicolor      0         33         2
  virginica       0          2        28

Overall Statistics
                                         
               Accuracy : 0.96           
                 95% CI : (0.9007, 0.989)
    No Information Rate : 0.35           
    P-Value
NIR] : < 2.2e-16      
                                         
                  Kappa : 0.9398         
                                         
 Mcnemar's Test P-Value : NA             

Statistics by Class:

                     Class: setosa Class: versicolor Class: virginica
Sensitivity                   1.00            0.9429           0.9333
Specificity                   1.00            0.9692           0.9714
Pos Pred Value                1.00            0.9429           0.9333
Neg Pred Value                1.00            0.9692           0.9714
Prevalence                    0.35            0.3500           0.3000
Detection Rate                0.35            0.3300           0.2800
Detection Prevalence          0.35            0.3500           0.3000
Balanced Accuracy             1.00            0.9560           0.9524

A saída mostra que as amostras foram classificadas corretamente com com um intervalo de confiança de 95%. Assim, podemos classificar as espécies de flores com razoável confiança usando um modelo de árvore de decisão.

Victor Santos

Físico, Universidade Federal do Ceará. Buracos negros, relatividade geral, computação e tudo mais.