Pertemuan 5

1 Persiapan Data

1.1 Package

library(tidyverse)
library(readxl)
library(janitor)
library(scales)
library(ggcorrplot)
library(knitr)

1.2 Import dan Pembersihan Ringkas

Langkah ini mengulang keputusan pada pertemuan 4: nama kolom dirapikan, dan 23 baris dengan visibility kosong dibuang (0.36% data, hampir seluruhnya di Papua Tengah).

path_data <- "C:/Users/user/OneDrive/backup/OneDrive/Desktop/5th/PSD/Kelompok 3 - Pekan 4/Dataset Kelompok 3 PSD (Data Scrapping).xlsx"

df_raw <- read_excel(path_data, sheet = 1)

df <- df_raw %>%
  mutate(tanggal = as.Date(tanggal)) %>%
  rename(
    suhu_c            = `suhu_C`,
    suhu_min_c        = `suhu_min_C`,
    suhu_max_c        = `suhu_max_C`,
    kelembapan_pct    = `kelembapan_%`,
    tekanan_hpa       = `tekanan_hPa`,
    kecepatan_angin   = `kecepatan_angin_m_s`,
    arah_angin        = `arah_angin_derajat`,
    tutupan_awan_pct  = `tutupan_awan_%`,
    visibility        = `visibility_m`
  ) %>%
  clean_names() %>%
  filter(!is.na(visibility))        

dim(df)

## [1] 6330   20

1.3 Membuat Peubah Turunan

Ada dua peubah kategorik baru yang dibuat dari data yang sudah ada:

• kategori: kategori_cuaca diberi label bahasa Indonesia (Cerah, Berawan, Hujan).
• pulau: 33 provinsi diringkas menjadi 7 kelompok pulau/wilayah. Provinsinya masih terlalu banyak (dan ada yang hanya berisi sedikit kecamatan) sehingga sulit dibaca dan tidak memenuhi syarat frekuensi harapan pada uji Chi-square.

df <- df %>%
  mutate(
    kategori = factor(kategori_cuaca,
                      levels = c("Clear", "Clouds", "Rain"),
                      labels = c("Cerah", "Berawan", "Hujan")),
    pulau = case_when(
      provinsi %in% c("Aceh", "Sumatera Utara", "Sumatera Barat", "Riau", "Kepulauan Riau",
                      "Jambi", "Sumatera Selatan", "Kepulauan Bangka Belitung",
                      "Bengkulu", "Lampung")                                   ~ "Sumatera",
      provinsi %in% c("DKI Jakarta", "Banten", "Jawa Barat", "Jawa Tengah",
                      "DI Yogyakarta", "Jawa Timur")                           ~ "Jawa",
      provinsi %in% c("Bali", "Nusa Tenggara Barat", "Nusa Tenggara Timur")    ~ "Bali & Nusa Tenggara",
      provinsi %in% c("Kalimantan Barat", "Kalimantan Tengah", "Kalimantan Selatan",
                      "Kalimantan Timur", "Kalimantan Utara")                  ~ "Kalimantan",
      provinsi %in% c("Sulawesi Barat", "Sulawesi Selatan", "Sulawesi Tengah",
                      "Sulawesi Tenggara", "Gorontalo")                        ~ "Sulawesi",
      provinsi %in% c("Maluku", "Maluku Utara")                                ~ "Maluku",
      provinsi %in% c("Papua Barat", "Papua Tengah")                           ~ "Papua"
    ),
    pulau = factor(pulau)
  )

sum(is.na(df$pulau))

## [1] 0

df %>% count(pulau, name = "jumlah_kecamatan") %>% arrange(desc(jumlah_kecamatan))

1.4 Tema dan Warna

warna_kategori <- c("Cerah" = "#F2A900", "Berawan" = "#7F8FA0", "Hujan" = "#1F6FB2")

warna_pulau <- c("Sumatera" = "#E4572E", "Jawa" = "#F3A712", "Bali & Nusa Tenggara" = "#8DB580",
                 "Kalimantan" = "#2E933C", "Sulawesi" = "#4F9DC9", "Maluku" = "#7B4B94",
                 "Papua" = "#264653")

theme_set(
  theme_minimal(base_size = 12) +
    theme(
      plot.title       = element_text(face = "bold", size = 14),
      plot.subtitle    = element_text(color = "grey35"),
      plot.caption     = element_text(color = "grey45", size = 9, hjust = 0),
      panel.grid.minor = element_blank(),
      legend.position  = "bottom"
    )
)

sumber <- paste0("Sumber: data scraping cuaca kecamatan, 29 Agustus 2026 pukul 12:00 (n = ",
                 format(nrow(df), big.mark = "."), " kecamatan)")

1.5 Peubah yang Dipakai dan yang Sengaja Tidak Dipakai

