Verificando existência de NaN’s

sum(is.na(data)) 

## [1] 0

sum(is.na(data_c)) 

## [1] 0

knitr::kable(summary(data))

knitr::kable(summary(data_c))

Podemos notar que a base não possui NaN’s. Além disso, apartir da ultima tabela gerada acima, podemos notar que não existe nenhum valor, em sua respectiva variavel, que aparenta estar fora do seu range de valores. Esse ultimo ponto de verificar se existe alguma variavel fora do range da sua natureza, foi utilizado com o intuito para verificar se existiam NaN’s tratados como 999 ou -1 por exemplo, como sabemos que muitas situações é comum tratarem observações com numeros bem discripantes da natureza da sua variavel para representar dados faltantes.

Definindo variável target

Vamos definir a nossa variável target para o nosso problema de classificação

target <- as.factor(data$Target)
levels(target)

## [1] "Dropout"  "Enrolled" "Graduate"

Vamos analisar a distribuição pelas classes da variável target, versus a variavel Unemployment.rate, Inflation.rate e Age.at.enrollment:

• Unemployment.rate:

boxplot(Unemployment.rate ~ target)

Acima podemos notar os grupos Dropout e Graduate aparentam ter comportamento semelhante com respetio a variável Unemployment.rate, porém ambos tem comportamento distinto da classe Enrolled, sendo esta ultima apresentando um menor intervalo interquartilico, aparentando ter uma variação menor dentre os membros da sua classe, com respeito a taxa de desemprego. Vale ressaltar também que a mediana manteve-se igual para os 3 grupos.

boxplot(Inflation.rate ~ target)

Podemos notar que para o cenário acima, agora a classe Graduate apresenta uma mudança quando comparado as outras duas classes, evidenciando uma mediana inferior quando comparado aos outros grupos, com respeito a variável inflation.rate.

names(data)

##  [1] "ï..Marital.status"           "Application.mode"           
##  [3] "Application.order"           "Course"                     
##  [5] "Daytime.evening.attendance." "Previous.qualification"     
##  [7] "Nacionality"                 "Mother.s.qualification"     
##  [9] "Father.s.qualification"      "Mother.s.occupation"        
## [11] "Father.s.occupation"         "Displaced"                  
## [13] "Educational.special.needs"   "Debtor"                     
## [15] "Tuition.fees.up.to.date"     "Gender"                     
## [17] "Scholarship.holder"          "Age.at.enrollment"          
## [19] "International"               "Target"

boxplot(Age.at.enrollment ~ target)

Podemos notar que agora com respeito a variável Age.at.enrollment um comportamento diferente entre as 3 classes, notando que o menor valor encontra-se na classe Graduate, e também a presença maior de outliers, analogo podemos notar que a classe Dropout possui em sua predominância valores mais altos. Esta variável referencia a idade na hora de inscrição do curso, podemos notar que estudantes que se inscrevem com idades menores se enquadram na classe de Graduate, ou seja, se graduam de fato, porém jovens com idades maiores ou mais avançadas que se inscrevem no curso, tendem a realizar o Dropout do curso, ou seja, desistem ou abandonam o curso.

• Trocando o nome das variáveis:

Vamos substituir os nomes das classes da variável target por valores numericos 0, 1 e 2, respectivamente

levels(target) <- c(0,1,2)
levels(target)

## [1] "0" "1" "2"

Proporção das classes da variável target

round((table(target)/length(target))*100,2)

## target
##     0     1     2 
## 32.12 17.95 49.93

Podemos notar que existe uma desproporção entre as classes, não estando balanceadas, sendo a classe majoritaria a classe 2, 49.93% dos dados, Graduate, e a minoritária é a classe 1, Enrolled, com 17.95% dos dados.

Padronizando variaveis continuas

Agora vamos padronizar as variáveis continuas na base de dados, por apresentarem escalas distintas.

