The PCA method is one of the methods in multivariate analysis used to reduce variables into certain dimensions or components (Jolliffe, 2014). From various dimensions or variables in assessing welfare, this study aims to determine the main components formed from variables in the Better Life Index dimension according to the OECD in the case of G20 member countries. This study also aims to determine the position of Indonesia in the cluster.

# library
library(knitr)
library(REdaS)
library(factoextra) 
library(clValid)
library(tidyverse)
library(cluster)
library(ggplot2)
library(tidyverse)
library(grDevices)
library(dplyr)
library(fclust)
library(psych)

Analysis begins with reading and getting to know the data. This process is also called descriptive analysis.

Data and Variable Sources

library(readxl)
dataset <- read_excel("STATA.xlsx", sheet = "Without HC_With BS",
                      col_names = TRUE) # dataset analisis
dataset <- as.data.frame(dataset)
rownames(dataset) <- dataset$Negara
dataset <- dataset[,-1]
colnames(dataset)=c("Basic Sanitation", "Median Income",
                    "Unemployment", "Years of Schooling", "Air Pollution",
                    "Regulatory Quality", "Life Expectancy", "Happiness Index",
                    "Homicide Rate", "Working Hours")
library(gt)
library(gtExtras)

Variabel <- c("Basic Sanitation", "Median Income",
              "Unemployment", "Years of Schooling", "Air Pollution",
              "Regulatory Quality", "Life Expectancy", "Happiness Index",
              "Homicide Rate", "Working Hours")
Unit <- c("percent (%)", "USD, constant 2015", "percent (%)", "Years",
         "Micrograms per cubic meter", "Score", "Years", "Score", "per 100,000 people",
         "hours per week")
Source <- c("World Bank", "World Bank", "World Bank", "UNDP", "World Bank",
            "World Bank", "World Health Organization", "World Happiness Report",
            "World Health Organization", "ILO")
df1 <- data.frame(Variabel, Unit, Source)

df1 %>% 
  gt() %>% 
  gt_theme_538() # generate tabel
Variabel Unit Source
Basic Sanitation percent (%) World Bank
Median Income USD, constant 2015 World Bank
Unemployment percent (%) World Bank
Years of Schooling Years UNDP
Air Pollution Micrograms per cubic meter World Bank
Regulatory Quality Score World Bank
Life Expectancy Years World Health Organization
Happiness Index Score World Happiness Report
Homicide Rate per 100,000 people World Health Organization
Working Hours hours per week ILO
# Visualisasi
library(corrplot)
library(viridis)
C <- cor(dataset)
library(RColorBrewer)
corrplot(C, method = "color", outline = T, addgrid.col = "darkgray", order="alphabet",
         tl.col = "black", addCoef.col = "black", number.digits = 1,
         number.cex =0.75, col =colorRampPalette(c("#8B0000", "white", "grey"))(100))

boxplot(dataset, col = colorRampPalette(c("#8B0000","white","grey"))(100))

When there are several variables with different unit sizes, it is necessary to normalize the data first.

# Normalisasi
datascaled <- scale(dataset)
datascaled <- as.data.frame(datascaled)
boxplot(datascaled, col = colorRampPalette(c("#8B0000","white","grey"))(100))

Factor Analysis

Regarding observation cluster analysis, factor analysis is required first. Factor analysis consists of data sufficiency analysis and multicollinearity analysis. Data sufficiency analysis can be tested with Kaiser-Meyer Olkin (KMO) by paying attention to the Measure of Sampling Adequacy (MSA) value of each variable. Data will be processed to clustering analysis if the KMO and MSA values are greater than 0.5. Meanwhile, multicollinearity can be detected through the Bartlett-Sphericity test. If the p-value < α, then there is multicollinearity in the data. PCA will be performed if there is multicollinearity in the data (Luthfi & Wijayanto, 2021; Afira & Wijayanto, 2021).

Data Sufficiency Analysis