Kolom	Keputusan	Alasan
suhu_min_c, suhu_max_c	Tidak dipakai	Nyaris sama kayaksuhu_c (r ≈ 0.99–1.00), jadi cuma menambah redundansi
arah_angin	Tidak dipakai di korelasi	Bentuknya derajat kompas : 1° dan 359° sebenarnya berdekatan, tetapi korelasi Pearson/Spearman ngeanggapnya berjauhan
cuaca_numerik, kategori_cuaca_numerik	Tidak dipakai	Cuma kode label dari kategori (misalnya Clouds = 1, Clear = 2, Rain = 3), bukan angka yang bermakna. Menghitung korelasi pada kode ini menyesatkan
id_kecamatan, tanggal, waktu	Tidak dipakai	Pengenal / konstan (satu tanggal dan satu jam)
latitude, longitude	Dipakai	Posisi geografis, bisa menjelaskan pola antar wilayah

num_vars <- c("suhu_c", "kelembapan_pct", "tekanan_hpa", "kecepatan_angin",
              "tutupan_awan_pct", "visibility", "suhu_min_c", "suhu_max_c")

label_num <- c(suhu_c = "Suhu", kelembapan_pct = "Kelembapan", tekanan_hpa = "Tekanan",
               kecepatan_angin = "Kec. angin", tutupan_awan_pct = "Tutupan awan",
               visibility = "Visibility", suhu_min_c= "Suhu Min", suhu_max_c = " Suhu Max")
kor_spearman <- cor(df[num_vars], method = "spearman", use = "complete.obs")
colnames(kor_spearman) <- label_num[colnames(kor_spearman)]
rownames(kor_spearman) <- label_num[rownames(kor_spearman)]

2 Q1: Pasangan Numerik yang Paling Berhubungan

ggcorrplot(kor_spearman, 
           type = "lower", 
           show.diag = TRUE,    
           lab = TRUE, 
           lab_size = 3.6,
           hc.order = FALSE,     
           colors = c("#1F6FB2", "white", "#C44E52"),
           outline.color = "white", 
           ggtheme = theme_minimal(base_size = 12),
           legend.title = "Spearman",
           title = "Matriks Korelasi Spearman antar Peubah Numerik") +
  labs(caption = sumber) +
  theme(
    panel.grid = element_blank(),
    axis.text.x = element_text(angle = 45, hjust = 1)
  )

tabel_kor <- combn(num_vars, 2, simplify = FALSE) %>%
  map_dfr(function(p) {
    tibble(
      peubah_1 = label_num[[p[1]]],
      peubah_2 = label_num[[p[2]]],
      pearson  = cor(df[[p[1]]], df[[p[2]]], method = "pearson", use = "complete.obs"),
      spearman = cor(df[[p[1]]], df[[p[2]]], method = "spearman", use = "complete.obs")
    )
  }) %>%
  mutate(
    arah = if_else(spearman > 0, "positif", "negatif"),
    kekuatan = case_when(
      abs(spearman) >= 0.80 ~ "sangat kuat",
      abs(spearman) >= 0.60 ~ "kuat",
      abs(spearman) >= 0.40 ~ "sedang",
      abs(spearman) >= 0.20 ~ "lemah",
      TRUE                  ~ "sangat lemah"
    )
  ) %>%
  arrange(desc(abs(spearman)))

tabel_kor %>%
  slice_head(n = 10) %>%
  kable(
    digits = 2, 
    caption = "10 Pasangan Peubah dengan Korelasi Terkuat (Urut |Spearman|)",
    col.names = c("Peubah 1", "Peubah 2", "Pearson", "Spearman", "Arah", "Kekuatan")
  )

10 Pasangan Peubah dengan Korelasi Terkuat (Urut |Spearman|)

Peubah 1	Peubah 2	Pearson	Spearman	Arah	Kekuatan
Suhu	Suhu Max	1.00	1.00	positif	sangat kuat
Suhu	Suhu Min	1.00	1.00	positif	sangat kuat
Suhu Min	Suhu Max	0.99	0.99	positif	sangat kuat
Kelembapan	Suhu Min	-0.70	-0.79	negatif	kuat
Suhu	Kelembapan	-0.70	-0.78	negatif	kuat
Kelembapan	Suhu Max	-0.70	-0.78	negatif	kuat
Tekanan	Suhu Max	-0.59	-0.56	negatif	sedang
Kec. angin	Suhu Max	0.46	0.56	positif	sedang
Suhu	Kec. angin	0.46	0.55	positif	sedang
Suhu	Tekanan	-0.58	-0.55	negatif	sedang

Interpretasi.