data_c_std = data.frame(scale(data_c))
knitr::kable(head(data_c_std))

round(colMeans(data_c_std),4)

##  PrevQualifiGrade    AdmissionGrade      X1stCredited      X1stEnrolled 
##                 0                 0                 0                 0 
##   X1stEvaluations      X1stApproved         X1stGrade    X1stWithoutEva 
##                 0                 0                 0                 0 
##      X2ndCredited      X2ndEnrolled   X2ndEvaluations      X2ndApproved 
##                 0                 0                 0                 0 
##         X2ndGrade    X2ndWithoutEva Unemployment.rate    Inflation.rate 
##                 0                 0                 0                 0 
##               GDP 
##                 0

Selecionando K pelo Elbow Plot

set.seed(13)
p1 <- elbow.plot(data_c_std)
p1

set.seed(13)
p2 <- elbow.plot(data_c_std, alg = "cmeans")
p2

set.seed(13)
p3 <- elbow.plot(data_c_std, alg = "hclust")
p3

n_cluster_p1 = p1$data$k[p1$data$d==max(p1$data$d)]
n_cluster_p2 = p2$data$k[p2$data$d==max(p2$data$d)]
n_cluster_p3 = p3$data$k[p3$data$d==max(p3$data$d)]
n_clusters_elbow = data.frame('kmeans'=n_cluster_p1,
                              'fuzzy.kmeans'=n_cluster_p2,
                              'hierarchical'=n_cluster_p3,
                              row.names = 'elbow')
knitr::kable(n_clusters_elbow)

grid.arrange(p1, p2, p3, nrow = 2)

Selecionando K pela maximização da medida de Silhueta

(clusters_km_sil = fviz_nbclust(data_c_std, kmeans, method = "silhouette"))

n_c_km_sil = clusters_km_sil$data$clusters[clusters_km_sil$data$y == max(clusters_km_sil$data$y)]

(clusters_fuzzy_sil = fviz_nbclust(data_c_std, cmeans, method = "silhouette"))

n_c_fuzzy_sil = clusters_fuzzy_sil$data$clusters[clusters_fuzzy_sil$data$y == max(clusters_fuzzy_sil$data$y)]

(clusters_hier_sil = fviz_nbclust(data_c_std, hcut, method = "silhouette"))

n_c_hier_sil = clusters_hier_sil$data$clusters[clusters_hier_sil$data$y == max(clusters_hier_sil$data$y)]

n_clusters_silhuete = data.frame('kmeans'=n_c_km_sil,
                                 'fuzzy.kmeans'=n_c_fuzzy_sil, 
                                 'hierarchical'=n_c_hier_sil,
                                 row.names = 'silhuete')
knitr::kable(n_clusters_silhuete)

grid.arrange(clusters_km_sil, clusters_fuzzy_sil, clusters_hier_sil, nrow = 2)

K-Means

set.seed(13)
fit.kmeans <- kmeans(data_c_std, as.integer(n_clusters_silhuete$kmeans)) # Ajustar K de acordo com o gráfico
#fit.kmeans
g1 <- cluster_viz(data_c_std, as.factor(fit.kmeans$cluster), geom = "point",
                  palette= c("#00AFBB", "#E7B800","#FC4E07"),
                  main = "Cluster plot for the K-means method")

Fuzzy-KMeans

set.seed(13)
fit.cmeans <- cmeans(data_c_std, as.integer(n_clusters_silhuete$fuzzy.kmeans)) # Ajustar K de acordo com o gráfico
#fit.cmeans
g2 <- cluster_viz(data_c_std, as.factor(fit.cmeans$cluster), geom = "point",
                  palette= c("#00AFBB","#B266FF"),
                  main = "Cluster plot for the Fuzzy K-means method")

Hierarchical