library(REdaS)
KMO_res <- KMOS(datascaled)
as.data.frame(KMO_res$MSA)
##                    KMO_res$MSA
## Basic Sanitation     0.7175191
## Median Income        0.7399388
## Unemployment         0.5236364
## Years of Schooling   0.6941534
## Air Pollution        0.7796938
## Regulatory Quality   0.7080481
## Life Expectancy      0.6806854
## Happiness Index      0.7675867
## Homicide Rate        0.6369711
## Working Hours        0.8330062
KMO_res$KMO
## [1] 0.7171577

Multicollinearity Analysis

bs_res <- bart_spher(datascaled) # metode Barlett test of sphericity
bs_res
##  Bartlett's Test of Sphericity
## 
## Call: bart_spher(x = datascaled)
## 
##      X2 = 284.066
##      df = 45
## p-value < 2.22e-16
Test <- c("Kaiser Meyer Olkin (KMO)", "Bartlett-Sphericity")
Value <- c(0.7171577, 2.22e-16)
factor_analysis <- data.frame(Test, Value)
factor_analysis %>% 
  gt() %>% 
  gt_theme_538()
Test Value
Kaiser Meyer Olkin (KMO) 7.171577e-01
Bartlett-Sphericity 2.220000e-16

The table above shows that there is multicollinearity in the data. This indicates the need for Principal Component Analysis.

# Principal Component Analysis
PCA <- prcomp(datascaled)
PCA
## Standard deviations (1, .., p=10):
##  [1] 2.2864420 1.2771166 0.8875965 0.8588440 0.6691782 0.6346937 0.6064621
##  [8] 0.4524031 0.3603092 0.2505839
## 
## Rotation (n x k) = (10 x 10):
##                           PC1          PC2         PC3         PC4         PC5
## Basic Sanitation   -0.3598144 -0.001805468  0.05567502 -0.37608011  0.18584605
## Median Income      -0.3935204  0.041013042  0.05616831 -0.02108076 -0.40615146
## Unemployment        0.1751425  0.418733048  0.79456401 -0.24063747 -0.01006975
## Years of Schooling -0.3242909  0.096243580  0.10826354  0.62879635  0.03963471
## Air Pollution       0.2714083 -0.424048623  0.11742210 -0.22388649 -0.58598845
## Regulatory Quality  0.2728626  0.382233032 -0.51365475 -0.30490900  0.16649947
## Life Expectancy    -0.3702465 -0.091097328 -0.00328146 -0.49486335  0.13967271
## Happiness Index    -0.3353680  0.295086018 -0.25279384 -0.08683493 -0.57326186
## Homicide Rate       0.2590961  0.559483297 -0.07463951  0.09149083 -0.25854900
## Working Hours       0.3391672 -0.281778481 -0.05980603  0.03054234 -0.11287244
##                             PC6          PC7          PC8         PC9
## Basic Sanitation   -0.547256744  0.384434069  0.246960780  0.17042837
## Median Income       0.296948914  0.007424408 -0.520619832 -0.05019177
## Unemployment        0.001315893 -0.002921670  0.001039330 -0.32186767
## Years of Schooling -0.474844955  0.172446217 -0.316092266 -0.23395310
## Air Pollution      -0.494000784 -0.179757285 -0.195202930  0.08391455
## Regulatory Quality -0.270265706 -0.118538548 -0.360165350 -0.41050763
## Life Expectancy     0.190250528  0.190236524 -0.369922944  0.09640055
## Happiness Index     0.029047961  0.019162640  0.494895712 -0.31473411
## Homicide Rate      -0.026842085  0.337233018 -0.135520717  0.63020766
## Working Hours       0.178919155  0.791051055 -0.009802288 -0.35920971
##                           PC10
## Basic Sanitation   -0.39260748
## Median Income      -0.55955638
## Unemployment        0.03221557
## Years of Schooling  0.25890539
## Air Pollution       0.13315201
## Regulatory Quality -0.09788061
## Life Expectancy     0.60966270
## Happiness Index     0.23492140
## Homicide Rate       0.09531128
## Working Hours      -0.03761905
eigen_value <- get_eigenvalue(PCA)
eigen_value %>% 
  gt() %>% 
  gt_theme_538()