• Hubungan terkuat adalah suhu dan kelembapan (Spearman ≈ −0.78, kuat mendekati sangat kuat): kecamatan yang lebih panas cenderung lebih kering. Suhu juga cenderung lebih tinggi di tempat dengan tekanan lebih rendah (≈ −0.55) dan angin lebih kencang (≈ +0.55), keduanya berhubungan sedang.
• Selisih antara Pearson dan Spearman terlihat pada pasangan suhu–kelembapan (−0.70 vs −0.78), tanda bahwa hubungannya menurun secara konsisten tetapi tidak berbentuk garis lurus
• Lintang berkorelasi sedang dengan kelembapan (≈ +0.49) dan kecepatan angin (≈ −0.47).
• Peubah visibility hanya berkorelasi lemah dengan peubah lain (|ρ| ≤ 0.20), wajar karena sekitar 97% nilainya sama (10.000 m) sehingga nyaris tidak bervariasi.
• Korelasi hanya menunjukkan kecenderungan bergerak bersama, bukan bahwa satu peubah menyebabkan peubah lain.

3 Q2: Suhu vs Kelembapan (Apakah Konsisten di Setiap Kategori Cuaca?)

3.1 Scatter plot keseluruhan

ggplot(df, aes(x = suhu_c, y = kelembapan_pct)) +
  geom_point(aes(color = kategori), alpha = 0.45, size = 1.5) +
  geom_smooth(method = "loess", se = FALSE, color = "black", linewidth = 0.9) +
  scale_color_manual(values = warna_kategori, name = "Kategori cuaca") +
  guides(color = guide_legend(override.aes = list(alpha = 1, size = 3))) +
  labs(title = "Suhu vs Kelembapan",
       subtitle = "Setiap titik = satu kecamatan; garis hitam = tren loess",
       x = "Suhu (\u00b0C)", y = "Kelembapan (%)", caption = sumber)

uji_p <- cor.test(df$suhu_c, df$kelembapan_pct, method = "pearson")
uji_s <- cor.test(df$suhu_c, df$kelembapan_pct, method = "spearman", exact = FALSE)

tibble(
  Metode      = c("Pearson (r)", "Spearman (rho)"),
  Koefisien   = c(unname(uji_p$estimate), unname(uji_s$estimate)),
  `CI 95%`    = c(sprintf("[%.2f ; %.2f]", uji_p$conf.int[1], uji_p$conf.int[2]), "-"),
  `p-value`   = c(format.pval(uji_p$p.value, eps = .Machine$double.eps),
                  format.pval(uji_s$p.value, eps = .Machine$double.eps))
) %>% kable(digits = 3, caption = "Korelasi suhu vs kelembapan (seluruh kecamatan)")

Korelasi suhu vs kelembapan (seluruh kecamatan)

Metode	Koefisien	CI 95%	p-value
Pearson (r)	-0.704	[-0.72 ; -0.69]	< 2.22e-16
Spearman (rho)	-0.784	-	< 2.22e-16

3.2 Pola dalam tiap kategori cuaca

kor_dalam <- df %>%
  group_by(kategori) %>%
  summarise(n = n(),
            pearson  = cor(suhu_c, kelembapan_pct),
            spearman = cor(suhu_c, kelembapan_pct, method = "spearman"),
            .groups = "drop")

ggplot(df, aes(x = suhu_c, y = kelembapan_pct, color = kategori)) +
  geom_point(alpha = 0.4, size = 1.3, show.legend = FALSE) +
  geom_smooth(method = "loess", se = FALSE, color = "black", linewidth = 0.8) +
  geom_text(data = kor_dalam,
            aes(x = Inf, y = Inf,
                label = paste0("\u03c1 = ", round(spearman, 2), "\nn = ", n)),
            hjust = 1.1, vjust = 1.3, color = "black", size = 3.8, inherit.aes = FALSE) +
  facet_wrap(~ kategori) +
  scale_color_manual(values = warna_kategori) +
  labs(title = "Hubungan Suhu dengan Kelembapan per Kategori Cuaca",
       x = "Suhu (\u00b0C)", y = "Kelembapan (%)", caption = sumber)

kor_dalam %>% kable(digits = 2, col.names = c("Kategori", "n", "Pearson (r)", "Spearman (rho)"),
                    caption = "Korelasi suhu vs kelembapan di dalam tiap kategori")

Korelasi suhu vs kelembapan di dalam tiap kategori

Kategori	n	Pearson (r)	Spearman (rho)
Cerah	1716	-0.72	-0.78
Berawan	3958	-0.78	-0.82
Hujan	656	0.00	-0.18

Untuk memahami mengapa korelasi di kategori Hujan berbeda, berikut asal provinsi kecamatan Hujan yang bersuhu dingin (< 20 °C):

df %>%
  filter(kategori == "Hujan", suhu_c < 20) %>%
  count(provinsi, sort = TRUE) %>%
  mutate(persen = round(100 * n / sum(n), 1)) %>%
  slice_head(n = 5) %>%
  kable(col.names = c("Provinsi", "Jumlah kecamatan", "% dari kecamatan Hujan < 20 \u00b0C"),
        caption = "Asal kecamatan Hujan bersuhu < 20 \u00b0C (5 provinsi teratas)")

