library(readxl)
library(tidyverse)
## āā Attaching core tidyverse packages āāāāāāāāāāāāāāāāāāāāāāāā tidyverse 2.0.0 āā
## ā dplyr 1.1.4 ā readr 2.1.5
## ā forcats 1.0.0 ā stringr 1.5.1
## ā ggplot2 4.0.0 ā tibble 3.3.0
## ā lubridate 1.9.4 ā tidyr 1.3.1
## ā purrr 1.1.0
## āā Conflicts āāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāā tidyverse_conflicts() āā
## ā dplyr::filter() masks stats::filter()
## ā dplyr::lag() masks stats::lag()
## ā¹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(plm)
##
## Attaching package: 'plm'
##
## The following objects are masked from 'package:dplyr':
##
## between, lag, lead
library(gplots)
##
## Attaching package: 'gplots'
##
## The following object is masked from 'package:stats':
##
## lowess
library(stargazer)
##
## Please cite as:
##
## Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
## R package version 5.2.3. https://CRAN.R-project.org/package=stargazer
library(dplyr)
library(ggplot2)
data = read_excel("D:\\Semester 6\\ADP\\DATA PANEL_ (1).xlsx")
## New names:
## ⢠`` -> `...10`
## ⢠`` -> `...11`
## ⢠`` -> `...12`
str(data)
## tibble [190 Ć 12] (S3: tbl_df/tbl/data.frame)
## $ Provinsi: chr [1:190] "Aceh" NA NA NA ...
## $ Tahun : num [1:190] 2020 2021 2022 2023 2024 ...
## $ Y : chr [1:190] "31.8" "30.96" "26.6" "23.68" ...
## $ X1 : num [1:190] 41.7 42.1 36.6 42.6 43.7 ...
## $ X2 : chr [1:190] "0.51400000000000001" "0.503" "0.504" "0.48899999999999999" ...
## $ X3 : chr [1:190] "9.1300000000000008" "9.2200000000000006" "9.24" "9.4" ...
## $ X4 : chr [1:190] "131580.97" "135274.04" "140971.72" "146932.42000000001" ...
## $ X5 : chr [1:190] "0.34" "0.35" "0.33" "0.33" ...
## $ X6 : chr [1:190] "0.74" "0.56999999999999995" "0.5" "0.38" ...
## $ ...10 : logi [1:190] NA NA NA NA NA NA ...
## $ ...11 : chr [1:190] NA NA NA NA ...
## $ ...12 : chr [1:190] NA NA NA NA ...
# isi kolom provinsi yang NA degan nilai di atasnya
data = data %>% fill(Provinsi, .direction = "down")
head(data)
## # A tibble: 6 Ć 12
## Provinsi Tahun Y X1 X2 X3 X4 X5 X6 ...10 ...11 ...12
## <chr> <dbl> <chr> <dbl> <chr> <chr> <chr> <chr> <chr> <lgl> <chr> <chr>
## 1 Aceh 2020 31.8 41.7 0.51⦠9.13⦠1315⦠0.34 0.74 NA <NA> <NA>
## 2 Aceh 2021 30.96 42.1 0.503 9.22⦠1352⦠0.35 0.56⦠NA <NA> <NA>
## 3 Aceh 2022 26.6 36.6 0.504 9.24 1409⦠0.33 0.5 NA <NA> <NA>
## 4 Aceh 2023 23.68 42.6 0.48⦠9.4 1469⦠0.33 0.38 NA <NA> <NA>
## 5 Aceh 2024 21.59 43.7 0.45⦠9.49 1537⦠0.33 0.42 NA <NA> <NA>
## 6 Sumatera Ut⦠2020 33.03 51.7 0.46⦠9.27⦠5337⦠0.56⦠1.09⦠NA <NA> <NA>
data = data %>%
select(-starts_with("...")) %>%
mutate(
`Y` = as.numeric(`Y`),
`X1` = as.numeric(`X1`),
`X2` = as.numeric(`X2`),
`X3` = as.numeric(`X3`),
`X4` = as.numeric(`X4`),
`X5` = as.numeric(`X5`),
`X6` = as.numeric(`X6`)
)
## Warning: There were 6 warnings in `mutate()`.
## The first warning was:
## ā¹ In argument: `Y = as.numeric(Y)`.
## Caused by warning:
## ! NAs introduced by coercion
## ā¹ Run `dplyr::last_dplyr_warnings()` to see the 5 remaining warnings.
str(data)
## tibble [190 Ć 9] (S3: tbl_df/tbl/data.frame)
## $ Provinsi: chr [1:190] "Aceh" "Aceh" "Aceh" "Aceh" ...
## $ Tahun : num [1:190] 2020 2021 2022 2023 2024 ...
## $ Y : num [1:190] 31.8 31 26.6 23.7 21.6 ...
## $ X1 : num [1:190] 41.7 42.1 36.6 42.6 43.7 ...
## $ X2 : num [1:190] 0.514 0.503 0.504 0.489 0.459 0.468 0.445 0.442 0.425 0.399 ...
## $ X3 : num [1:190] 9.13 9.22 9.24 9.4 9.49 9.28 9.33 9.5 9.63 9.77 ...
## $ X4 : num [1:190] 131581 135274 140972 146932 153780 ...
## $ X5 : num [1:190] 0.34 0.35 0.33 0.33 0.33 0.57 0.58 0.11 0 0 ...
## $ X6 : num [1:190] 0.74 0.57 0.5 0.38 0.42 1.09 0.93 0.82 0.71 1.06 ...
head(data)
## # A tibble: 6 Ć 9
## Provinsi Tahun Y X1 X2 X3 X4 X5 X6
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Aceh 2020 31.8 41.7 0.514 9.13 131581. 0.34 0.74
## 2 Aceh 2021 31.0 42.1 0.503 9.22 135274. 0.35 0.57
## 3 Aceh 2022 26.6 36.6 0.504 9.24 140972. 0.33 0.5
## 4 Aceh 2023 23.7 42.6 0.489 9.4 146932. 0.33 0.38
## 5 Aceh 2024 21.6 43.7 0.459 9.49 153780. 0.33 0.42
## 6 Sumatera Utara 2020 33.0 51.7 0.468 9.28 533746. 0.57 1.09
# Statistika Deskriptif
summary(data)
## Provinsi Tahun Y X1
## Length:190 Min. :2020 Min. :16.04 Min. :34.54
## Class :character 1st Qu.:2021 1st Qu.:29.72 1st Qu.:44.15
## Mode :character Median :2022 Median :33.62 Median :47.37
## Mean :2022 Mean :33.61 Mean :48.47
## 3rd Qu.:2023 3rd Qu.:38.29 3rd Qu.:51.73
## Max. :2024 Max. :46.28 Max. :81.34
## NA's :16 NA's :16
## X2 X3 X4 X5
## Min. :0.1420 Min. : 3.400 Min. : 13172 Min. :0.0000
## 1st Qu.:0.4235 1st Qu.: 7.872 1st Qu.: 64008 1st Qu.:0.3425
## Median :0.4745 Median : 8.550 Median : 141776 Median :0.4650
## Mean :0.4568 Mean : 8.523 Mean : 333834 Mean :0.5018
## 3rd Qu.:0.5218 3rd Qu.: 9.195 3rd Qu.: 356561 3rd Qu.:0.5900
## Max. :0.6650 Max. :11.190 Max. :2151041 Max. :1.3300
## NA's :16 NA's :16 NA's :12 NA's :20
## X6
## Min. :0.380
## 1st Qu.:0.930
## Median :1.160
## Mean :1.183
## 3rd Qu.:1.410
## Max. :2.360
## NA's :17
pdata = pdata.frame(data, index = c ("Provinsi", "Tahun"))
class(pdata)
## [1] "pdata.frame" "data.frame"
pdim(pdata)
## Balanced Panel: n = 38, T = 5, N = 190
# Cek Heterogentitas
plotmeans(Y ~ Provinsi, data = data,
xlab = "Provinsi", ylab = "Y",
main = "Heterogenitas Antar Provinsi",
bars = TRUE, p = 0.95, minbar = TRUE, mean.labels = FALSE,
xaxt = "n")
## Warning in qt((1 + p)/2, ns - 1): NaNs produced
# Label sumbu X dibuat miring agar terbaca
axis(1, at = 1:length(unique(data$Provinsi)),
labels = unique(data$Provinsi), las = 2, cex.axis = 0.7)
# Cek jumlah data valid (bukan NA) per Provinsi
cek_data = data %>%
filter(!is.na(Y)) %>% # Hanya ambil data Y yang ada isinya
group_by(Provinsi) %>%
summarise(
Jumlah_Data = n(), # Hitung berapa tahun datanya
Variasi = var(Y) # Hitung varians (keragaman)
) %>%
arrange(Jumlah_Data) # Urutkan dari yang paling sedikit
# Tampilkan hasil
print(cek_data, n = 34)
## # A tibble: 38 Ć 3
## Provinsi Jumlah_Data Variasi
## <chr> <int> <dbl>
## 1 Papua Barat Daya 1 NA
## 2 Papua Pegunungan 1 NA
## 3 Papua Selatan 1 NA
## 4 Papua Tengah 1 NA
## 5 Aceh 5 19.8
## 6 Bali 5 11.4
## 7 Banten 5 45.8
## 8 Bengkulu 5 14.6
## 9 DI Yogyakarta 5 5.17
## 10 DKI Jakarta 5 42.4
## 11 Gorontalo 5 21.6
## 12 Jambi 5 12.2
## 13 Jawa Barat 5 34.3
## 14 Jawa Tengah 5 7.11
## 15 Jawa Timur 5 5.73
## 16 Kalimantan Barat 5 13.4
## 17 Kalimantan Selatan 5 26.2
## 18 Kalimantan Tengah 5 24.3
## 19 Kalimantan Utara 5 27.4
## 20 Kalimantan timur 5 20.7
## 21 Kep. Bangka Belitung 5 28.6
## 22 Kep. Riau 5 10.9
## 23 Lampung 5 17.5
## 24 Maluku 5 12.6
## 25 Maluku Utara 5 19.1
## 26 NTB 5 23.9
## 27 NTT 5 11.8
## 28 Papua 5 91.8
## 29 Papua Barat 5 14.9
## 30 Riau 5 35.0
## 31 Sulawesi Barat 5 19.7
## 32 Sulawesi Selatan 5 16.3
## 33 Sulawesi Tengah 5 11.5
## 34 Sulawesi Tenggara 5 15.6
## # ā¹ 4 more rows
# Visualisasi Tren Waktu
ggplot(data, aes(x = Tahun, y = Y, group = Provinsi)) +
geom_line(color = "grey", alpha = 0.6) +
geom_point(color = "grey", size = 1, alpha = 0.6) +
stat_summary(fun = mean, geom = "line", aes(group = 1), color = "red", size = 1.2) +
theme_minimal() +
labs(
title = "Tren Persentase Pemuda Kawin (20220-2024)",
subtitle ="Garis Abu: per provinsi | Garis Merah: Rataan Nasional",
y = "Persentase Pemuda Kawin",
x = "Tahun"
)
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ā¹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning: Removed 16 rows containing non-finite outside the scale range
## (`stat_summary()`).
## Warning: Removed 16 rows containing missing values or values outside the scale range
## (`geom_line()`).
## Warning: Removed 16 rows containing missing values or values outside the scale range
## (`geom_point()`).
garis merah menandakan bahwa rata-rata presentase pemuda kawin di tiap provinsi di Indonesia cenderung menurun dari waktu ke waktu.
ada beberapa garis abu yang menaik tajam/menurun tajam, bisa jadi perhatian lebih untuk didalami.
# 1. Ubah struktur data X menjadi memanjang ke bawah (Long Format)
# Ini trik agar kita bisa mem-plot semuanya sekaligus
data_long = data %>%
select(Y, X1, X2, X3, X4, X5, X6) %>%
pivot_longer(
cols = starts_with("X"), # Pilih semua kolom yang diawali huruf X
names_to = "Variabel_X", # Nama kolom baru untuk label (X1, X2, dst)
values_to = "Nilai_X" # Nama kolom baru untuk angkanya
)
# 2. Buat Plot Gabungan (Facet Wrap)
ggplot(data_long, aes(x = Nilai_X, y = Y)) +
geom_point(alpha = 0.4, color = "darkblue", size = 1) + # Titik data
geom_smooth(method = "lm", color = "red", se = FALSE) + # Garis regresi linier
facet_wrap(~Variabel_X, scales = "free_x", ncol = 3) + # Memecah gambar per variabel X
theme_minimal() +
labs(
title = "Hubungan Y (Pemuda Kawin) dengan Semua Variabel X",
subtitle = "Garis Merah: Tren Regresi Linier",
y = "Persen Pemuda Kawin (Y)",
x = "Nilai Variabel X"
) +
theme(strip.text = element_text(face = "bold", size = 10))
## `geom_smooth()` using formula = 'y ~ x'
## Warning: Removed 101 rows containing non-finite outside the scale range
## (`stat_smooth()`).
## Warning: Removed 101 rows containing missing values or values outside the scale range
## (`geom_point()`).
Garis merah merupakan slope yang menunjukkan hubungan Y dengan variabel peubah
A. Hubungan Negatif
Y vs X1 = Y (Persen Pemuda Kawin) meningkat, maka X1 (TPAK Perempuan) menurun.
Y vs X3 = Y (Persen Pemuda Kawin) meningkat, maka X3 (RLS Perempuan) menurun.
Y vs X4 = Y (Persen Pemuda Kawin) meningkat, maka X4 (PDRB Konstan) menurun.
Y vs X5 = Y (Persen Pemuda Kawin) meningkat, maka X5 (Urbanisasi) menurun.
B. Hubungan Positif
Y vs X2 = Y (Persen Pemuda Kawin) menurun, maka X2 (IKG) meningkat.
Y vs X6 = Y (Persen Pemuda Kawin) menurun, maka X6 (Cerai) meningkat.
# Membuang Provinsi yang datanya < 2 Tahun
data_clean = data %>%
group_by(Provinsi) %>%
filter(n() > 1) %>% # hanya ambil provinsi dengan data > 1
ungroup()
table(data_clean$Provinsi)
##
## Aceh Bali Banten
## 5 5 5
## Bengkulu DI Yogyakarta DKI Jakarta
## 5 5 5
## Gorontalo Jambi Jawa Barat
## 5 5 5
## Jawa Tengah Jawa Timur Kalimantan Barat
## 5 5 5
## Kalimantan Selatan Kalimantan Tengah Kalimantan timur
## 5 5 5
## Kalimantan Utara Kep. Bangka Belitung Kep. Riau
## 5 5 5
## Lampung Maluku Maluku Utara
## 5 5 5
## NTB NTT Papua
## 5 5 5
## Papua Barat Papua Barat Daya Papua Pegunungan
## 5 5 5
## Papua Selatan Papua Tengah Riau
## 5 5 5
## Sulawesi Barat Sulawesi Selatan Sulawesi Tengah
## 5 5 5
## Sulawesi Tenggara Sulawesi Utara Sumatera Barat
## 5 5 5
## Sumatera Selatan Sumatera Utara
## 5 5
# Plot Korelasi: Melihat hubungan antar variabel X
library(corrplot)
## corrplot 0.95 loaded
# Ambil hanya kolom numerik
corr_data = data_clean %>%
select(Y, X1, X2, X3, X4, X5, X6) %>%
na.omit() # Buang NA sebentar untuk hitung korelasi
# Hitung matriks
M = cor(corr_data)
# Visualisasi yang cantik
corrplot(M, method = "number", type = "upper",
tl.col = "black", tl.srt = 45,
title = "Matriks Korelasi", mar=c(0,0,2,0))
# Buat ranking Provinsi
ranking_provinsi = data_clean %>%
group_by(Provinsi) %>%
summarise(Rata_Rata_Y = mean(Y, na.rm = TRUE)) %>%
arrange(desc(Rata_Rata_Y))
# Barplot 5 Provinsi dengan Persentase Pemuda Kawin tertinggi dan Terendah
top_bottom = rbind(
head(ranking_provinsi, 5) %>% mutate(Status = "Tertinggi"),
tail(ranking_provinsi, 5) %>% mutate(Status = "Terendah")
)
# 2. Buat Plot
ggplot(top_bottom, aes(x = reorder(Provinsi, Rata_Rata_Y), y = Rata_Rata_Y, fill = Status)) +
geom_col() +
coord_flip() +
scale_fill_manual(values = c("Tertinggi" = "red", "Terendah" = "darkgreen")) +
geom_text(aes(label = round(Rata_Rata_Y, 1)), hjust = -0.2) +
theme_minimal() +
labs(
title = "Kontras Ekstrem: Provinsi Tertinggi vs Terendah",
y = "Persentase Pemuda Kawin",
x = ""
)
5 Provinsi dengan Persentase Pemuda Kawin Tertinggi
5 Provinsi dengan Persentase Pemuda Kawin Terendah
# Plot Korelasi: Melihat hubungan antar variabel X
library(corrplot)
# Ambil hanya kolom numerik
corr_data = data_clean %>%
select(X1, X2, X3, X4, X5, X6) %>%
na.omit() # Buang NA sebentar untuk hitung korelasi
# Hitung matriks
M = cor(corr_data)
# Visualisasi yang cantik
corrplot(M, method = "number", type = "upper",
tl.col = "black", tl.srt = 45,
title = "Matriks Korelasi", mar=c(0,0,2,0))
TPAK Perempuan (\(X1\)) vs IKG (\(X2\))= Korelasi Negatif kuat
Pendidikan (\(X3\) ) vs IKG (\(X2\)) = Korelasi Negatif cukup kuat.
IKG (\(X2\)) vs Urbanisasi (\(X5\)) = Korelasi Negatif kuat.
Cerai (\(X6\)) vs Pendidikan (\(X3\)) = Korelasi Negatif cukup kuat.
Ekonomi (\(X4\)) vs Urbanisasi (\(X5\)) = Korelasi Positif cukup kuat.
Pendidikan (\(X3\)) vs Urbanisasi (\(X5\)) = Korelasi Positif cukup kuat.
Cerai (\(X6\)) vs IKG (\(X2\)) = Korelasi Positif cukup kuat.
# 1. Tentukan variabel mana saja yang mau dicek
# (X5 kita skip karena isinya 0, tidak bisa di-log)
vars_x = c("X1", "X2", "X3", "X4", "X5", "X6")
par(mfrow = c(1, 2)) # 1 Baris, 2 Kolom
for (var in vars_x) {
nilai = data[[var]]
# Histogram
hist(nilai,
main = paste("Distribusi ASLI:", var),
xlab = var,
col = "lightblue",
border = "white",
breaks = 15) # Breaks biar batangnya lebih detail
# Periksa pakah semua angkanya > 0? (Log butuh angka positif)
if (all(nilai > 0, na.rm = TRUE)) {
hist(log(nilai),
main = paste("Distribusi LOG:", var),
xlab = paste("Log(", var, ")"),
col = "lightgreen",
border = "white",
breaks = 15)
} else {
# Jika ada angka 0 atau negatif, tampilkan pesan error di gambar
plot.new()
text(0.5, 0.5, "Data mengandung 0/ Negatif\nTidak bisa di-Log", col = "red", cex = 1.2)
}
# Beri jeda sedikit biar kamu sempat lihat (opsional, bisa dihapus)
Sys.sleep(0.5)
}
# 4. Kembalikan tata letak gambar ke normal (1 gambar per layar)
par(mfrow = c(1, 1))
# Buat kolom baru "Wilayah"
data2 = data_clean %>%
mutate(Wilayah = ifelse(Provinsi %in% c("DKI Jakarta", "Jawa Barat", "Jawa Tengah",
"DI Yogyakarta", "Jawa Timur", "Banten"),
"Jawa", "Luar Jawa"))
# Plot Boxplot membandingkan Jawa vs Luar Jawa
ggplot(data2, aes(x = Wilayah, y = Y, fill = Wilayah)) +
geom_boxplot() +
labs(title = "Perbandingan Persentase Pemuda Kawin: Jawa vs Luar Jawa",
y = "Persentase Pemuda Kawin") +
theme_minimal()
## Warning: Removed 16 rows containing non-finite outside the scale range
## (`stat_boxplot()`).