eigenvalue variance.percent cumulative.variance.percent
5.2278172 52.278172 52.27817
1.6310268 16.310268 68.58844
0.7878276 7.878276 76.46672
0.7376130 7.376130 83.84285
0.4477994 4.477994 88.32084
0.4028361 4.028361 92.34920
0.3677963 3.677963 96.02716
0.2046686 2.046686 98.07385
0.1298227 1.298227 99.37208
0.0627923 0.627923 100.00000
summary(PCA)
## Importance of components:
##                           PC1    PC2     PC3     PC4     PC5     PC6     PC7
## Standard deviation     2.2864 1.2771 0.88760 0.85884 0.66918 0.63469 0.60646
## Proportion of Variance 0.5228 0.1631 0.07878 0.07376 0.04478 0.04028 0.03678
## Cumulative Proportion  0.5228 0.6859 0.76467 0.83843 0.88321 0.92349 0.96027
##                            PC8     PC9    PC10
## Standard deviation     0.45240 0.36031 0.25058
## Proportion of Variance 0.02047 0.01298 0.00628
## Cumulative Proportion  0.98074 0.99372 1.00000
# if using graph
fviz_eig(PCA, addlabels = TRUE, barfill = "grey", barcolor = "grey")

After performing PCA, based on eigenvalues greater than or equal to 1, it was found that there are two principal components in Better Life in G20 countries. Cumulatively, the two principal components explain 68.6 percent of the total variance.

PCA_df <- as.data.frame(PCA$rotation)
print(round(PCA_df[1:2],4))
##                        PC1     PC2
## Basic Sanitation   -0.3598 -0.0018
## Median Income      -0.3935  0.0410
## Unemployment        0.1751  0.4187
## Years of Schooling -0.3243  0.0962
## Air Pollution       0.2714 -0.4240
## Regulatory Quality  0.2729  0.3822
## Life Expectancy    -0.3702 -0.0911
## Happiness Index    -0.3354  0.2951
## Homicide Rate       0.2591  0.5595
## Working Hours       0.3392 -0.2818
PCA_fix=PCA$x[,1:2]

# Graph of the variables
fviz_pca_var(PCA, col.var = "black")

# contribution of each var
fviz_cos2(PCA, choice = "var", axes = 1:2, fill = "grey", color = "grey")

# combine
fviz_pca_var(PCA, col.var = "cos2",
            gradient.cols = "uchicago",
            repel = TRUE)

The figure above shows how the variables are distributed in the components and which variables have a high quality of representation (cos2) or have a large contribution to components 1 and 2. It can be seen that the Homicide Rate and Median Income are the two variables with the highest quality of representation, which are in the interval of more than 0.7. Meanwhile, the variable with the lowest quality of representation is Unemployment.

Through this result, we can see that the first component is explained by Basic Sanitation, Median Income, Years of Schooling, Life Expectancy, Happiness Index, and Working Hours. It can be said that the first component relates to basic human needs such as health and education which can be obtained using a certain amount of money derived from individual income. This is internal to each individual or household. Hence, component 1 is named Internal Aspects. Meanwhile, the second component is explained more by Unemployment, Air Pollution, Regulatory Quality, and Homicide Rate. The variables of the second component are more directed towards something external to the life of an individual or household. Therefore, component 2 is named External Aspects.

Validity Analysis

Internal Validation

In this section we will find out which clustering method fits the data.

library(factoextra)
library(fpc)
library(NbClust)
intern <- clValid(PCA_fix, 2:10, clMethods = c("hierarchical","kmeans", "pam"), 
                  validation = "internal") # 3 metode