Asal kecamatan Hujan bersuhu < 20 °C (5 provinsi teratas)

Provinsi	Jumlah kecamatan	% dari kecamatan Hujan < 20 °C
Papua Tengah	84	71.2
Sumatera Utara	10	8.5
Aceh	6	5.1
Jambi	4	3.4
Maluku Utara	3	2.5

Interpretasi.

• Secara keseluruhan hubungan suhu–kelembapan negatif dan kuat (Pearson ≈ −0.70, Spearman ≈ −0.78) tetapi tidak berbentuk garis lurus: di bawah ±25 °C kelembapan tinggi dan relatif datar (sekitar 90%), sedangkan di atas ±25 °C kelembapan turun tajam sampai sekitar 45% pada suhu 33 °C. Karena itu Spearman lebih mewakili dibanding Pearson.
• Hubungan kuat di kategori Cerah dan Berawan (ρ ≈ −0.78 dan −0.82), tetapi hampir hilang di kategori Hujan (ρ ≈ −0.18; Pearson ≈ 0).
• Pada kategori Hujan, sebagian besar titik berkumpul di sekitar 23–27 °C dengan kelembapan tinggi, sementara titik bersuhu dingin justru tersebar pada kelembapan 50–85%. Tabel di atas menunjukkan titik dingin ini sebagian besar berasal dari Papua Tengah, wilayah yang di pertemuan 4 sudah kita identifikasi sebagai pegunungan. Jadi kategori Hujan tampaknya mencampur dua kelompok lokasi yang berbeda, dan itu yang meredam korelasinya. Ini dugaan berdasarkan pola data, belum diuji.
• Ini contoh bahwa satu angka korelasi global bisa menutupi perbedaan antar subkelompok.
• Kehati-hatian: ini bukan bukti bahwa suhu menentukan kelembapan. Kelembapan relatif memang secara fisik bergantung pada suhu, tetapi pada data ini ada faktor lain yang ikut berperan (ketinggian tempat, wilayah, jam lokal saat data direkam) dan tidak tercatat.

4 Q3 : Apakah Karakteristik Cuaca Berbeda Antar Kategori Cuaca?

df %>%
  select(kategori, suhu_c, kelembapan_pct, kecepatan_angin) %>%
  pivot_longer(-kategori, names_to = "peubah", values_to = "nilai") %>%
  mutate(peubah = factor(peubah,
                         levels = c("suhu_c", "kelembapan_pct", "kecepatan_angin"),
                         labels = c("Suhu (\u00b0C)", "Kelembapan (%)", "Kecepatan angin (m/s)"))) %>%
  ggplot(aes(x = kategori, y = nilai, fill = kategori)) +
  geom_violin(alpha = 0.35, color = NA, show.legend = FALSE) +
  geom_boxplot(width = 0.22, outlier.size = 0.6, outlier.alpha = 0.35, show.legend = FALSE) +
  facet_wrap(~ peubah, scales = "free_y", nrow = 1) +
  scale_fill_manual(values = warna_kategori) +
  labs(title = "Karakteristik Cuaca menurut Kategori Cuaca",
       x = NULL, y = NULL, caption = sumber)

df %>%
  group_by(kategori) %>%
  summarise(n = n(),
            suhu_median       = median(suhu_c),
            kelembapan_median = median(kelembapan_pct),
            angin_median      = median(kecepatan_angin),
            .groups = "drop") %>%
  kable(digits = 2, caption = "Median tiap kategori cuaca",
        col.names = c("Kategori", "n", "Suhu (C)", "Kelembapan (%)", "Kec. angin (m/s)"))

Median tiap kategori cuaca

Kategori	n	Suhu (C)	Kelembapan (%)	Kec. angin (m/s)
Cerah	1716	27.24	72	2.87
Berawan	3958	27.20	71	2.59
Hujan	656	24.54	90	1.24

Lanjut pakai Uji Kruskal-Wallis (bukan ANOVA) karena sebaran peubah tidak normal dan ragamnya tidak homogen. Ukuran efek yang dilaporkan adalah \(\eta^2_H = (H - k + 1)/(n - k)\) dengan patokan 0.01 = kecil, 0.06 = sedang, 0.14 = besar.

uji_kw <- function(data, peubah, kelompok) {
  x <- data[[peubah]]
  g <- droplevels(data[[kelompok]])
  kw <- kruskal.test(x ~ g)
  k <- nlevels(g); n <- length(x); H <- unname(kw$statistic)
  tibble(peubah = peubah, H = H, db = k - 1, p_value = kw$p.value,
         eta2_H = (H - k + 1) / (n - k))
}

