Exemplo do livro Fávero e Belfiore (2017), Tabela 9.11. O exercício se baseará no método do encadeamento único.
Obs.: Para visualização das referências, ver ementa da disciplina.
# Abrindo a base de dados
library(haven)
exemp3 <- read_dta("D:/dados/exemp3.dta")
View(exemp3)
attach(exemp3)
# Deixando os nomes como label
row.names(exemp3) <- exemp3$estudante
## Warning: Setting row names on a tibble is deprecated.
exemp3[1] <- NULL
# Pacotes
library(dendextend)
##
## ---------------------
## Welcome to dendextend version 1.10.0
## Type citation('dendextend') for how to cite the package.
##
## Type browseVignettes(package = 'dendextend') for the package vignette.
## The github page is: https://github.com/talgalili/dendextend/
##
## Suggestions and bug-reports can be submitted at: https://github.com/talgalili/dendextend/issues
## Or contact: <tal.galili@gmail.com>
##
## To suppress this message use: suppressPackageStartupMessages(library(dendextend))
## ---------------------
##
## Attaching package: 'dendextend'
## The following object is masked from 'package:stats':
##
## cutree
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(ggplot2)
# Base de dados padronizada
exemp3_scaled<-data.frame(scale(exemp3))
# Distância
dist_m<-dist(exemp3, method = "euclidean")
# Método de aglomeração hierárquica- Único
hc_m<-hclust(dist_m,method="single")
# Dendograma
plot(hc_m, hang=-1, cex=0.5, main=" Método Único")
#Ou
dend_d <- as.dendrogram(hc_m)
plot(dend_d, hang=-1)
## Warning in plot.window(...): "hang" não é um parâmetro gráfico
## Warning in plot.xy(xy, type, ...): "hang" não é um parâmetro gráfico
## Warning in axis(side = side, at = at, labels = labels, ...): "hang" não é
## um parâmetro gráfico
## Warning in axis(side = side, at = at, labels = labels, ...): "hang" não é
## um parâmetro gráfico
## Warning in title(...): "hang" não é um parâmetro gráfico
dend_20 <- color_branches(dend_d, h = 3)
plot(dend_20)
plot(dend_20)
# Número de Cluster
library(NbClust)
devAskNewPage(ask=TRUE)
nc<-NbClust(exemp3_scaled, distance="euclidean", min.nc=2, max.nc=4, method="single", index="ch")
nc
## $All.index
## 2 3 4
## 2.8970 4.0886 5.2366
##
## $Best.nc
## Number_clusters Value_Index
## 4.0000 5.2366
##
## $Best.partition
## Gabriela Luiz Felipe Patrícia Ovídio Leonor
## 1 2 3 4 1
# Modelo K-means- iter.max=0 no default
model_exemp3 <- kmeans(exemp3, centers = 3)
print(model_exemp3)
## K-means clustering with 3 clusters of sizes 1, 1, 3
##
## Cluster means:
## matematica fisica quimica
## 1 7.8 8.0 1.5
## 2 8.9 1.0 2.7
## 3 4.7 1.9 7.7
##
## Clustering vector:
## Gabriela Luiz Felipe Patrícia Ovídio Leonor
## 3 1 2 3 3
##
## Within cluster sum of squares by cluster:
## [1] 0.00 0.00 20.38
## (between_SS / total_SS = 81.1 %)
##
## Available components:
##
## [1] "cluster" "centers" "totss" "withinss"
## [5] "tot.withinss" "betweenss" "size" "iter"
## [9] "ifault"
# Extraindo o vetor dos cluster com base no modelo kmeans
clust_exemp3 <- model_exemp3$cluster
# Nova base de dados com os dados originais mais o cluster
lineup_exemp3 <- mutate(exemp3, cluster = clust_exemp3)
# Gráficos de dispersão
graph1<-ggplot(lineup_exemp3, aes(matematica,
fisica, colour=factor(cluster))) + geom_point()
graph1
graph2<-ggplot(lineup_exemp3, aes(matematica,
quimica, colour=factor(cluster))) + geom_point()
graph2
graph3<-ggplot(lineup_exemp3, aes(fisica,
quimica, colour=factor(cluster))) + geom_point()
graph3
# Médias por cluster
library(stats)
# Média de matemática por cluster
tapply(lineup_exemp3$matematica, lineup_exemp3$cluster, summary)
## $`1`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 7.8 7.8 7.8 7.8 7.8 7.8
##
## $`2`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 8.9 8.9 8.9 8.9 8.9 8.9
##
## $`3`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 3.40 3.55 3.70 4.70 5.35 7.00
# Média de física por cluster
tapply(lineup_exemp3$fisica, lineup_exemp3$cluster, summary)
## $`1`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 8 8 8 8 8 8
##
## $`2`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1 1 1 1 1 1
##
## $`3`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.00 1.50 2.00 1.90 2.35 2.70
# Média de química por cluster
tapply(lineup_exemp3$quimica, lineup_exemp3$cluster, summary)
## $`1`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.5 1.5 1.5 1.5 1.5 1.5
##
## $`2`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 2.7 2.7 2.7 2.7 2.7 2.7
##
## $`3`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 5.00 7.00 9.00 7.70 9.05 9.10
# Anova
anova1 = anova(lm(formula = matematica ~ factor(cluster), data=lineup_exemp3))
anova1
## Analysis of Variance Table
##
## Response: matematica
## Df Sum Sq Mean Sq F value Pr(>F)
## factor(cluster) 2 16.592 8.296 2.0792 0.3248
## Residuals 2 7.980 3.990
anova2 = anova(lm(formula = fisica ~ factor(cluster), data=lineup_exemp3))
anova2
## Analysis of Variance Table
##
## Response: fisica
## Df Sum Sq Mean Sq F value Pr(>F)
## factor(cluster) 2 32.612 16.306 22.337 0.04285 *
## Residuals 2 1.460 0.730
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova3 = anova(lm(formula = quimica ~ factor(cluster), data=lineup_exemp3))
anova3
## Analysis of Variance Table
##
## Response: quimica
## Df Sum Sq Mean Sq F value Pr(>F)
## factor(cluster) 2 38.352 19.176 3.5057 0.2219
## Residuals 2 10.940 5.470
### Número relevante de clusters no conjunto de dados
#nc<-NbClust(exemp3, min.nc=2, max.nc=4, method="kmeans")
#nc