set.seed(13)
fit.hclust <- hcut(data_c_std, as.integer(n_clusters_silhuete$hierarchical)) # Ajustar K de acordo com o gráfico
#fit.hclust
g3 <- cluster_viz(data_c_std, as.factor(fit.hclust$cluster), geom = "point",
                  palette= c("#00AFBB", "#E7B800","#FC4E07","#B266FF"),
                  main = "Cluster plot for the Hierarchical method")

Comparando todos

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

Acima, considerando o metodo KMeans, temos uma separação siginficante entre os 3 grupos, porém ainda com uma leve interseção entre os poligos da região dos grupos. Para o modelo Fuzzy, temos dois grupos muito bem separados, porém ainda com uma leve interseção entre seus poligonos. Por fim, temos o metodo hierarquico que foram encontrados 4 grupos, apresentando aparentemente uma grande interseção entre todos os grupos, aparentando não estar tão bem separado. Vale ressaltar que os pontos mais distantes, com respeito a componente 1, foram agrupados em um grupo diferente para o metodo hierarquico, quando comparamos o mesmo com os outros dois metodos Kmeans, fuzzy e normal.

Gerando PCA

pca1 =  prcomp(data_c_std);P = pca1$rotation ;sdev = pca1$sdev 

fviz_pca_ind(pca1,
             col.ind = data$Target,
             palette = c("#00AFBB", "#E7B800", "#FC4E07"),
             addEllipses = TRUE, 
             ellipse.type = "confidence",
             legend.title = "Engine shape",
             repel = TRUE
)

Acima podemos rever um gráfico feito na primeira prova, aonde temos as duas primeiras componentes obtidas pelo PCA, e separado pelos grupos de interesse, Dropout, Enrolled e Graduate. Dessa forma, como visto nas outras avaliações, notamos que o grupo enrolled é um grupo “indeciso”, que aparentemente apresenta caracteristicas similares em ambos os grupos, tornando uma tarefa de agrupamento bem complicada para separar esses 3 grupos. Para tanto, vamos mudar o objetivo, vamos mesclar as classes Enrolled e Graduate em uma classe “Others” e mantermos a classe Dropout, dessa forma teremos dois grupos, e assim rever o PCA acima separado por esses 2 grandes grupos.

Y = factor(data$Target)
levels(Y) = c("Dropout","Other","Other")
pca_fviz = fviz_pca_ind(pca1,
             col.ind = Y,
             palette = c("red","#B266FF"),
             addEllipses = TRUE, 
             ellipse.type = "confidence",
             legend.title = "Engine shape",
             repel = TRUE
)
pca_fviz

Agora com 2 grupos bem separados, vamos comparar os resultados obtidos pelos outros modelos com os rotulos pre-existentes.

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

## Warning: ggrepel: 4424 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

#expand.grid(pca_fviz, g1,g2,g3)

Visto que agora temos dois rotulos reais na base, “Others” e “Dropout”, temos dois grupos bem separados, porém os modelos Hierarquico e K-means classico identificaram mais do que 2 grupos. Dessa forma, afim de comparação graficas e para metricas futuras, iremos retreinar os modelos e considerar 2 clusters.

Reajustando modelos

set.seed(13)
fit.kmeans.f <- kmeans(data_c_std, 2) # Ajustar K de acordo com o gráfico
#fit.kmeans
g1.f <- cluster_viz(data_c_std, as.factor(fit.kmeans.f$cluster), geom = "point",
                  palette= c("red","#B266FF"),
                  main = "Cluster plot for the K-means method")

set.seed(13)
fit.cmeans.f <- cmeans(data_c_std, 2) # Ajustar K de acordo com o gráfico
#fit.cmeans
g2.f <- cluster_viz(data_c_std, as.factor(fit.cmeans.f$cluster), geom = "point",
                  palette= c("red","#B266FF"),
                  main = "Cluster plot for the Fuzzy K-means method")