hasil_kw_kategori <- c("suhu_c", "kelembapan_pct", "kecepatan_angin", "tekanan_hpa", "tutupan_awan_pct") %>%
  map_dfr(~ uji_kw(df, .x, "kategori")) %>%
  mutate(
    peubah  = recode(peubah, suhu_c = "Suhu", kelembapan_pct = "Kelembapan",
                     kecepatan_angin = "Kec. angin", tekanan_hpa = "Tekanan",
                     tutupan_awan_pct = "Tutupan awan (*)"),
    p_value = format.pval(p_value, digits = 3, eps = 1e-16),
    efek    = case_when(eta2_H >= 0.14 ~ "besar", eta2_H >= 0.06 ~ "sedang",
                        eta2_H >= 0.01 ~ "kecil", TRUE ~ "diabaikan")
  )

hasil_kw_kategori %>%
  kable(digits = 3, caption = "Kruskal-Wallis: peubah numerik ~ kategori cuaca",
        col.names = c("Peubah", "H", "db", "p-value", "eta^2_H", "Ukuran efek"))

Kruskal-Wallis: peubah numerik ~ kategori cuaca

Peubah	H	db	p-value	eta^2_H	Ukuran efek
Suhu	388.540	2	<1e-16	0.061	sedang
Kelembapan	331.741	2	<1e-16	0.052	kecil
Kec. angin	435.138	2	<1e-16	0.068	sedang
Tekanan	250.364	2	<1e-16	0.039	kecil
Tutupan awan (*)	3753.579	2	<1e-16	0.593	besar

(*) tutupan_awan_pct ikut ditampilkan hanya sebagai pembanding. Kategori cuaca (Cerah/Berawan/Hujan) memang ditentukan dari tutupan awan, sehingga hubungannya bersifat definisional dan bukan temuan.

Karena Kruskal-Wallis hanya menjawab “ada beda atau tidak”, uji lanjutan (Wilcoxon berpasangan dengan koreksi Holm) dipakai untuk melihat kategori mana yang berbeda pada peubah suhu:

pairwise.wilcox.test(df$suhu_c, df$kategori, p.adjust.method = "holm", exact = FALSE)

## 
##  Pairwise comparisons using Wilcoxon rank sum test with continuity correction 
## 
## data:  df$suhu_c and df$kategori 
## 
##         Cerah  Berawan
## Berawan 0.14   -      
## Hujan   <2e-16 <2e-16 
## 
## P value adjustment method: holm

Interpretasi.

• Kecamatan Hujan paling menonjol: suhu lebih rendah (median ≈ 24 °C dibanding ≈ 27 °C), kelembapan jauh lebih tinggi (median ≈ 90%), dan angin lebih tenang (median ≈ 1.2 m/s dibanding ≈ 2.6–2.9 m/s).
• Cerah dan Berawan hampir tidak berbeda untuk suhu (median 27.24 vs 27.20 °C; uji lanjutan p = 0.14, tidak nyata) dan kelembapan. Perbedaan di uji Kruskal-Wallis terutama disumbang oleh kategori Hujan.
• Semua p-value sangat kecil, tetapi itu karena jumlah data besar (n > 6.000). Yang lebih informatif adalah ukuran efek: η² suhu ≈ 0.06, kelembapan ≈ 0.05, dan angin ≈ 0.07, yaitu kecil sampai sedang. Jadi kategori cuaca hanya menjelaskan sebagian kecil variasi suhu, kelembapan, dan angin.
• Jangan tafsirkan sebagai “hujan menurunkan suhu”. Kecamatan hujan terkonsentrasi di wilayah tertentu (Papua dan Sumatera), dan titik hujan bersuhu dingin sebagian besar berasal dari Papua Tengah (ada di poin suhu vs kelembapan). Jadi suhu rendah bisa mencerminkan lokasi (misalnya ketinggian), bukan efek hujan.

5 Q4 : Apakah Proporsi Kategori Cuaca Berbeda Antar Pulau?

Karena kedua peubah kategorik, jadi pakai tabel kontingensi, stacked bar 100%

tab_pulau <- table(df$pulau, df$kategori)

addmargins(tab_pulau) %>%
  as.data.frame.matrix() %>%
  kable(caption = "Tabel kontingensi: jumlah kecamatan menurut pulau dan kategori cuaca")

Tabel kontingensi: jumlah kecamatan menurut pulau dan kategori cuaca

	Cerah	Berawan	Hujan	Sum
Bali & Nusa Tenggara	93	235	0	328
Jawa	673	1466	4	2143
Kalimantan	312	297	9	618
Maluku	35	188	9	232
Papua	9	276	209	494
Sulawesi	179	392	13	584
Sumatera	415	1104	412	1931
Sum	1716	3958	656	6330

round(100 * prop.table(tab_pulau, margin = 1), 1) %>%
  as.data.frame.matrix() %>%
  kable(caption = "Persentase per baris (% kategori cuaca di dalam tiap pulau)")