summary(intern)
## 
## Clustering Methods:
##  hierarchical kmeans pam 
## 
## Cluster sizes:
##  2 3 4 5 6 7 8 9 10 
## 
## Validation Measures:
##                                  2       3       4       5       6       7       8       9      10
##                                                                                                   
## hierarchical Connectivity   2.9290  6.7869 13.8464 15.7988 24.1865 26.1865 28.7444 31.6968 36.0516
##              Dunn           0.4335  0.6345  0.2788  0.2788  0.2268  0.2268  0.2581  0.2581  0.3130
##              Silhouette     0.6246  0.6453  0.5426  0.4855  0.4230  0.4023  0.3361  0.3252  0.3477
## kmeans       Connectivity   7.8365 10.7806 14.6940 19.9111 26.8258 28.8258 33.9302 37.6444 38.3143
##              Dunn           0.1075  0.1423  0.2430  0.0905  0.1312  0.1312  0.1860  0.1889  0.1939
##              Silhouette     0.5716  0.5912  0.5532  0.4714  0.4524  0.4252  0.3950  0.3955  0.3941
## pam          Connectivity   7.8365 16.4286 18.6548 23.1190 25.5690 34.8460 37.9790 39.9790 41.3734
##              Dunn           0.1075  0.0355  0.0517  0.0862  0.1507  0.0756  0.0952  0.0952  0.0952
##              Silhouette     0.5716  0.4015  0.4534  0.4053  0.4167  0.3939  0.4044  0.3772  0.3780
## 
## Optimal Scores:
## 
##              Score  Method       Clusters
## Connectivity 2.9290 hierarchical 2       
## Dunn         0.6345 hierarchical 3       
## Silhouette   0.6453 hierarchical 3

Cluster Analysis

After determining the clustering method, we will look at the cluster member states formed from the previously obtained components.

# mengukur similaritas
res_dist <- dist(datascaled, method = "euclidean") # default measurement
m <- c("single","complete","ward", "average")
names(m) <- c("single","complete","ward", "average") 
ac <- function(x){
  agnes(PCA_fix,method=x)$ac
}
Nilai_AC <- map_dbl(m,ac) # melihar hasil ac
as.data.frame(Nilai_AC)
##           Nilai_AC
## single   0.8604325
## complete 0.9265089
## ward     0.9595709
## average  0.9101382
Metode <- c("Single", "Complete", "Ward", "Average")
AC_Scores <- data.frame(Metode, Nilai_AC)
AC_Scores %>% 
  gt() %>% 
  gt_theme_538()
Metode Nilai_AC
Single 0.8604325
Complete 0.9265089
Ward 0.9595709
Average 0.9101382
library(grDevices)
data.hcc2 <- PCA_fix %>%      
  dist(method = "euclidean") %>% 
  hclust(method = "ward.D2")
fviz_dend(data.hcc2, k = 3, 
          cex = 0.6, palette = "uchicago",
          rect = TRUE, rect_fill = TRUE, rect_border = "uchicago",
          main = "Cluster Dendogram (Ward's Method)")

Referring to the previous two main components, namely Internal Aspects and External Aspects, it can be seen that Cluster 1 dominates the variable values in Internal Aspects (PC1) in this analysis. Focusing on the two variables with the highest quality of representation (cos2), namely Median Income and Homicide Rate, the summary above illustrates that Cluster 1 is a collection of countries with a relatively high median income and a relatively high level of security compared to the other two clusters. Meanwhile, Cluster 2 has a relatively poor median component score compared to the other clusters. The summary above illustrates that Cluster 2 is a group of countries with relatively low-income populations with relatively low levels of security compared to the other two clusters. Finally, the cluster to which Indonesia belongs, Cluster 3, has a relatively moderate median component value. The insight that can be gained is that Cluster 3 is a group of countries with a population that has a moderate income and a moderate level of security.

Afira, N., & Wijayanto, A. W. (2021). Analisis Cluster dengan Metode Partitioning dan Hierarki pada Data Informasi Kemiskinan Provinsi di Indonesia Tahun 2019. Komputika: Jurnal Sistem Komputer, 10(2), 101-109. https://doi.org/10.34010/komputika.v10i2.4317

Jolliffe, I. (2014). Principal Component Analysis. International Encyclopedia of Statistical Science. https://doi.org/10.1007/978-3-642-04898-2_455

Luthfi, E., & Wijayanto, A. W. (2021). Analisis perbandingan metode hirearchical, k - means, dan k - medoids clustering dalam pengelompokkan indeks pembangunan manusia Indonesia. INOVASI, 17(4), 761. https://journal.feb.unmul.ac.id/index.php/INOVASI/article/view/10106/1513