set.seed(13)
fit.hclust.f <- hcut(data_c_std, 2) # Ajustar K de acordo com o gráfico
#fit.hclust
g3.f <- cluster_viz(data_c_std, as.factor(fit.hclust.f$cluster), geom = "point",
                  palette= c("red","#B266FF"),
                  main = "Cluster plot for the Hierarchical method")

Comparando modelos reajustados

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

## Warning: ggrepel: 4424 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

Acima podemos notar a separação dos dois grupos, para os rotulos reais e para os agrupamentos de cada modelo. Graficamente, podemos notar que o pior modelo, quando comparado aos rotulos reais, é o modelo hierarquico, super agrupando a classe 2, referente ao grupo “Others”, e sub agrupando o grupo 1, “Dropout”. Aparentemente, o modelo que melhor separou os dois grupos foi o modelo Fuzzy K-Means.

Vamos agora, primeiramente para o “K” otimo encontrado pelos modelos e depois para os dois grupos finais, “Others” e “Dropout”, verificar qual dos modelos melhor agrupou os dados.

Conectividade

c1 <- connectivity(clusters = fit.kmeans$cluster, Data = data_c_std)
c2 <- connectivity(clusters = fit.cmeans$cluster, Data = data_c_std)
c3 <- connectivity(clusters = fit.hclust$cluster, Data = data_c_std)

c1.f <- connectivity(clusters = fit.kmeans.f$cluster, Data = data_c_std)
c2.f <- connectivity(clusters = fit.cmeans.f$cluster, Data = data_c_std)
c3.f <- connectivity(clusters = fit.hclust.f$cluster, Data = data_c_std)

Dunn

d1 <- dunn(clusters = fit.kmeans$cluster, Data = data_c_std)
d2 <- dunn(clusters = fit.cmeans$cluster, Data = data_c_std)
d3 <- dunn(clusters = fit.hclust$cluster, Data = data_c_std)

d1.f <- dunn(clusters = fit.kmeans.f$cluster, Data = data_c_std)
d2.f <- dunn(clusters = fit.cmeans.f$cluster, Data = data_c_std)
d3.f <- dunn(clusters = fit.hclust.f$cluster, Data = data_c_std)

Silhueta

s.kmeans <- silhouette(fit.kmeans$cluster, dist = dist(data_c_std)^2)
s.cmeans <- silhouette(fit.cmeans$cluster, dist = dist(data_c_std)^2)
s.hclust <- silhouette(fit.hclust$cluster, dist = dist(data_c_std)^2)

s1 <- summary(s.kmeans)$avg.width
s2 <- summary(s.cmeans)$avg.width
s3 <- summary(s.hclust)$avg.width

s.kmeans.f <- silhouette(fit.kmeans.f$cluster, dist = dist(data_c_std)^2)
s.cmeans.f <- silhouette(fit.cmeans.f$cluster, dist = dist(data_c_std)^2)
s.hclust.f <- silhouette(fit.hclust.f$cluster, dist = dist(data_c_std)^2)

s1.f <- summary(s.kmeans.f)$avg.width
s2.f <- summary(s.cmeans.f)$avg.width
s3.f <- summary(s.hclust.f)$avg.width

Indice de Rand

Como temos conhecimento dos rotulos verdadeiros dos dados, vamos utilizar do indice de rand também a nivel de comparação da similaridade entre os grupos, verdadeiro e agrupado pelo respectivo modelo.

Y_f = Y
levels(Y_f)=c(1,2)
Y_f = as.integer(Y_f)
kmeans_rand = rand.index(Y_f, fit.kmeans.f$cluster)
cmeans_rand = rand.index(Y_f, fit.cmeans.f$cluster)
hclust_rand = rand.index(Y_f, fit.hclust.f$cluster)

rand_results = t(data.frame('kmeans'=kmeans_rand,
                          'fuzzy'=cmeans_rand,
                          'hierarchical'=hclust_rand,row.names = 'rand'))
knitr::kable(rand_results)