Persentase per baris (% kategori cuaca di dalam tiap pulau)

	Cerah	Berawan	Hujan
Bali & Nusa Tenggara	28.4	71.6	0.0
Jawa	31.4	68.4	0.2
Kalimantan	50.5	48.1	1.5
Maluku	15.1	81.0	3.9
Papua	1.8	55.9	42.3
Sulawesi	30.7	67.1	2.2
Sumatera	21.5	57.2	21.3

df_pulau <- df %>%
  count(pulau, kategori, .drop = FALSE) %>%
  group_by(pulau) %>%
  mutate(n_pulau = sum(n), pct = n / n_pulau) %>%
  ungroup()

urutan_pulau <- df_pulau %>%
  filter(kategori == "Hujan") %>%
  arrange(pct) %>%
  pull(pulau) %>%
  as.character()

df_pulau <- df_pulau %>%
  mutate(pulau = factor(pulau, levels = urutan_pulau),
         label_pulau = paste0(pulau, "  (n = ", format(n_pulau, big.mark = ".", trim = TRUE), ")"),
         label_pulau = factor(label_pulau, levels = unique(label_pulau[order(pulau)])))

ggplot(df_pulau, aes(x = label_pulau, y = pct, fill = kategori)) +
  geom_col(width = 0.72, color = "white", linewidth = 0.4,
           position = position_stack(reverse = TRUE)) +
  geom_text(aes(label = ifelse(pct >= 0.04, percent(pct, accuracy = 1), ""), color = kategori),
            position = position_stack(vjust = 0.5, reverse = TRUE),
            fontface = "bold", size = 3.6, show.legend = FALSE) +
  coord_flip() +
  scale_y_continuous(labels = percent, expand = expansion(mult = c(0, 0.03))) +
  scale_fill_manual(values = warna_kategori, name = "Kategori cuaca") +
  scale_color_manual(values = c("Cerah" = "grey15", "Berawan" = "white", "Hujan" = "white")) +
  labs(title = "Komposisi Kategori Cuaca menurut Pulau",
       subtitle = "Diurutkan dari pulau dengan proporsi hujan terkecil ke terbesar",
       x = NULL, y = "Proporsi kecamatan", caption = sumber)

uji_chi   <- chisq.test(tab_pulau)
cramers_v <- sqrt(unname(uji_chi$statistic) / (sum(tab_pulau) * (min(dim(tab_pulau)) - 1)))

tibble(
  `Chi-square`        = unname(uji_chi$statistic),
  db                  = unname(uji_chi$parameter),
  `p-value`           = format.pval(uji_chi$p.value, eps = .Machine$double.eps),
  `Cramer's V`        = cramers_v,
  `Frekuensi harapan minimum` = min(uji_chi$expected)
) %>% kable(digits = 3, caption = "Uji Chi-square kemandirian pulau dan kategori cuaca")

Uji Chi-square kemandirian pulau dan kategori cuaca

Chi-square	db	p-value	Cramer’s V	Frekuensi harapan minimum
1410.324	12	< 2.22e-16	0.334	24.043

round(uji_chi$stdres, 1) %>%
  as.data.frame.matrix() %>%
  kable(caption = "Residual terstandarisasi (nilai di luar -2 s.d. +2 = berbeda nyata dari harapan)")

Residual terstandarisasi (nilai di luar -2 s.d. +2 = berbeda nyata dari harapan)

	Cerah	Berawan	Hujan
Bali & Nusa Tenggara	0.5	3.5	-6.3
Jawa	5.5	6.9	-19.0
Kalimantan	13.8	-7.8	-7.6
Maluku	-4.2	5.9	-3.3
Papua	-13.2	-3.2	24.3
Sulawesi	2.0	2.4	-6.8
Sumatera	-6.7	-5.8	19.0

Patokan Cramér’s V untuk tabel dengan min(baris, kolom) - 1 = 2: 0.07 kecil, 0.21 sedang, 0.35 besar (Cohen).

Interpretasi.

• Proporsi kategori cuaca berbeda jelas antar pulau (Chi-square sangat signifikan, Cramér’s V ≈ 0.34, mendekati efek besar).
• Hujan pada snapshot ini terkonsentrasi di dua wilayah: Papua (sekitar 42% kecamatan berhujan) dan Sumatera (sekitar 21%). Pulau lain nyaris kering: Jawa dan Bali–Nusa Tenggara hampir 0%, Kalimantan, Sulawesi, dan Maluku kurang dari 4%. Residual terstandarisasi mengonfirmasi bahwa Hujan di Papua dan Sumatera jauh di atas harapan, dan di Jawa serta Bali–Nusa Tenggara jauh di bawah harapan.
• Kalimantan unik karena paling banyak Cerah (sekitar 50%), sedangkan Maluku paling didominasi Berawan (sekitar 81%).
• Syarat uji Chi-square terpenuhi: frekuensi harapan terkecil ≈ 24, jauh di atas batas minimum 5.
• Kehati-hatian: ini potret satu hari, satu jam, bukan iklim. Lokasi yang kebetulan berhujan pada 29 Agustus 2026 pukul 12:00 belum tentu mencerminkan pola hujan sepanjang tahun. Kecamatan bertetangga juga cenderung bercuaca sama, sehingga data tidak sepenuhnya independen. # Q5: Apakah Suhu Berbeda Antar Pulau?

