Tugas Individu TPG
Mhd Haekal Hakim
Library
library(tidyverse)
library(lubridate)
library(janitor)
library(imputeTS)
library(zoo)
library(TSclust)
library(cluster)
library(factoextra)
library(purrr)
library(sf)
library(dtwclust)
library(proxy)
library(dplyr)
library(ggplot2)
library(stringr)Import Data
df <- read.csv("C:/Users/mhdha/Downloads/0_df_long.csv")Cleaning Data
Memisahkan Kolom
df_new <- read_delim("C:/Users/mhdha/Downloads/0_df_long.csv", delim = ";") %>% clean_names()
#Cek kolom awal
str(df_new)## spc_tbl_ [340,620 × 9] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ kab_kot : chr [1:340620] "Kab. Aceh Barat" "Kab. Aceh Barat" "Kab. Aceh Barat" "Kab. Aceh Barat" ...
## $ tahun : num [1:340620] 2022 2022 2022 2022 2022 ...
## $ bulan : chr [1:340620] "Januari" "Januari" "Januari" "Januari" ...
## $ produk : chr [1:340620] "Beras Premium" "Beras Medium" "Bawang Merah" "Bawang Putih (Bonggol)" ...
## $ harga : num [1:340620] 11429 9979 31000 28636 19409 ...
## $ kategori : chr [1:340620] "Beras" "Beras" "Bawang" "Bawang" ...
## $ kode_bps : num [1:340620] 11.1 11.1 11.1 11.1 11.1 ...
## $ kode_prov: num [1:340620] 11 11 11 11 11 11 11 11 11 11 ...
## $ nama_prov: chr [1:340620] "Aceh" "Aceh" "Aceh" "Aceh" ...
## - attr(*, "spec")=
## .. cols(
## .. KabKot = col_character(),
## .. Tahun = col_double(),
## .. Bulan = col_character(),
## .. Produk = col_character(),
## .. Harga = col_double(),
## .. Kategori = col_character(),
## .. KodeBPS = col_double(),
## .. KodeProv = col_double(),
## .. NamaProv = col_character()
## .. )
## - attr(*, "problems")=<externalptr>head(df_new)## # A tibble: 6 × 9
## kab_kot tahun bulan produk harga kategori kode_bps kode_prov nama_prov
## <chr> <dbl> <chr> <chr> <dbl> <chr> <dbl> <dbl> <chr>
## 1 Kab. Aceh Barat 2022 Janu… Beras… 11429 Beras 11.0 11 Aceh
## 2 Kab. Aceh Barat 2022 Janu… Beras… 9979 Beras 11.0 11 Aceh
## 3 Kab. Aceh Barat 2022 Janu… Bawan… 31000 Bawang 11.0 11 Aceh
## 4 Kab. Aceh Barat 2022 Janu… Bawan… 28636 Bawang 11.0 11 Aceh
## 5 Kab. Aceh Barat 2022 Janu… Cabai… 19409 Cabai 11.0 11 Aceh
## 6 Kab. Aceh Barat 2022 Janu… Cabai… 45706 Cabai 11.0 11 Acehread_delim berguna untuk memisahkan data yang terjadi seperti data kita, yang dimana pemisah datanya itu berupa tanda “;”, sehingga kita jadi lebih mudah membaca datanya
Mengambil Provinsi Sumatera Barat
df_sb <- df_new %>% filter(nama_prov == "Sumatera Barat")
head(df_sb)## # A tibble: 6 × 9
## kab_kot tahun bulan produk harga kategori kode_bps kode_prov nama_prov
## <chr> <dbl> <chr> <chr> <dbl> <chr> <dbl> <dbl> <chr>
## 1 Kab. Agam 2022 Januari Beras Pre… 13000 Beras 13.1 13 Sumatera…
## 2 Kab. Agam 2022 Januari Beras Med… 12000 Beras 13.1 13 Sumatera…
## 3 Kab. Agam 2022 Januari Bawang Me… 24833 Bawang 13.1 13 Sumatera…
## 4 Kab. Agam 2022 Januari Bawang Pu… 24833 Bawang 13.1 13 Sumatera…
## 5 Kab. Agam 2022 Januari Cabai Mer… 30167 Cabai 13.1 13 Sumatera…
## 6 Kab. Agam 2022 Januari Cabai Raw… NA Cabai 13.1 13 Sumatera…Bertujuan untuk memfilter data, yang di mana kita ingin melihat data dari Provinsi Sumatera Barat
Bentuk Data Baru
bulan_map <- c(
"Januari" = 1, "Februari" = 2, "Maret" = 3, "April" = 4,
"Mei" = 5, "Juni" = 6, "Juli" = 7, "Agustus" = 8,
"September" = 9, "Oktober" = 10, "November" = 11, "Desember" = 12
)
df_sb_new <- df_sb %>%
mutate(
tahun = as.numeric(tahun),
bulan = unname(bulan_map[bulan]),
tanggal = make_date(tahun, bulan, 1)
) %>%
select(
kab_kot = kab_kot,
tanggal,
produk,
harga,
kategori,
kode_bps
)
head(df_sb_new)## # A tibble: 6 × 6
## kab_kot tanggal produk harga kategori kode_bps
## <chr> <date> <chr> <dbl> <chr> <dbl>
## 1 Kab. Agam 2022-01-01 Beras Premium 13000 Beras 13.1
## 2 Kab. Agam 2022-01-01 Beras Medium 12000 Beras 13.1
## 3 Kab. Agam 2022-01-01 Bawang Merah 24833 Bawang 13.1
## 4 Kab. Agam 2022-01-01 Bawang Putih (Bonggol) 24833 Bawang 13.1
## 5 Kab. Agam 2022-01-01 Cabai Merah Keriting 30167 Cabai 13.1
## 6 Kab. Agam 2022-01-01 Cabai Rawit Merah NA Cabai 13.1Mengubah bentuk data, yang awalnya terdiri dari kolom tahun dan bulan, kita gabung menjadi kolom baru yang bernama “tanggal” yang terdiri dari kolom tahun dan bulan, serta ditambahkan tanggal agar kita bisa melihat pola waktu nya
Mengisi Data yang Kosong
df_impute <- df_sb_new %>%
group_by(kab_kot, kategori) %>%
arrange(tanggal) %>%
mutate(
harga_filled = na.locf(harga, na.rm = FALSE),
harga_filled = na.locf(harga_filled, fromLast = TRUE),
harga = harga_filled
) %>%
select(-harga_filled) %>%
ungroup()
head(df_impute)## # A tibble: 6 × 6
## kab_kot tanggal produk harga kategori kode_bps
## <chr> <date> <chr> <dbl> <chr> <dbl>
## 1 Kab. Agam 2022-01-01 Beras Premium 13000 Beras 13.1
## 2 Kab. Agam 2022-01-01 Beras Medium 12000 Beras 13.1
## 3 Kab. Agam 2022-01-01 Bawang Merah 24833 Bawang 13.1
## 4 Kab. Agam 2022-01-01 Bawang Putih (Bonggol) 24833 Bawang 13.1
## 5 Kab. Agam 2022-01-01 Cabai Merah Keriting 30167 Cabai 13.1
## 6 Kab. Agam 2022-01-01 Cabai Rawit Merah 30167 Cabai 13.1Bertujuan untuk mengisi data yang kosong menggunakan tetangga terdekatnya berdasarkan dari kategori dan dari kab/kotanya
Eksplorasi Data
Boxplot
ggplot(df_impute, aes(kategori, harga)) +
geom_boxplot() +
labs(title = "Distribusi Harga Komoditas Pangan di Sumatera Barat") +
theme_minimal()
Menunjukkan bagaimana sebaran harga berdasarkan kategori yang ada
Plot Data Time Series
df_month <- df_impute %>%
filter(kab_kot == "Kota Pariaman", kategori == "Beras") %>%
mutate(bulan = floor_date(tanggal, "month")) %>%
group_by(bulan) %>%
summarise(harga = mean(harga))
head(df_month)## # A tibble: 6 × 2
## bulan harga
## <date> <dbl>
## 1 2022-01-01 14000
## 2 2022-02-01 14000
## 3 2022-03-01 14000
## 4 2022-04-01 14000
## 5 2022-05-01 14010.
## 6 2022-06-01 13488ggplot(df_month, aes(bulan, harga)) +
geom_line(color = "blue") +
theme_minimal()
Melihat bagaimana harga dari komoditas beras dari Januari 2022 sampai Oktober 2025
Agregasi Data
df_wide <- df_impute %>%
group_by(kab_kot, tanggal, kategori) %>%
summarise(harga = mean(harga, na.rm = TRUE), .groups = 'drop') %>%
arrange(kab_kot, tanggal) %>%
pivot_wider(names_from = kategori, values_from = harga)
head(df_wide)## # A tibble: 6 × 9
## kab_kot tanggal Bawang Beras Cabai Gula Jagung MinyakGoreng Protein
## <chr> <date> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Kab. Agam 2022-01-01 24833 12500 30167 14000 5400 20000 63572.
## 2 Kab. Agam 2022-02-01 28500 13037 38705 14741 5400 17008 58444.
## 3 Kab. Agam 2022-03-01 28403 13250 36813 14484 5400 19206 63721.
## 4 Kab. Agam 2022-04-01 33224. 13008. 32138 15000 5400 17733. 64813.
## 5 Kab. Agam 2022-05-01 35278 13000 38530. 14815 5400 17852. 68012.
## 6 Kab. Agam 2022-06-01 38429 12654. 76361. 14938 5400 17840. 64939.# Pastikan semua kolom numeric kecuali kab & tanggal
df_wide <- df_wide %>%
mutate(across(-c(kab_kot, tanggal), as.numeric))
# Split menjadi list time-series per kab/kot
list_ts <- df_wide %>%
split(.$kab_kot) %>%
map(~ .x %>%
select(-kab_kot, -tanggal) %>%
as.matrix()
)dist_matrix <- proxy::dist(list_ts, method = "dtw_basic", window.size = 5)
dist_mat <- as.matrix(dist_matrix)Bertujuan untuk mengubah data menjadi format time-series wide, dan memecahnya menjadi data time series sper kabupaten/kota. Selain itu juga, berguna untuk mengukur kemiripan pola harga antar kab/kota sehingga menghasilkan distance matriks Dynamic Time Warping (DTW)
Penentuan Cluster yang Terbaik
# Silhouette
sil_scores <- sapply(2:8, function(k){
cl <- cutree(hclust(dist_matrix, method = "ward.D2"), k)
ss <- silhouette(cl, dist_matrix)
summary(ss)$avg.width
})
best_k <- which.max(sil_scores) + 1
best_k## [1] 2# Elbow WSS
wss <- sapply(2:8, function(k){
cl <- cutree(hclust(dist_matrix, method = "ward.D2"), k)
sum(sapply(unique(cl), function(ci){
members <- which(cl == ci)
if(length(members) > 1) sum(dist_mat[members, members]) else 0
}))
})
plot(2:8, wss, type = "b", pch = 19, xlab = "Jumlah Cluster", ylab = "Within-Cluster Distance",
main = "Elbow Method untuk K Optimal")
sil_scores## [1] 0.39969945 0.11536792 0.10946844 0.10788082 0.10360626 0.07359473 0.07479442best_k <- which.max(sil_scores) + 1
best_k## [1] 2Menentukan jumlah cluster terbaik digunakan menggunakan Silhouette dan Elbow
Dendrogram Clustering
hc <- hclust(dist_matrix, method = "ward.D2")
cluster_assign <- tibble(
kab_kot = names(list_ts),
cluster = cutree(hc, k = best_k)
)plot(hc, labels = names(list_ts), main = "Dendrogram Clustering Harga Pangan (DTW)")
rect.hclust(hc, k = best_k, border = 2:3)
Melihat bagaimana pola clustering yang terbentuk berdasarkan kab/kota
Pemetaan
Membersihkan Data di Kode BPS
cluster_assign <- df_impute %>%
distinct(kab_kot, kode_bps) %>%
left_join(cluster_assign, by = "kab_kot") %>%
mutate(
kode_bps = as.character(kode_bps),
# Mapping manual jika perlu
kode_bps = case_when(
kode_bps == "13.1" ~ "1310",
TRUE ~ str_replace(kode_bps, "\\.", "") %>% str_pad(width = 4, pad = "0")
)
)
cluster_assign## # A tibble: 19 × 3
## kab_kot kode_bps cluster
## <chr> <chr> <int>
## 1 Kab. Agam 1306 1
## 2 Kab. Dharmasraya 1310 1
## 3 Kab. Kepulauan Mentawai 1309 2
## 4 Kab. Lima Puluh Kota 1307 1
## 5 Kab. Padang Pariaman 1305 1
## 6 Kab. Pasaman 1308 1
## 7 Kab. Pasaman Barat 1312 1
## 8 Kab. Pesisir Selatan 1301 1
## 9 Kab. Sijunjung 1303 1
## 10 Kab. Solok 1302 1
## 11 Kab. Solok Selatan 1311 1
## 12 Kab. Tanah Datar 1304 1
## 13 Kota Bukittinggi 1375 1
## 14 Kota Padang 1371 1
## 15 Kota Padang Panjang 1374 1
## 16 Kota Pariaman 1377 1
## 17 Kota Payakumbuh 1376 1
## 18 Kota Sawah Lunto 1373 1
## 19 Kota Solok 1372 1Membersihakn beberapa data yang miss terlewat
Membaca File JSON
map_idn <- st_read("C:/Users/mhdha/Downloads/gadm41_IDN_2.json/gadm41_IDN_2.json")## Reading layer `gadm41_IDN_2' from data source
## `C:\Users\mhdha\Downloads\gadm41_IDN_2.json\gadm41_IDN_2.json'
## using driver `GeoJSON'
## Simple feature collection with 502 features and 13 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 95.0102 ymin: -11.0075 xmax: 141.0194 ymax: 6.0766
## Geodetic CRS: WGS 84colnames(map_idn)## [1] "GID_2" "GID_0" "COUNTRY" "GID_1" "NAME_1" "NL_NAME_1"
## [7] "NAME_2" "VARNAME_2" "NL_NAME_2" "TYPE_2" "ENGTYPE_2" "CC_2"
## [13] "HASC_2" "geometry"map_sumbar <- map_idn %>%
filter(NAME_1 == "SumateraBarat") %>%
select(KabKot = NAME_2, geometry)
view(map_sumbar)Mengubah Nama Kab/Kota
cluster_assign <- cluster_assign %>%
mutate(
kab_kot = case_when(
kab_kot == "Kab. Agam" ~ "Agam",
kab_kot == "Kab. Dharmasraya" ~ "Dharmasraya",
kab_kot == "Kab. Kepulauan Mentawai" ~ "KepulauanMentawai",
kab_kot == "Kab. Lima Puluh Kota" ~ "LimaPuluhKota",
kab_kot == "Kab. Padang Pariaman" ~ "PadangPariaman",
kab_kot == "Kab. Pasaman" ~ "Pasaman",
kab_kot == "Kab. Pasaman Barat" ~ "PasamanBarat",
kab_kot == "Kab. Pesisir Selatan" ~ "PesisirSelatan",
kab_kot == "Kab. Sijunjung" ~ "Sijunjung",
kab_kot == "Kab. Solok" ~ "Solok",
kab_kot == "Kab. Solok Selatan" ~ "SolokSelatan",
kab_kot == "Kab. Tanah Datar" ~ "TanahDatar",
kab_kot == "Kota Bukittinggi" ~ "Bukittinggi",
kab_kot == "Kota Padang" ~ "Padang",
kab_kot == "Kota Padang Panjang" ~ "PadangPanjang",
kab_kot == "Kota Pariaman" ~ "Pariaman",
kab_kot == "Kota Payakumbuh" ~ "Payakumbuh",
kab_kot == "Kota Sawah Lunto" ~ "Sawahlunto",
kab_kot == "Kota Solok" ~ "KotaSolok",
)
)
cluster_assign## # A tibble: 19 × 3
## kab_kot kode_bps cluster
## <chr> <chr> <int>
## 1 Agam 1306 1
## 2 Dharmasraya 1310 1
## 3 KepulauanMentawai 1309 2
## 4 LimaPuluhKota 1307 1
## 5 PadangPariaman 1305 1
## 6 Pasaman 1308 1
## 7 PasamanBarat 1312 1
## 8 PesisirSelatan 1301 1
## 9 Sijunjung 1303 1
## 10 Solok 1302 1
## 11 SolokSelatan 1311 1
## 12 TanahDatar 1304 1
## 13 Bukittinggi 1375 1
## 14 Padang 1371 1
## 15 PadangPanjang 1374 1
## 16 Pariaman 1377 1
## 17 Payakumbuh 1376 1
## 18 Sawahlunto 1373 1
## 19 KotaSolok 1372 1Menyamakan nama kab/kota di data dengan file map json
Pengecekan Kesesuaian
# Cek intersect lagi
intersect(map_sumbar$KabKot, cluster_assign$kab_kot)## [1] "Agam" "Bukittinggi" "Dharmasraya"
## [4] "KepulauanMentawai" "KotaSolok" "LimaPuluhKota"
## [7] "Padang" "PadangPanjang" "PadangPariaman"
## [10] "Pariaman" "Pasaman" "PasamanBarat"
## [13] "Payakumbuh" "PesisirSelatan" "Sawahlunto"
## [16] "Sijunjung" "Solok" "SolokSelatan"
## [19] "TanahDatar"setdiff(map_sumbar$KabKot, cluster_assign$kab_kot)## [1] "Danau"setdiff(cluster_assign$kab_kot, map_sumbar$KabKot)## character(0)Pengecekan apakah kedua data sudah memiliki kecocokan yang sama atau tidak
Penggabungan Data
# Join cluster ke peta
map_cluster <- map_sumbar %>%
left_join(cluster_assign %>% select(kab_kot, cluster),
by = c("KabKot" = "kab_kot"))Menggabungkan data dengan json sehingga nanti hasil akhirnya bisa ditampilkan dalam bentuk peta
Cloropeth MAP
# Plot peta
ggplot(map_cluster) +
geom_sf(aes(fill = factor(cluster)), color = "black") +
scale_fill_brewer(palette = "Set2", name = "Cluster") +
theme_minimal() +
labs(
title = "Cluster Kabupaten/Kota Sumatera Barat",
subtitle = "Hasil clustering berdasarkan data",
caption = "Sumber: df_impute & GeoJSON GADM"
) +
theme(
axis.text = element_blank(),
axis.ticks = element_blank()
)
Hasil clustering dari yang telah dilakukan
Hijau - Pola harga komoditas realtif mirip antar daerah daratan Sumatera Barat, berfluktuasi konsisten, memiliki musiman yang sama, dan mengikuti tren pasar utama. Kendali distribusi pangan antarwilyaha cenderung terhubung, sehingga pola harganya mirip
Oranye - Pola harga lebih fluktuatif dibanding daerah daratan, harga komoditas memungkinkan lebih mahal, berubah dengan cepat, dan dipengaruhi dari cuaca laut serta transportasi laut. Selain itu juga rantai pasoknya lebih terpencil dibandingkan dengan daerah daratan.
NA - Danau Singkarak