Utilizando o dataset Iris

Os passos contidos nesse tutorial foram retirados da seguinte página:

[https://rstudio-pubs-static.s3.amazonaws.com/574220_2ddd7eb0d872432886602576d6f9801a.html]

Os dados consistem de 50 unidades amostrais de três espécies (setosa, virginica, versicolor) de íris (uma espécie de planta). Temos então um total de 150 unidades amostrais.

head(iris)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa
summary(iris)
##   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  
##                 
##                 
## 

Se forem observar a função head informa os primeiros valores da tabela de dados. Percebe-se que temos 4 variáveis que foram analisadas:

Vamos explorar e conhecer um pouco mais sobre esses dados?

table(iris$Species)
## 
##     setosa versicolor  virginica 
##         50         50         50

Podemos verificar que temos 50 amostras por espécie da planta. Isso foi apresentado no início desse estudo, mas agora com a função “table” podemos comprovar na prática.

Vamos investigar um pouco mais?? Inicialmente com a fução “summary” obtivemos a estatística descritiva de cada variável independentemente de cada espécie. Agora vamos obter os dados de maneira mais detalhada.

tapply(X = iris$Sepal.Length, INDEX = list(iris$Species), FUN = mean)
##     setosa versicolor  virginica 
##      5.006      5.936      6.588

Utilizamos acima a função “tapply”. Vamos agora proceder o mesmo cálculo, porém com outra função. Veja como fica:

aggregate(Sepal.Length ~ Species, data = iris, mean)
##      Species Sepal.Length
## 1     setosa        5.006
## 2 versicolor        5.936
## 3  virginica        6.588

Vamos agora utilizar esse mesmo comando para o cálculo da média das outras variáveis:

aggregate(Sepal.Length ~ Species, data = iris, mean)
##      Species Sepal.Length
## 1     setosa        5.006
## 2 versicolor        5.936
## 3  virginica        6.588
aggregate(Sepal.Width ~ Species, data = iris, mean)
##      Species Sepal.Width
## 1     setosa       3.428
## 2 versicolor       2.770
## 3  virginica       2.974
aggregate(Petal.Length ~ Species, data = iris, mean)
##      Species Petal.Length
## 1     setosa        1.462
## 2 versicolor        4.260
## 3  virginica        5.552

Vamos agora utilizar dessa mesma função “aggregate” para calcular o desvio padrão:

aggregate(Sepal.Length ~ Species, data = iris, sd)
##      Species Sepal.Length
## 1     setosa    0.3524897
## 2 versicolor    0.5161711
## 3  virginica    0.6358796
aggregate(Sepal.Width ~ Species, data = iris, sd)
##      Species Sepal.Width
## 1     setosa   0.3790644
## 2 versicolor   0.3137983
## 3  virginica   0.3224966
aggregate(Petal.Length ~ Species, data = iris, sd)
##      Species Petal.Length
## 1     setosa    0.1736640
## 2 versicolor    0.4699110
## 3  virginica    0.5518947

Interessante né? Mas vamos tentar otimizar e colocar esses resultados de médias e desvio padrão numa tabela só:

medias <- matrix(NA, ncol = 3, nrow = 4)
colnames(medias) <- unique(iris$Species)
rownames(medias) <- names(iris)[-5]
for (i in 1:4){medias[i,] <- tapply(iris[,i], iris$Species, mean)}
medias
##              setosa versicolor virginica
## Sepal.Length  5.006      5.936     6.588
## Sepal.Width   3.428      2.770     2.974
## Petal.Length  1.462      4.260     5.552
## Petal.Width   0.246      1.326     2.026

Agora que aprendemos alguns conceitos básicos de manipulação no R, vamos falar da Estatística descritiva.

Noções sobre estatística descritiva

Vamos falar um pouco sobre:

Parâmetro Descrição Função no R
Média Média aritmética mean()
Mediana valor central median()
Moda Valor mais frequente sort(table(), decreasing = TRUE)[1]
Desvio padrão variação em torno da média sd()
Quantis pontos de corte dividindo uma distribuição de probabilidade quantile()

Média

A média é um dos fundamentos para obtenção de medidas de tendência central.

\(\overline{x}=\frac{\sum_{i=1}^{n}\times _{i}}{n}\)

vars <- iris[, -5]
apply(vars, 2, mean)
## Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
##     5.843333     3.057333     3.758000     1.199333

Mediana

Vamos ao cálculo da mediana. Como dito anteriormente, é o valor central.

apply(vars, 2, median)
## Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
##         5.80         3.00         4.35         1.30

Moda

Calculando agora a moda da minha amostra:

moda <- sort(table(iris$Petal.Length), decreasing = TRUE)
moda[1]
## 1.4 
##  13

Medidas de dispersão

Agora que vimos alguns cálculos básicos de medidas de tendência central é necessário verificar quanto desses valores estão espalhados em torno da média. Isso chama-se na linguagem estatística de medidas de dispersão.

Variância da amostra

A variância mede o quanto os dados estão dispersos em torno da média, ou seja, do valor esperado.

Monitorar a variância é essencial para as indústrias de manufatura e qualidade porque a redução da variância dos processos aumenta a precisão e diminui o número de defeitos dos produtos manufaturados.

\(s^{2}=\sum_{i=1}^{n}\frac{(x_{i}-x)^{2}}{n-1}\)

apply(vars, 2, var)
## Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
##    0.6856935    0.1899794    3.1162779    0.5810063

Desvio Padrão

O desvio padrão é uma medida que indica a dispersão dos dados dentro de uma amostra com relação à média.

\(s=\sqrt{\sum_{i=1}^{n}\frac{(x_{i}-x)^{2}}{n-1}}\)

sd01 <- apply(vars, 2, sd)
sd01
## Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
##    0.8280661    0.4358663    1.7652982    0.7622377

Correlação

Podemos avaliar a correlação entre as variáveis.

Tipos de Correlação

cor(vars)
##              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

Métodos Gráficos

Vamos agora falar dos métodos gráficos:

Tipo de gráfico Descrição
barras cada barra representa a frequência de observações de um grupo
histograma cada barra representa a frequência de observalçoes em um dado conjunto de valores
densidade estimativa da distribuição estatística baseada nos dados
quantil-quantil comparação dos dados em relação a uma distribuição normal
box-plot representação visual da média, quartis, simetria e outliers
dispersão display gráfico de variáveis contínuas nos eixos x e y

O R em seu pacote base possui funções para execução dos gráficos acima. Contudo, pode ser instalado o pacote “ggplot2” que é uma ferramenta poderosa para ambiente gráfico.

Na tabela abaixo segue um compilado das funções do Rbase e do pacote ggplot2 para obtenção dos gráficos apresentados anteriormente.

Tipo de gráfico Rbase ggplot2
barras barplot() geom_bar()
histograma histogram() geom_histogram()
densidade plot(density()) geom_density()
quantil-quantil qqnorm() geom_qq()
box-plot boxplot() geom_boxplot()
dispersão plot() geom_point()

Gráfico de barras

barplot(table(iris$Species))

barplot(table(iris$Sepal.Width))

Histograma

par(mfrow=c(2, 2))
hist(iris$Sepal.Length)
hist(iris$Sepal.Width)
hist(iris$Petal.Length)
hist(iris$Petal.Length)

par(mfrow=c(1, 2))
hist(iris$Sepal.Width)
hist(iris$Sepal.Width, breaks = 4)

Densidade

A curva de densidade mostra a probabilidade de observar determinado valor. Em comparação ao histograma, no eixo y, ao invés de termos a frequência, temos a densidade probabilística.

par(mfrow=c(1, 2))
hist(iris$Sepal.Width)
hist(iris$Sepal.Width, freq = FALSE)

Vamos utilizar a função density.

par(mfrow=c(1, 2))
# plot da curva de densidade
plot(density(iris$Sepal.Width))
# plot da curva de densidade sobre o histograma de densidade
hist(iris$Sepal.Width, freq = FALSE)
lines(density(iris$Sepal.Width), col="blue")

Box-plot

O boxplot nos fornece uma análise visual da posição, dispersão, simetria, caudas e valores discrepantes (outliers) do conjunto de dados.

  • Posição: Em relação à posição dos dados, observa-se a linha central do retângulo (a mediana ou segundo quartil).

  • Dispersão: A dispersão dos dados pode ser representada pelo intervalo interquatilico que é a diferença entre o terceiro quartil e o primeiro quartil (tamanho da caixa), ou ainda pela amplitude que é calculada da seguinte maneira: valor máximo – valor mínimo. Embora a amplitude seja de fácil entendimento, o intervalo interqualitico é uma estatística mais robusta para medir variabilidade uma vez que não sofre influência de outliers.

  • Simetria: Um conjunto de dados que tem uma distribuição simétrica, terá a linha da mediana no centro do retângulo. Quando a linha da mediana está próxima ao primeiro quartil, os dados são assimétricos positivos e quando a posição da linha da mediana é próxima ao terceiro quartil, os dados são assimétricos negativos. Vale ressaltar que a mediana é a medida de tendência central mais indicada quando os dados possuem distribuição assimétrica, uma vez que a média aritmética é influenciada pelos valores extremos.

  • Caudas: As linhas que vão do retângulo até aos outliers podem fornecer o comprimento das caudas da distribuição.

  • Outliers: Já os outliers indicam possíveis valores discrepantes. No boxplot, as observações são consideradas outliers quando estão abaixo ou acima do limite de detecção de outliers. O limite de detecção de outliers é construído utilizando o intervalo interquartílico, dado pela distância entre o primeiro e o terceiro quartil. Sendo assim, os limites inferior e superior de detecção de outlier são dados por:

  • Limite Inferior = Primeiro Quartil – 1,5 * (Terceiro Quartil – Primeiro Quartil)

  • Limite Superior = Terceiro Quartil + 1,5 * (Terceiro Quartil – Primeiro Quartil)

[Oliveira, Bruno. “O que é Desvio Padrão e Erro Padrão? - Blog da Oper”. Oper, 14 de agosto de 2019, https://operdata.com.br/blog/desvio-padrao-e-erro-padrao/.]

Vamos elaborar os box-plots de nossas variáveis:

boxplot(iris$Sepal.Length)

boxplot(iris$Sepal.Width)

boxplot(iris$Petal.Length)

boxplot(iris$Petal.Width)

Vamos rodar por espécie:

boxplot(Sepal.Length ~ Species, data = iris)

boxplot(Sepal.Width ~ Species, data = iris)

boxplot(Petal.Length ~ Species, data = iris)

boxplot(Petal.Width ~ Species, data = iris)

Podemos perceber a existência de alguns outliers. Vamos usar a própria função boxplot para indentificá-los.

boxplot(iris$Sepal.Width)

my_boxplot <- boxplot(iris$Sepal.Width, plot=FALSE)
my_boxplot
## $stats
##      [,1]
## [1,]  2.2
## [2,]  2.8
## [3,]  3.0
## [4,]  3.3
## [5,]  4.0
## 
## $n
## [1] 150
## 
## $conf
##          [,1]
## [1,] 2.935497
## [2,] 3.064503
## 
## $out
## [1] 4.4 4.1 4.2 2.0
## 
## $group
## [1] 1 1 1 1
## 
## $names
## [1] "1"
# o objeto é uma lista e os valores outliers estão guardados no elemento $out da lista
outliers <- my_boxplot$out
#qual a posicao dos outliers
which(iris$Sepal.Width %in% outliers)
## [1] 16 33 34 61
# vamos usar a posicao para indexar o objeto
iris[which(iris$Sepal.Width %in% outliers), c("Sepal.Width", "Species")]
##    Sepal.Width    Species
## 16         4.4     setosa
## 33         4.1     setosa
## 34         4.2     setosa
## 61         2.0 versicolor

No caso anterior consideramos outliers em relação à distribuição da variável para todas as espécies juntas. É razoável assumir que cada espécie tenha um padrão morfométrico distinto de modo que poderíamos identificar outliers de maneira específica por espécie.

boxplot(Sepal.Width ~ Species, data = iris)

my_boxplot2 <- boxplot(Sepal.Width ~ Species, data=iris, plot=FALSE)
my_boxplot2
## $stats
##      [,1] [,2] [,3]
## [1,]  2.9  2.0  2.2
## [2,]  3.2  2.5  2.8
## [3,]  3.4  2.8  3.0
## [4,]  3.7  3.0  3.2
## [5,]  4.4  3.4  3.8
## 
## $n
## [1] 50 50 50
## 
## $conf
##          [,1]     [,2]     [,3]
## [1,] 3.288277 2.688277 2.910622
## [2,] 3.511723 2.911723 3.089378
## 
## $out
## [1] 2.3
## 
## $group
## [1] 1
## 
## $names
## [1] "setosa"     "versicolor" "virginica"
# o objeto é uma lista e os valores outliers estão guardados no elemento $out da lista
outliers2 <- my_boxplot2$out
# neste caso, queremos apenas os outliers da especie setosa
# vamos usar a posicao para indexar o objeto
iris[iris$Sepal.Width %in% outliers2 & 
       iris$Species == "setosa", 
     c("Sepal.Width", "Species")]
##    Sepal.Width Species
## 42         2.3  setosa

Distribuição dos dados

A análise da distribuição dos dados é fundamental. Os testes paramétricos por exemplo só podem ser aplicados em dados que apresentem uma distribuição normal. Observe a tabela abaixo com algumas distribuições:

Distribuição Tipo E (X) \(\sigma^{2}(X)\) Uso Exemplo
normal contínua \(\sigma^{2}\) A distribuição normal é simétrica em torno da média o que implica que e média, a mediana e a moda são todas coincidentes. distribuição de tamanho
binomial discreta \(np\) \(np(1-p)\) número de sucessos em \(n\) tentativas presença ou ausência de espécies
Poisson discreta \(\lambda\) \(\lambda\) modelar o número de ocorrências de um evento por um certo período de tempo ou por um certo volume ou por uma certa área. número de chegas em um estabelecimento bancário.
Log-normal contínua \(log(\mu)\) \(log(\sigma^{2})\) curva assimétrica Distribuição de abundância de espécies

Vamos verificar a normalidade através dos seguintes gráficos:

qqnorm e qqline

par(mfrow = c(1,3))
qqnorm(iris$Sepal.Length[iris$Species == "setosa"], 
       main = "setosa")
qqline(iris$Sepal.Length[iris$Species == "setosa"])
qqnorm(iris$Sepal.Length[iris$Species == "versicolor"], 
       main = "versicolor")
qqline(iris$Sepal.Length[iris$Species == "versicolor"])
qqnorm(iris$Sepal.Length[iris$Species == "virginica"], 
       main = "virginica")
qqline(iris$Sepal.Length[iris$Species == "virginica"])

Relação entre as variáveis

Esse gráfico mostra as correelações anteriormente calculadas. Podemos observar uma forte correlação entre as variáveis Petalwidth e Petallenght.

pairs(vars)

library("GGally")
## Loading required package: ggplot2
## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2
ggpairs(vars)