set.seed(2026)

df %>%
  mutate(pulau = fct_reorder(pulau, suhu_c, .fun = median)) %>%
  ggplot(aes(x = pulau, y = suhu_c, fill = pulau)) +
  geom_jitter(width = 0.18, alpha = 0.12, size = 0.8, color = "grey30", show.legend = FALSE) +
  geom_boxplot(alpha = 0.75, outlier.shape = NA, show.legend = FALSE) +
  coord_flip() +
  scale_fill_manual(values = warna_pulau) +
  labs(title = "Sebaran Suhu menurut Pulau",
       subtitle = "Diurutkan menurut median suhu (tertinggi di atas)",
       x = NULL, y = "Suhu (\u00b0C)", caption = sumber)

df %>%
  group_by(pulau) %>%
  summarise(n = n(),
            median = median(suhu_c),
            Q1 = quantile(suhu_c, 0.25),
            Q3 = quantile(suhu_c, 0.75),
            min = min(suhu_c), max = max(suhu_c),
            .groups = "drop") %>%
  arrange(desc(median)) %>%
  kable(digits = 1, caption = "Ringkasan suhu (\u00b0C) per pulau, diurutkan dari median tertinggi")

Ringkasan suhu (°C) per pulau, diurutkan dari median tertinggi

pulau	n	median	Q1	Q3	min	max
Jawa	2143	29.4	27.6	30.6	14.1	36.4
Sulawesi	584	29.2	27.5	30.6	18.2	33.3
Papua	494	26.9	21.2	29.7	5.6	34.5
Sumatera	1931	25.7	24.0	27.0	12.1	29.8
Maluku	232	23.6	22.4	24.8	14.3	28.1
Kalimantan	618	22.9	21.7	24.8	14.6	29.2
Bali & Nusa Tenggara	328	22.0	19.4	25.4	11.7	28.0

uji_kw(df, "suhu_c", "pulau") %>%
  mutate(p_value = format.pval(p_value, digits = 3, eps = 1e-16)) %>%
  kable(digits = 3, caption = "Kruskal-Wallis: suhu ~ pulau",
        col.names = c("Peubah", "H", "db", "p-value", "eta^2_H"))

Kruskal-Wallis: suhu ~ pulau

Peubah	H	db	p-value	eta^2_H
suhu_c	3191.743	6	<1e-16	0.504

Untuk memeriksa kejanggalan suhu rendah pada beberapa wilayah, berikut provinsi dengan median suhu terendah (minimal 30 kecamatan):

df %>%
  group_by(provinsi) %>%
  summarise(n = n(),
            median_suhu   = median(suhu_c),
            IQR_suhu      = IQR(suhu_c),
            median_kelembapan = median(kelembapan_pct),
            median_angin  = median(kecepatan_angin),
            .groups = "drop") %>%
  filter(n >= 30) %>%
  arrange(median_suhu) %>%
  slice_head(n = 6) %>%
  kable(digits = 2, caption = "6 provinsi dengan median suhu terendah")

6 provinsi dengan median suhu terendah

provinsi	n	median_suhu	IQR_suhu	median_kelembapan	median_angin
Nusa Tenggara Timur	227	21.33	5.35	83	2.17
Kalimantan Tengah	136	21.61	0.90	95	0.74
Kalimantan Selatan	152	22.66	1.41	93	1.20
Kalimantan Utara	53	23.13	6.24	95	0.70
Kalimantan Timur	103	23.27	2.68	90	0.86
Maluku Utara	115	23.30	2.49	85	2.55

Interpretasi.

• Pulau menjelaskan sebagian besar variasi suhu: Kruskal-Wallis sangat signifikan dengan ukuran efek η² ≈ 0.50 (besar, jauh di atas ambang 0.14). Artinya pulau jauh lebih berpengaruh dibanding kategori cuaca (Pertanyaan 3, η² ≈ 0.06) dalam membedakan suhu.
• Jawa dan Sulawesi paling panas (median ≈ 29 °C, sebaran sempit), sedangkan Bali–Nusa Tenggara (≈ 22 °C) dan Kalimantan (≈ 23 °C) paling dingin.
• Papua memiliki sebaran terlebar (dari sekitar 6 °C hingga di atas 34 °C) karena wilayahnya mencakup pegunungan tinggi dan pesisir.
• Kehati-hatian (perlu dicek sebelum masuk infografis): tabel di atas menunjukkan beberapa provinsi Kalimantan bersuhu rendah untuk ukuran siang hari, dengan kelembapan 90% ke atas dan angin sangat tenang (< 1.3 m/s). Kalimantan Tengah bahkan sangat seragam: median 21.6 °C dengan IQR hanya sekitar 0.9 °C di 136 kecamatan. Pola ini mirip kondisi dini hari, sehingga ada kemungkinan data tidak direkam serempak menurut jam lokal, atau ada pengaruh faktor lain yang tidak kita miliki (misalnya ketinggian). Ini perlu dikonfirmasi ke sumber/proses scraping, karena bisa mencampur perbedaan iklim dengan perbedaan jam pengukuran.

