Eksplorasi Data: Hubungan Stunting dengan
Early Childhood Development Index (ECDI)
Sapta Hastho Ponco
2025-11-21
Eksplorasi Data atau Exploratory Data Analysis (EDA) merupakan langkah awal yang sangat penting sebelum melakukan analisis statistika yang lebih lanjut. Exploratory Data Analysis dilakukan untuk memahami struktur data, menemukan pola penting, mendeteksi anomali, serta memeriksa asumsi statistik sebelum dilakukan analisis lebih lanjut (Tukey, 1977).
EDA juga membantu memastikan bahwa metode analisis yang dipilih sesuai dengan karakteristik data. Misalnya disyaratkan distribusi data yang normal atau tidak boleh mengandung outlier ekstrem. Tanpa pemeriksaan awal ini, hasil analisis dapat menjadi bias atau menyesatkan.
Dalam kajian ini, dilakukan proses eksplorasi data:
ECDI_24_29: ECDI kelompok umur 24-29 bulan
ECDI_36_41: ECDI kelompok umur 36-41 bulan
ECDI_42_47: ECDI kelompok umur 42-47 bulan
ECDI_48_59: ECDI kelompok umur 48-59 bulan
ECDI_total: ECDI kelompok umur 24-59 bulan
Variabel X: Balita stunting 2024.
pacman::p_load(readxl, psych, dplyr, tidyr, survey, janitor, magrittr, car, ggplot2)1 Import Data
setwd ("D:\\Kesra ECDI 2025")
df <- read_excel("ECDI SSN MODUL 2025 - weighted + Stunting.xlsx")
str(df)## tibble [38 × 11] (S3: tbl_df/tbl/data.frame)
## $ No : num [1:38] 1 2 3 4 5 6 7 8 9 10 ...
## $ Provinsi : chr [1:38] "Aceh" "Sumatera Utara" "Sumatera Barat" "Riau" ...
## $ Kode_prov : num [1:38] 11 12 13 14 15 16 17 18 19 21 ...
## $ ECDI_24_29: num [1:38] 85.8 81.1 88.8 83.7 85.8 ...
## $ ECDI_30_35: num [1:38] 86.1 78.9 87.7 91.9 76 ...
## $ ECDI_36_41: num [1:38] 80.6 78 96.2 90.5 87.4 ...
## $ ECDI_42_47: num [1:38] 76.7 69.3 72.8 69.7 63.6 ...
## $ ECDI_48_59: num [1:38] 76.8 68.6 86.5 73.1 74.5 ...
## $ ECDI_total: num [1:38] 80.2 73.7 86.7 78.9 77.4 ...
## $ VarX : num [1:38] 28.6 22 24.9 20.1 17.1 15.9 18.8 15.9 20.1 15 ...
## $ RSE : num [1:38] 2.07 2.51 2.4 3.23 4.21 ...2 Statistik Deskriptif
describe(df %>% select(ECDI_24_29,
ECDI_30_35,
ECDI_36_41,
ECDI_42_47,
ECDI_48_59,
ECDI_total))## vars n mean sd median trimmed mad min max range skew
## ECDI_24_29 1 38 80.83 12.92 83.99 82.15 6.08 42.34 100.00 57.66 -1.15
## ECDI_30_35 2 38 81.73 8.90 83.33 82.17 9.88 59.77 97.04 37.27 -0.41
## ECDI_36_41 3 38 82.13 11.58 83.61 83.05 9.84 42.51 100.00 57.49 -1.06
## ECDI_42_47 4 38 71.40 14.71 72.58 72.05 13.06 15.01 100.00 84.99 -1.17
## ECDI_48_59 5 38 71.93 12.11 73.40 73.20 8.44 35.35 89.48 54.14 -1.01
## ECDI_total 6 38 76.25 10.38 77.57 77.21 10.32 40.83 91.17 50.35 -1.15
## kurtosis se
## ECDI_24_29 0.86 2.10
## ECDI_30_35 -0.24 1.44
## ECDI_36_41 1.80 1.88
## ECDI_42_47 3.83 2.39
## ECDI_48_59 1.01 1.96
## ECDI_total 1.86 1.68describe(df$VarX)## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 38 23.01 6.66 22.55 22.61 6.82 8.7 40 31.3 0.49 0.02 1.08Output fungsi describe(df$VarX) memberikan rangkuman
statistik deskriptif untuk Variabel X.
Dari hasil diperoleh:
n = 38 menunjukkan jumlah observasi yang tersedia dalam Variabel X.
Mean = 23.01 Rata-rata nilai Variabel X berada pada kisaran 23.
SD = 6.66. Standar deviasi sebesar 6.66 berarti bahwa sebaran nilai X1 cukup bervariasi dari nilai rata-rata. Semakin besar SD, semakin besar variasi dalam data.
Median = 22.55. Median adalah nilai tengah data setelah diurutkan. Karena median 22.55, dan mean 23.01, maka datanya sedikit miring ke kanan (right-skewed), yang dikonfirmasi oleh nilai skew di bawah.
Trimmed Mean = 22.61. Rata-rata yang dihitung setelah menghilangkan sebagian data ekstrem (biasanya 10% teratas dan terbawah). Nilainya mendekati mean, jadi outlier tidak terlalu memengaruhi data.
MAD = 6.82 Median absolute deviation menunjukkan penyebaran yang lebih robust terhadap outlier. Nilainya sekitar 6.82, mirip dengan SD, menandakan tidak ada outlier ekstrem.
Min = 8.7; Max = 40.0; Range = 31.3
Skew = 0.49. Ada beberapa nilai yang lebih tinggi yang menarik mean sedikit ke atas. Kurva menceng ke kanan.
Kurtosis = 0.02. Nilai 0.02, mendekati 0, berarti distribusi mendekati normal (mesokurtic).
3 Presisi Estimasi
Kualitas presisi hasil estimasi suatu survei bisa diamati dari nilai RSE yang dihasilkan. Kesalahan sampling dari beberapa estimasi harus digunakan secara hati-hati. Untuk estimasi yang berdasarkan jumlah kasus yang kecil, kesalahan relatif cenderung besar.
Batasan nilai RSE ( Australian Bureau Statistics):
Jika \(RSE \le 25\%\), estimasi bersifat presisi.
Jika \(25\% < RSE \le 50\%\), estimasi perlu dilakukan dengan hati-hati.
Jika \(RSE > 50\%\), estimasi dianggap sangat tidak presisi.
df <- df %>%
mutate(
RSE_kategori = case_when(
RSE < 10 ~ "< 10%",
RSE >= 10 & RSE <= 25 ~ "10–25%",
RSE > 25 & RSE <= 50 ~ "26–50%",
RSE > 50 ~ "> 50%",
TRUE ~ NA_character_
)
)
table(df$RSE_kategori)##
## < 10%
## 38Dari hasil di atas, tampak semua data diestimasi dengan nilai RSE kurang dari 10%.
Secara visual dapat digambarkan sebagai berikut:
ggplot(df, aes(x = Provinsi, y = RSE, group = 1, color = RSE_kategori)) +
geom_line(linewidth = 1.4) +
geom_point(size = 4) +
geom_text(aes(label = round(RSE,1)), vjust = -0.8, size = 3.5, color = "black") +
geom_hline(yintercept = 25, linetype = "dashed", linewidth = 1, color = "red") +
scale_color_manual(
name = "Kategori RSE",
values = c(
"< 10%" = "#7570B3",
"10–25%" = "#FF8C00",
"26–50%" = "#D95F02",
"> 50%" = "#1B9E77"
)
) +
labs(
title = "RSE per Provinsi",
subtitle = "Termasuk garis batas kualitas survei (RSE ≤ 25%)",
x = NULL,
y = "Nilai RSE (%)"
) +
theme_minimal(base_size = 13) +
theme(
plot.title = element_text(face = "bold", size = 15),
legend.position = "bottom",
axis.text.x = element_text(angle = 40, vjust = 1, hjust = 1)
)
4 Korelasi Pearson
Korelasi adalah ukuran statistik yang digunakan untuk menggambarkan derajat hubungan linier antara dua variabel. Nilai korelasi menunjukkan seberapa kuat kedua variabel, apakah meningkat secara bersamaan, bergerak berlawanan arah, atau tidak memiliki hubungan sama sekali.
Menurut Field (2013), korelasi Pearson digunakan untuk mengukur hubungan linear antara dua variabel numerik dan memiliki nilai antara -1 hingga +1.
Hipotesis:
\(H_0 : r = 0\)
\(H_1 : r \neq 0\)
Rumus koefisien korelasi Pearson:
\[ r = \frac{\sum (X_i - \bar{X})(Y_i - \bar{Y})} {\sqrt{\sum (X_i - \bar{X})^2 \cdot \sum (Y_i - \bar{Y})^2}} \]
cor_24_29 <- cor.test(df$ECDI_24_29, df$VarX, method = "pearson")
cor_30_35 <- cor.test(df$ECDI_30_35, df$VarX, method = "pearson")
cor_36_41 <- cor.test(df$ECDI_36_41, df$VarX, method = "pearson")
cor_42_47 <- cor.test(df$ECDI_42_47, df$VarX, method = "pearson")
cor_48_59 <- cor.test(df$ECDI_48_59, df$VarX, method = "pearson")
cor_total <- cor.test(df$ECDI_total, df$VarX, method = "pearson")
#cor_24_29; cor_30_35; cor_36_41; cor_42_47; cor_48_59; cor_totalMatrik Korelasi yang dihasilkan adalah:
library(corrplot)
df_corr <- df %>%
select(ECDI_24_29, ECDI_30_35, ECDI_36_41, ECDI_42_47, ECDI_48_59, ECDI_total, VarX)
corr_matrix <- cor(df_corr, use = "pairwise.complete.obs")cor_matrix <- cor(df %>%
select(ECDI_24_29,
ECDI_30_35,
ECDI_36_41,
ECDI_42_47,
ECDI_48_59,
ECDI_total,
VarX),
use = "pairwise.complete.obs")
cor_matrix## ECDI_24_29 ECDI_30_35 ECDI_36_41 ECDI_42_47 ECDI_48_59 ECDI_total
## ECDI_24_29 1.0000000 0.4628475 0.6802541 0.5142396 0.6866327 0.8100854
## ECDI_30_35 0.4628475 1.0000000 0.4703781 0.5372181 0.7173902 0.7612357
## ECDI_36_41 0.6802541 0.4703781 1.0000000 0.3190239 0.5695030 0.7012494
## ECDI_42_47 0.5142396 0.5372181 0.3190239 1.0000000 0.7057030 0.7950495
## ECDI_48_59 0.6866327 0.7173902 0.5695030 0.7057030 1.0000000 0.9372175
## ECDI_total 0.8100854 0.7612357 0.7012494 0.7950495 0.9372175 1.0000000
## VarX -0.5401172 -0.4715659 -0.4694895 -0.6448003 -0.6393618 -0.7046309
## VarX
## ECDI_24_29 -0.5401172
## ECDI_30_35 -0.4715659
## ECDI_36_41 -0.4694895
## ECDI_42_47 -0.6448003
## ECDI_48_59 -0.6393618
## ECDI_total -0.7046309
## VarX 1.0000000Secara keseluruhan, matriks ini menunjukkan bahwa VarX memiliki hubungan negatif yang kuat dan konsisten dengan seluruh indikator ECDI. Hubungan negatif tersebut tampak stabil pada semua kelompok usia, namun cenderung menguat pada rentang usia 42 bulan ke atas serta pada skor total ECDI, sehingga mengindikasikan bahwa semakin tinggi nilai VarX, semakin rendah capaian perkembangan anak usia dini.
Visualisasi Korelogram
corrplot(corr_matrix,
method = "color",
type = "upper",
addCoef.col = "black",
tl.col = "black",
tl.srt = 45,
col = colorRampPalette(c("red", "white", "blue"))(200))
Signifikansi korelasi
hasil_korelasi <- data.frame(
Variabel = c("ECDI_24_29", "ECDI_30_35", "ECDI_36_41", "ECDI_42_47", "ECDI_48_59", "ECDI_total"),
Korelasi = c(cor_24_29$estimate,
cor_30_35$estimate,
cor_36_41$estimate,
cor_42_47$estimate,
cor_48_59$estimate,
cor_total$estimate),
p_value = c(cor_24_29$p.value,
cor_30_35$p.value,
cor_36_41$p.value,
cor_42_47$p.value,
cor_48_59$p.value,
cor_total$p.value)
)
print(hasil_korelasi)## Variabel Korelasi p_value
## 1 ECDI_24_29 -0.5401172 4.650078e-04
## 2 ECDI_30_35 -0.4715659 2.801977e-03
## 3 ECDI_36_41 -0.4694895 2.942193e-03
## 4 ECDI_42_47 -0.6448003 1.244632e-05
## 5 ECDI_48_59 -0.6393618 1.553137e-05
## 6 ECDI_total -0.7046309 7.902989e-07Secara keseluruhan, nilai koefisien korelasi berada pada kategori sedang hingga kuat (−0.47 hingga −0.70), dan seluruhnya signifikan. Pola ini menunjukkan adanya hubungan yang konsisten dan bermakna antara kondisi yang diukur oleh VarX dengan indikator perkembangan anak (ECDI). Semakin tinggi nilai VarX, semakin rendah skor perkembangan anak, sehingga VarX dapat dipandang sebagai faktor yang berpotensi menghambat perkembangan anak usia dini.
5 Normality Test
Dilakukan uji normalitas Variabel X menggunakan uji Shapiro–Wilk. Uji ini digunakan untuk mengevaluasi apakah distribusi Variabel X mengikuti distribusi normal. Nilai p < 0.05 menunjukkan bahwa Variabel X tidak berdistribusi normal, sedangkan p ≥ 0.05 menunjukkan bahwa data mendekati distribusi normal.
Hipotesis
\(H_{0}\): Data berdistribusi Normal
\(H_{1}\): Data tidak berdistribusi Normal
shapiro.test(df$VarX)##
## Shapiro-Wilk normality test
##
## data: df$VarX
## W = 0.97131, p-value = 0.4281Hasil uji Shapiro–Wilk menunjukkan nilai statistik W = 0.97131 dengan p-value = 0.4281, yang berada jauh di atas batas signifikansi umum (α = 0.05).
Karena p-value lebih besar dari 0.05, tidak ada alasan untuk menolak hipotesis nol (H₀) yang menyatakan bahwa data berdistribusi normal. Dengan demikian dapat dikatakan Variabel X mengikuti distribusi normal.
Visualisasi
par(mfrow = c(1, 2))
# QQ Plot
qqnorm(df$VarX,
main = "QQ Plot Variabel X")
qqline(df$VarX, col = "red", lwd = 2)
# Histogram
hist(df$VarX,
breaks = 25,
probability = TRUE,
main = "Histogram Variabel X",
col = "skyblue",
xlab = "Var X")
lines(density(df$VarX), col = "red", lwd = 2)
par(mfrow = c(1, 1))6 Outlier Test
Untuk mengidentifikasi apakah terdapat nilai ekstrem (outlier) pada Variabel X, digunakan pendekatan Interquartile Range (IQR). Metode IQR digunakan untuk menentukan batas nilai yang dianggap wajar dalam sebuah himpunan data berdasarkan nilai kuartil pertama (Q1) dan kuartil ketiga (Q3).
Nilai yang berada di luar rentang:
\[ Q_1 - 1.5 \times IQR \quad \text{dan} \quad Q_3 + 1.5 \times IQR \]
dikategorikan sebagai outlier atau nilai pencilan karena dianggap menyimpang dari distribusi data utama (Tukey, 1977; Bruce, Miles, & Field, 2020).
boxplot(df$VarX,
main = "Boxplot Variabel X",
ylab = "Nilai Variabel X",
col = "lightblue")
Nilai yang berada di bawah 4.475 atau di atas 39.075 dianggap sebagai outlier karena berada di luar rentang variasi data. Batas-batas ini kemudian digunakan untuk mengidentifikasi provinsi mana saja yang memiliki nilai Variabel X yang menyimpang secara signifikan dari pola umum.
Pemeriksaan Outlier
Q1 <- quantile(df$VarX, 0.25, na.rm = TRUE)
Q3 <- quantile(df$VarX, 0.75, na.rm = TRUE)
IQR <- Q3 - Q1
# Batas bawah dan atas
lower_bound <- Q1 - 1.5 * IQR
upper_bound <- Q3 + 1.5 * IQR
lower_bound## 25%
## 4.475upper_bound## 75%
## 39.075Mengidentifikasi outlier
outliers_iqr <- df %>% filter(VarX < lower_bound | VarX > upper_bound)
outliers_iqr## # A tibble: 1 × 12
## No Provinsi Kode_prov ECDI_24_29 ECDI_30_35 ECDI_36_41 ECDI_42_47
## <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 38 Papua Pegunungan 97 42.3 64.1 74.1 15.0
## # ℹ 5 more variables: ECDI_48_59 <dbl>, ECDI_total <dbl>, VarX <dbl>,
## # RSE <dbl>, RSE_kategori <chr>Nilai estimasi yang dianggap outlier:
bp <- boxplot(df$VarX, plot = FALSE)
bp$out## [1] 40Posisi data dari nilai estimasi yang outlier:
detect_outlier <- function(x){
Q1 <- quantile(x, 0.25, na.rm = TRUE)
Q3 <- quantile(x, 0.75, na.rm = TRUE)
IQR <- Q3 - Q1
lower <- Q1 - 1.5 * IQR
upper <- Q3 + 1.5 * IQR
return(which(x < lower | x > upper))
}
detect_outlier(df$VarX)## [1] 38Berdasarkan perhitungan tersebut, terdapat satu provinsi yang teridentifikasi sebagai outlier, yaitu Provinsi Papua Pegunungan (kode 97). Provinsi ini memiliki nilai VarX sebesar 40, yang berada di luar batas atas rentang wajar, sehingga dikategorikan sebagai observasi yang menyimpang secara signifikan dari pola umum.
7 References
[1] Bruce, P., Miles, J., & Field, Z. (2020). Practical Statistics for Data Scientists: 50+ Essential Concepts Using R and Python. O’Reilly Media.
[2] Field, A. (2013). Discovering statistics using IBM SPSS statistics (4th ed.). SAGE Publications.
[3] Confidentiality and relative standard error, Australian Bureau Statistics. https://www.abs.gov.au/statistics/microdata-tablebuilder/tablebuilder/confidentiality-and-relative-standard-error
[4] Tukey, J. W. (1977). Exploratory Data Analysis. Reading, MA: Addison-Wesley.
Direktorat Statistik Kesejahteraan Rakyat, BPS, saptahas@bps.go.id