6 Q6: Apakah Semakin ke Utara, Kelembapan Semakin Tinggi?

Lintang (latitude) bernilai negatif di selatan khatulistiwa dan positif di utaranya.

ggplot(df, aes(x = latitude, y = kelembapan_pct, color = pulau)) +
  geom_vline(xintercept = 0, linetype = "dashed", color = "grey40") +
  geom_point(alpha = 0.55, size = 1.5) +
  annotate("text", x = 0.2, y = 34, label = "Khatulistiwa", hjust = 0, size = 3.5, color = "grey30") +
  scale_color_manual(values = warna_pulau, name = "Pulau") +
  guides(color = guide_legend(override.aes = list(alpha = 1, size = 3), nrow = 1)) +
  labs(title = "Lintang vs Kelembapan Udara",
       subtitle = "Setiap titik = satu kecamatan",
       x = "Lintang (\u00b0, negatif = selatan)", y = "Kelembapan (%)", caption = sumber)

rho_lat <- cor(df$latitude, df$kelembapan_pct, method = "spearman")

df %>%
  group_by(pulau) %>%
  summarise(n = n(),
            `rho lintang-kelembapan` = cor(latitude, kelembapan_pct, method = "spearman"),
            .groups = "drop") %>%
  bind_rows(tibble(pulau = "Semua (gabungan)", n = nrow(df), `rho lintang-kelembapan` = rho_lat)) %>%
  kable(digits = 2, caption = "Korelasi Spearman lintang vs kelembapan: ")

Korelasi Spearman lintang vs kelembapan:

pulau	n	rho lintang-kelembapan
Bali & Nusa Tenggara	328	-0.36
Jawa	2143	-0.10
Kalimantan	618	-0.17
Maluku	232	0.19
Papua	494	-0.15
Sulawesi	584	0.18
Sumatera	1931	0.29
Semua (gabungan)	6330	0.49

Interpretasi.

• Secara gabungan, lintang dan kelembapan berkorelasi sedang positif (ρ ≈ +0.49): kecamatan yang lebih ke utara cenderung lebih lembap.
• Namun di dalam tiap pulau, hubungannya lemah dan arahnya tidak konsisten (dari sekitar −0.36 di Bali–Nusa Tenggara sampai +0.29 di Sumatera). Jadi korelasi gabungan terutama berasal dari perbedaan antar wilayah (mis. Sumatera dan Kalimantan lebih lembap dibanding Jawa dan Sulawesi), bukan dari efek “garis lintang” itu sendiri.
• Lintang tidak “menyebabkan” udara lembap. Ada faktor lain yang mungkin bisa menjadi penyebabnya, misalnya musim (Agustus umumnya kemarau di sebagian besar wilayah selatan khatulistiwa), ketinggian tempat, dan pola angin monsun. Faktor-faktor ini belum tercatat di data, sehingga hanya bisa disebut sebagai dugaan.

7 Kesimpulan

No	Pertanyaan	Temuan	Ukuran hubungan
1	Pasangan numerik terkuat	Suhu–kelembapan paling kuat (negatif); suhu juga berkaitan sedang dengan tekanan dan angin	Spearman −0.78; −0.56; +0.55
2	Suhu vs kelembapan	Negatif kuat dan tidak lurus; hampir hilang di kategori Hujan (campuran lokasi)	ρ −0.78 (Cerah −0.78, Berawan −0.82, Hujan −0.18)
3	Karakteristik per kategori cuaca	Hujan lebih dingin, lembap, dan tenang; Cerah ≈ Berawan	η² kecil–sedang (≈ 0.05–0.07 untuk suhu, kelembapan, angin)
4	Kategori cuaca per pulau	Hujan terkonsentrasi di Papua dan Sumatera; Jawa dan Nusa Tenggara kering	Cramér’s V ≈ 0.34
5	Suhu per pulau	Jawa dan Sulawesi terpanas; Bali–Nusa Tenggara dan Kalimantan terdingin	η² ≈ 0.50 (besar)
6	Lintang vs kelembapan	Positif sedang secara gabungan, tetapi tidak konsisten di dalam pulau	ρ +0.49 (gabungan)