1 Pendahuluan

1.1 Tentang Studi Kasus

Project menerapkan teknik eksplorasi data lingkungan yang dibahas pada Pertemuan 2 — statistik deskriptif, visualisasi, deteksi outlier, dan autokorelasi (temporal & spasial), pada studi kasus luas kebakaran hutan dan lahan (karhutla) dan curah hujan di 34 provinsi Indonesia.

Dataset utama (Data_Project) berisi 190 baris (34 dan atau 38 provinsi × 5 tahun, 2022–2026), memuat rata-rata curah hujan tahunan (rfh, mm), luas karhutla tahunan (luas_karhutla_ha), serta koordinat (lat, long) tiap provinsi. Data disusun dari dua sumber mentah pada file yang sama: sheet Data Karhutla (rekap luas karhutla per provinsi per tahun) dan sheet Curah Hujan 2022-2026 (data curah hujan dekadal/10-harian per wilayah administrasi, format rfh/rfq gaya WFP-FEWS NET), yang digabungkan dengan koordinat pada sheet Bujur lintang.

1.2 Rumusan Pertanyaan Analisis

  1. Bagaimana karakteristik statistik deskriptif luas karhutla dan curah hujan secara nasional maupun per provinsi?
  2. Bagaimana pola distribusi, tren tahunan, dan sebaran spasial karhutla dapat digambarkan melalui visualisasi?
  3. Apakah terdapat provinsi/tahun dengan luas karhutla yang tergolong outlier?
  4. Apakah terdapat indikasi autokorelasi temporal (antarwaktu) pada curah hujan dan autokorelasi spasial (antarprovinsi berdekatan) pada karhutla?

Keempat pertanyaan ini dijawab bertahap melalui delapan langkah analisis di bawah ini.

2 Persiapan

2.1 Memuat Package

library(readxl)         # membaca file Excel
library(dplyr)           # manipulasi data
library(tidyr)            # merapikan data
library(ggplot2)         # visualisasi
library(sf)                # objek spasial
library(tmap)             # peta tematik
library(spdep)            # autokorelasi spasial (Moran's I)
library(rnaturalearth)  # garis batas negara Indonesia (basemap, tanpa perlu shapefile manual)

2.2 Langkah 1 — Import dan Struktur Data

# Input data
data_karhutla <- read_excel("C:/Users/user/OneDrive - untirta.ac.id/Datacoba1 (1).xlsx", sheet = "Data_Project")

str(data_karhutla)
## tibble [190 × 8] (S3: tbl_df/tbl/data.frame)
##  $ kode_provinsi   : chr [1:190] "11" "11" "11" "11" ...
##  $ provinsi        : chr [1:190] "Aceh" "Aceh" "Aceh" "Aceh" ...
##  $ tahun           : num [1:190] 2022 2023 2024 2025 2026 ...
##  $ rfh             : chr [1:190] "87.065249777777808" "76.088721000000007" "72.041313486111122" "84.924350111111096" ...
##  $ rfq             : chr [1:190] "119.20335880555558" "101.73461622222224" "98.849901000000003" "110.10899866666668" ...
##  $ luas_karhutla_ha: num [1:190] 3716 1937 7257 8925 4070 ...
##  $ lat             : chr [1:190] "4.3685499999999999" "4.3685499999999999" "4.3685499999999999" "4.3685499999999999" ...
##  $ long            : chr [1:190] "97.025300000000001" "97.025300000000001" "97.025300000000001" "97.025300000000001" ...
data_karhutla <- data_karhutla %>%
  dplyr::mutate(dplyr::across(c(rfh, rfq, lat, long), as.numeric))

str(data_karhutla)
## tibble [190 × 8] (S3: tbl_df/tbl/data.frame)
##  $ kode_provinsi   : chr [1:190] "11" "11" "11" "11" ...
##  $ provinsi        : chr [1:190] "Aceh" "Aceh" "Aceh" "Aceh" ...
##  $ tahun           : num [1:190] 2022 2023 2024 2025 2026 ...
##  $ rfh             : num [1:190] 87.1 76.1 72 84.9 67.1 ...
##  $ rfq             : num [1:190] 119.2 101.7 98.8 110.1 103.9 ...
##  $ luas_karhutla_ha: num [1:190] 3716 1937 7257 8925 4070 ...
##  $ lat             : num [1:190] 4.37 4.37 4.37 4.37 4.37 ...
##  $ long            : num [1:190] 97 97 97 97 97 ...
dim(data_karhutla)
## [1] 190   8
# Memeriksa data hilang per kolom
colSums(is.na(data_karhutla))
##    kode_provinsi         provinsi            tahun              rfh 
##                0                0                0               20 
##              rfq luas_karhutla_ha              lat             long 
##               20                0               20               20
data_karhutla %>%
  dplyr::filter(is.na(lat)) %>%
  dplyr::distinct(provinsi)
## # A tibble: 4 × 1
##   provinsi        
##   <chr>           
## 1 Papua Barat     
## 2 Papua Barat Daya
## 3 Papua Pegunungan
## 4 Papua Selatan

2.2.1 Interpretasi

Dataset berisi 190 baris dan 8 kolom, sesuai rancangan 38 provinsi × 5 tahun (2022–2026). Kolom luas_karhutla_ha dan rfh sudah bertipe numerik dan provinsi bertipe karakter.

Dataset tidak sepenuhnya lengkap: kolom rfh, rfq, lat, dan long masing-masing memiliki 20 nilai hilang (20 = 4 provinsi × 5 tahun), seluruhnya berasal dari empat provinsi hasil pemekaran Papua 2022 (Papua Barat, Papua Barat Daya, Papua Pegunungan, Papua Selatan) yang belum tersedia data curah hujan maupun koordinatnya pada sumber data ini kemungkinan karena sumber WFP-FEWS NET dan basemap yang dipakai belum memperbarui pembagian wilayah administratifnya. Sebaliknya, luas_karhutla_ha tetap lengkap untuk keempat provinsi tersebut karena berasal dari sumber data yang berbeda (sheet Data Karhutla). Ketidaklengkapan ini ditangani dengan mengecualikan baris terkait hanya dari analisis yang membutuhkan rfh atau koordinat, bukan menghapusnya dari keseluruhan dataset

3 Statistik Deskriptif

3.1 Langkah 2 — Statistik Deskriptif

# Statistik deskriptif keseluruhan — luas karhutla
summary(data_karhutla$luas_karhutla_ha)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
##      0.0    384.4   2404.6  12129.1   8288.4 190394.6
sd(data_karhutla$luas_karhutla_ha)
## [1] 28932.59
IQR(data_karhutla$luas_karhutla_ha)
## [1] 7903.957
# Statistik deskriptif keseluruhan — curah hujan (rfh)
summary(data_karhutla$rfh)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.     NAs 
##   42.21   64.20   78.28   78.92   94.23  122.16      20
sd(data_karhutla$rfh, na.rm = TRUE)
## [1] 18.59752

Rata-rata luas karhutla nasional adalah 12.129 ha per provinsi per tahun, jauh di atas mediannya (2.405 ha), dengan simpangan baku yang sangat besar (28.933 ha) hampir 2,4 kali nilai rata-ratanya sendiri. Selisih mean–median yang lebar ini adalah indikasi kuat distribusi right-skewed: sebagian kecil provinsi/tahun dengan karhutla ekstrem menarik rata-rata jauh di atas kondisi yang tipikal, sehingga median lebih representatif sebagai ukuran pemusatan untuk data ini.

# Statistik deskriptif per provinsi (diurutkan dari rerata karhutla tertinggi)
rekap_provinsi <- data_karhutla %>%
  group_by(provinsi) %>%
  summarise(
    n       = n(),
    rerata  = mean(luas_karhutla_ha),
    median  = median(luas_karhutla_ha),
    sd      = sd(luas_karhutla_ha),
    min     = min(luas_karhutla_ha),
    maks    = max(luas_karhutla_ha)
  ) %>%
  arrange(desc(rerata))

rekap_provinsi
## # A tibble: 38 × 7
##    provinsi                n rerata  median     sd    min    maks
##    <chr>               <int>  <dbl>   <dbl>  <dbl>  <dbl>   <dbl>
##  1 Nusa Tenggara Timur     5 84848. 102537. 40929. 19565. 118409.
##  2 Kalimantan Barat        5 44571.  26703. 38139. 21836  111848.
##  3 Kalimantan Selatan      5 43147.   8020. 82421.   429  190395.
##  4 Papua Selatan           5 42388.  25839. 61711.     0  150813.
##  5 Kalimantan Tengah       5 38270.   7634. 71439.  1554  165896.
##  6 Nusa Tenggara Barat     5 37743.  30567  19197. 17406.  66716.
##  7 Sumatera Selatan        5 31819.   5940. 56290.  1926. 132083.
##  8 Maluku                  5 22194.  14954  13720. 13405.  45999.
##  9 Jawa Timur              5 18868.  10835. 18345.  2380   49498.
## 10 Kalimantan Timur        5 14945.   9573. 16214.   373   39494.
## # ℹ 28 more rows

Rekap per provinsi menunjukkan adanya variasi spasial yang cukup jelas dalam luas kebakaran hutan dan lahan (karhutla). Provinsi-provinsi dengan rerata luas karhutla tertinggi didominasi oleh wilayah Kalimantan (Kalimantan Barat, Kalimantan Selatan, dan Kalimantan Tengah), Sumatera Selatan, Nusa Tenggara Timur, serta Papua Selatan. Sebaliknya, beberapa provinsi dengan karakteristik wilayah yang lebih urban atau memiliki area berhutan dan bergambut yang relatif lebih terbatas menunjukkan rerata luas karhutla yang lebih rendah.

# Rerata per tahun — melihat tahun mana yang paling parah secara nasional
rekap_tahun <- data_karhutla %>%
  group_by(tahun) %>%
  summarise(
    rerata_karhutla = mean(luas_karhutla_ha),
    rerata_hujan    = mean(rfh, na.rm = TRUE)
  )
rekap_tahun
## # A tibble: 5 × 3
##   tahun rerata_karhutla rerata_hujan
##   <dbl>           <dbl>        <dbl>
## 1  2022           5392.         90.8
## 2  2023          30558.         66.5
## 3  2024           9916.         79.6
## 4  2025           9464.         92.7
## 5  2026           5316.         65.0

3.1.1 Interpretasi

Rekap per tahun menunjukkan 2023 sebagai tahun terparah secara nasional (rerata karhutla 30.558 ha, lebih dari 3 kali tahun lainnya), justru bersamaan dengan curah hujan rata-rata terendah kedua (66,5 mm).

4 Visualisasi Data

4.1 Langkah 3 — Histogram dan Boxplot

# Histogram distribusi luas karhutla keseluruhan
ggplot(data_karhutla, aes(x = luas_karhutla_ha)) +
  geom_histogram(bins = 25, fill = "#147D64", color = "white") +
  labs(title = "Distribusi Luas Karhutla — Indonesia 2022–2026",
       x = "Luas Karhutla (ha)", y = "Frekuensi") +
  theme_minimal()

4.1.1 Interpretasi

Histogram memperlihatkan distribusi luas karhutla yang sangat menceng ke kanan (positively skewed). Sebagian besar observasi berada pada rentang rendah (<20.000 ha), sedangkan hanya sedikit observasi yang memiliki luas karhutla sangat besar hingga mendekati 190.000 ha. Kondisi ini menunjukkan bahwa kejadian karhutla berskala besar relatif jarang terjadi, namun memberikan kontribusi yang signifikan terhadap variasi data secara keseluruhan. Oleh karena itu, nilai rata-rata cenderung lebih tinggi dibandingkan median akibat pengaruh pengamatan ekstrem pada ekor kanan distribusi.

# Boxplot luas karhutla per provinsi (diurutkan berdasarkan median)
ggplot(data_karhutla, aes(x = reorder(provinsi, luas_karhutla_ha, median),
                            y = luas_karhutla_ha)) +
  geom_boxplot(fill = "#FFC200") +
  coord_flip() +
  labs(title = "Sebaran Luas Karhutla per Provinsi (2022–2026)",
       x = "Provinsi", y = "Luas Karhutla (ha)") +
  theme_minimal()

4.1.2 Interpretasi

Boxplot per provinsi (diurutkan dari median terendah ke tertinggi, ditampilkan horizontal karena banyaknya kategori) mempertegas pola regional: provinsi di Kalimantan, Sumatera bagian selatan, dan Nusa Tenggara memiliki kotak (IQR) dan whisker jauh lebih lebar dibanding provinsi Jawa dan kepulauan kecil artinya karhutla di wilayah tersebut tidak hanya lebih tinggi rata-ratanya, tetapi juga jauh lebih variatif antartahun. Titik-titik terpisah di atas whisker atas pada Kalimantan Selatan, Kalimantan Tengah, dan beberapa provinsi lain adalah representasi visual dari nilai ekstrem yang akan diuji secara formal pada Langkah 6.

4.2 Langkah 4 — Diagram Deret Waktu

# Fokus pada 8 provinsi dengan rerata karhutla tertinggi agar grafik tetap terbaca
provinsi_top8 <- rekap_provinsi$provinsi[1:8]

data_karhutla %>%
  dplyr::filter(provinsi %in% provinsi_top8) %>%
  ggplot(aes(x = tahun, y = luas_karhutla_ha, color = provinsi)) +
  geom_line(linewidth = 1) +
  geom_point(size = 2) +
  scale_x_continuous(breaks = 2022:2026) +
  labs(title = "Tren Tahunan Luas Karhutla — 8 Provinsi Tertinggi",
       x = "Tahun", y = "Luas Karhutla (ha)", color = "Provinsi") +
  theme_minimal()

4.2.1 Interpretasi

Grafik deret waktu pada 8 provinsi ber-karhutla tertinggi memperlihatkan lonjakan tajam dan serentak pada 2023 di hampir semua provinsi tersebut paling ekstrem di Kalimantan Selatan, Kalimantan Tengah, dan Papua Selatan — sebelum kembali menurun tajam pada 2024–2026 mendekati level 2022.

# Scatter plot: rerata curah hujan vs rerata karhutla per tahun (agregat nasional)
ggplot(rekap_tahun, aes(x = rerata_hujan, y = rerata_karhutla)) +
  geom_point(color = "blue", size = 3) +
  geom_text(aes(label = tahun), vjust = -1) +
  geom_smooth(method = "lm", color = "red", se = TRUE) +
  labs(title = "Rerata Curah Hujan vs Rerata Luas Karhutla per Tahun (Nasional)",
       x = "Rerata Curah Hujan (mm)", y = "Rerata Luas Karhutla (ha)") +
  theme_minimal()

cor(rekap_tahun$rerata_hujan, rekap_tahun$rerata_karhutla)
## [1] -0.4604231

4.2.2 Interpretasi

Scatter plot menunjukkan kecenderungan hubungan negatif antara rerata curah hujan dan rerata luas karhutla pada tingkat nasional selama periode 2022–2026. Hal ini ditunjukkan oleh garis regresi yang menurun, yang mengindikasikan bahwa tahun dengan curah hujan lebih tinggi cenderung memiliki rerata luas karhutla yang lebih rendah. Sebaliknya, ketika curah hujan relatif rendah, luas karhutla cenderung lebih besar.

Namun demikian, hubungan tersebut belum terlihat kuat karena jumlah pengamatan terbatas dan pita kepercayaan di sekitar garis regresi masih cukup lebar. Selain itu, terdapat variasi yang cukup besar antar tahun, seperti pada tahun 2023 yang mencatat luas karhutla jauh lebih tinggi dibandingkan tahun-tahun lainnya. Oleh karena itu, grafik ini sebagai indikasi awal adanya hubungan negatif antara curah hujan dan luas karhutla daripada sebagai bukti hubungan yang pasti.

Korelasi antara rerata curah hujan dan rerata luas karhutla pada agregat tahunan bernilai negatif sedang (r = -0,46) tahun dengan curah hujan lebih rendah cenderung memiliki luas karhutla lebih besar, sejalan dengan mekanisme musim kemarau panjang yang mengeringkan bahan bakar vegetasi dan meningkatkan kerentanan kebakaran. Meski demikian, korelasi ini hanya dihitung dari 5 titik data tahunan, sehingga nilainya bersifat indikatif dan rentan terhadap pengaruh titik individual (khususnya 2023) — bukan bukti hubungan sebab-akibat yang konklusif, dan idealnya divalidasi dengan data bulanan/dekadal atau periode tahun yang lebih panjang.

4.3 Langkah 5

# Garis batas Indonesia diambil otomatis (tidak perlu shapefile provinsi manual)

indo <- ne_countries(country = "Indonesia", scale = "medium", returnclass = "sf")

# Titik rerata karhutla per provinsi (2022-2026), provinsi tanpa koordinat otomatis dikecualikan


titik_provinsi <- data_karhutla %>%
  dplyr::filter(!is.na(lat), !is.na(long)) %>%
  tidyr::drop_na(lat, long) %>%   # jaring pengaman kedua
  dplyr::group_by(provinsi, lat, long) %>%
  dplyr::summarise(rerata_karhutla = mean(luas_karhutla_ha), .groups = "drop")

# Sanity check: pastikan tidak ada NA tersisa sebelum dibuat objek spasial
stopifnot(sum(is.na(titik_provinsi$lat)) == 0, sum(is.na(titik_provinsi$long)) == 0)

titik_provinsi <- titik_provinsi %>%
  st_as_sf(coords = c("long", "lat"), crs = 4326, remove = FALSE)

tm_shape(indo) +
  tm_borders(col = "grey50") +
tm_shape(titik_provinsi) +
  tm_symbols(size = "rerata_karhutla", col = "rerata_karhutla",
             palette = "YlOrRd", scale = 2,
             title.col = "Rerata Karhutla (ha)",
             title.size = "Rerata Karhutla (ha)") +
  tm_layout(title = "Sebaran Rerata Luas Karhutla per Provinsi (2022-2026)",
            legend.outside = TRUE)

4.4 Langkah 5 Peta Tematik Choropleth

## Menggunakan shapefile provinsi asli (indonesia.shp/.shx/.dbf/.prj),

peta_provinsi <- sf::st_read("C:/Users/user/OneDrive - untirta.ac.id/Documents/Semester 5/shp indo", quiet = TRUE)

# Cek nama provinsi di shapefile vs di data
sort(peta_provinsi$name)
##  [1] "Aceh"                "Bali"                "Bangka-Belitung"    
##  [4] "Banten"              "Bengkulu"            "Gorontalo"          
##  [7] "Jakarta Raya"        "Jambi"               "Jawa Barat"         
## [10] "Jawa Tengah"         "Jawa Timur"          "Kalimantan Barat"   
## [13] "Kalimantan Selatan"  "Kalimantan Tengah"   "Kalimantan Timur"   
## [16] "Kepulauan Riau"      "Lampung"             "Maluku"             
## [19] "Maluku Utara"        "North Kalimantan"    "Nusa Tenggara Barat"
## [22] "Nusa Tenggara Timur" "Papua"               "Papua Barat"        
## [25] "Riau"                "Sulawesi Barat"      "Sulawesi Selatan"   
## [28] "Sulawesi Tengah"     "Sulawesi Tenggara"   "Sulawesi Utara"     
## [31] "Sumatera Barat"      "Sumatera Selatan"    "Sumatera Utara"     
## [34] "Yogyakarta"
# --- Harmonisasi nama provinsi ---
# Shapefile ini masih memakai batas 34 provinsi (SEBELUM pemekaran Papua 2022),
# jadi 5 provinsi baru di data kita (Papua, Papua Tengah, Papua Pegunungan,
# Papua Selatan, Papua Barat Daya) tidak punya poligon sendiri.
# Solusi: agregasi data provinsi baru ke wilayah induknya agar tetap bisa dipetakan.
crosswalk_papua <- c(
  "Papua"              = "Papua",
  "Papua Tengah"        = "Papua",
  "Papua Pegunungan"    = "Papua",
  "Papua Selatan"       = "Papua",
  "Papua Barat"         = "Papua Barat",
  "Papua Barat Daya"    = "Papua Barat"
)

titik_provinsi_peta <- data_karhutla %>%
  dplyr::mutate(
    provinsi_peta = dplyr::case_when(
      provinsi %in% names(crosswalk_papua) ~ crosswalk_papua[provinsi],
      provinsi == "DKI Jakarta"       ~ "Jakarta Raya",       # beda penulisan di shapefile
      provinsi == "Bangka Belitung"   ~ "Bangka-Belitung",    # beda penulisan di shapefile
      provinsi == "Kalimantan Utara"  ~ "North Kalimantan",   # beda penulisan di shapefile
      TRUE ~ provinsi
    )
  ) %>%
  dplyr::group_by(provinsi_peta) %>%
  dplyr::summarise(rerata_karhutla = mean(luas_karhutla_ha), .groups = "drop")

# Gabungkan data dengan poligon provinsi
peta_karhutla <- peta_provinsi %>%
  dplyr::left_join(titik_provinsi_peta, by = c("name" = "provinsi_peta"))

# Peta choropleth: seluruh wilayah provinsi diwarnai, bukan hanya titik pusatnya
tm_shape(peta_karhutla) +
  tm_polygons(col = "rerata_karhutla",
              palette = "YlOrRd", style = "jenks",
              title = "Rerata Karhutla (ha)",
              border.col = "grey40", colorNA = "grey90",
              textNA = "Data tidak tersedia") +
  tm_layout(title = "Sebaran Rerata Luas Karhutla per Provinsi (2022-2026)",
            legend.outside = TRUE)

4.4.1 Interpretasi

Peta tematik memperlihatkan pola spasial yang cukup jelas: simbol besar dan berwarna gelap terkonsentrasi di Kalimantan dan Nusa Tenggara Timur, sedangkan Pulau Jawa hampir seluruhnya bersimbol kecil dan terang. Pola “provinsi berdekatan cenderung mirip” ini menjadi petunjuk visual awal adanya autokorelasi spasial positif yang akan diuji formal pada Langkah 8.

5 Deteksi Outlier

5.1 Langkah 6 — Deteksi Outlier

# Metode IQR
Q1 <- quantile(data_karhutla$luas_karhutla_ha, 0.25)
Q3 <- quantile(data_karhutla$luas_karhutla_ha, 0.75)
IQR_val <- Q3 - Q1
batas_bawah <- Q1 - 1.5 * IQR_val
batas_atas  <- Q3 + 1.5 * IQR_val

cat("Batas bawah:", round(batas_bawah, 1), " | Batas atas:", round(batas_atas, 1), "\n")
## Batas bawah: -11471.5  | Batas atas: 20144.3
outlier_iqr <- data_karhutla %>%
  dplyr::filter(luas_karhutla_ha < batas_bawah | luas_karhutla_ha > batas_atas) %>%
  dplyr::arrange(desc(luas_karhutla_ha))
outlier_iqr %>% select(provinsi, tahun, luas_karhutla_ha)
## # A tibble: 27 × 3
##    provinsi            tahun luas_karhutla_ha
##    <chr>               <dbl>            <dbl>
##  1 Kalimantan Selatan   2023          190395.
##  2 Kalimantan Tengah    2023          165896.
##  3 Papua Selatan        2023          150813.
##  4 Sumatera Selatan     2023          132083.
##  5 Nusa Tenggara Timur  2025          118409.
##  6 Nusa Tenggara Timur  2024          113090.
##  7 Kalimantan Barat     2023          111848.
##  8 Nusa Tenggara Timur  2023          102537.
##  9 Nusa Tenggara Timur  2022           70637 
## 10 Nusa Tenggara Barat  2023           66716.
## # ℹ 17 more rows
# Metode Z-score
data_karhutla$z_karhutla <- as.numeric(scale(data_karhutla$luas_karhutla_ha))
outlier_z <- data_karhutla %>%
  dplyr::filter(abs(z_karhutla) > 3) %>%
  dplyr::arrange(desc(z_karhutla))
outlier_z %>% select(provinsi, tahun, luas_karhutla_ha, z_karhutla)
## # A tibble: 8 × 4
##   provinsi            tahun luas_karhutla_ha z_karhutla
##   <chr>               <dbl>            <dbl>      <dbl>
## 1 Kalimantan Selatan   2023          190395.       6.16
## 2 Kalimantan Tengah    2023          165896.       5.31
## 3 Papua Selatan        2023          150813.       4.79
## 4 Sumatera Selatan     2023          132083.       4.15
## 5 Nusa Tenggara Timur  2025          118409.       3.67
## 6 Nusa Tenggara Timur  2024          113090.       3.49
## 7 Kalimantan Barat     2023          111848.       3.45
## 8 Nusa Tenggara Timur  2023          102537.       3.12

5.1.1 Interpretasi

Dengan aturan IQR, batas atas kewajaran adalah 20.144,3 ha; sebanyak 27 dari 190 observasi (≈14,2%) melampauinya. Proporsi ini terlihat besar untuk standar IQR konvensional, namun wajar mengingat data memang sangat right-skewed — pada distribusi seperti ini, aturan 1,5×IQR cenderung menandai banyak titik karena IQR-nya sendiri relatif sempit dibanding rentang datanya. Metode Z-score yang lebih konservatif (|z| > 3) menandai hanya 8 observasi paling ekstrem.

Kedua metode sepakat bahwa observasi paling ekstrem berasal dari tahun 2023, tersebar di beberapa provinsi rawan karhutla (Kalimantan Selatan, Kalimantan Tengah, Papua Selatan, Sumatera Selatan, Nusa Tenggara Timur, Kalimantan Barat).

6 Autokorelasi

6.1 Langkah 7 — Autokorelasi Temporal (ACF)

Deret waktu tahunan karhutla per provinsi hanya memiliki 5 titik (2022–2026) terlalu pendek untuk ACF yang andal. Sebagai gantinya, autokorelasi temporal didemonstrasikan pada data curah hujan dekadal (10-harian) dari sheet Curah Hujan 2022-2026, yang tersedia untuk 34 provinsi dengan ±168 titik waktu per provinsi (Januari 2022 – Agustus 2026) — jauh lebih memadai untuk melihat pola autokorelasi/musiman secara meyakinkan.

# Curah hujan dekadal level provinsi (adm_level == 1), contoh: Nusa Tenggara Timur (PCODE "ID53")
# —provinsi dengan rerata karhutla tertinggi pada rekap_provinsi di Langkah 2
data_hujan_dekadal <- read_excel("C:/Users/user/OneDrive - untirta.ac.id/Datacoba1 (1).xlsx", sheet = "Curah Hujan 2022-2026")

ntt_dekadal <- data_hujan_dekadal %>%
  dplyr::filter(adm_level == 1, PCODE == "ID53") %>%
  dplyr::arrange(date)

ts_ntt <- ts(ntt_dekadal$rfh, start = c(2022, 1), frequency = 36)  # 36 dekade per tahun

acf(ts_ntt, lag.max = 40, main = "ACF Curah Hujan Dekadal — Nusa Tenggara Timur")

6.1.1 Interpretasi

Grafik ACF menunjukkan bahwa pola curah hujan di Nusa Tenggara Timur sangat dipengaruhi oleh musim. Nilai autokorelasi yang tinggi pada lag awal menunjukkan bahwa jika suatu periode memiliki curah hujan tinggi, maka periode-periode setelahnya juga cenderung memiliki curah hujan tinggi. Sebaliknya, pada lag menengah autokorelasi menjadi negatif, yang mencerminkan perbedaan antara musim hujan dan musim kemarau. Setelah mendekati satu tahun, autokorelasi kembali positif, menandakan bahwa pola curah hujan tahun ini cenderung berulang pada periode yang sama di tahun berikutnya. Dengan kata lain, curah hujan di NTT memiliki siklus musiman tahunan yang kuat dan relatif konsisten dari tahun ke tahun.

6.2 Langkah 8 — Autokorelasi Spasial (Moran’s I)

# Matriks pembobot spasial: k = 4 tetangga terdekat, dari rerata karhutla per provinsi (Langkah 5)
# (provinsi tanpa koordinat otomatis dikecualikan — lihat catatan pada Langkah 5)
koordinat <- cbind(titik_provinsi$long, titik_provinsi$lat)
knn4 <- knearneigh(koordinat, k = 4)
wij <- nb2listw(knn2nb(knn4), style = "W")

moran.test(titik_provinsi$rerata_karhutla, listw = wij)
## 
##  Moran I test under randomisation
## 
## data:  titik_provinsi$rerata_karhutla  
## weights: wij    
## 
## Moran I statistic standard deviate = 1.8432, p-value = 0.03265
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.140049158      -0.030303030       0.008542134

6.2.1 Interpretasi

Hasil uji Moran’s I menghasilkan nilai sebesar 0,140 dengan p-value 0,0327. Karena p-value lebih kecil dari 0,05, terdapat autokorelasi spasial positif yang signifikan pada data rerata luas karhutla antarprovinsi di Indonesia.

Autokorelasi spasial positif berarti provinsi yang berdekatan secara geografis cenderung memiliki tingkat karhutla yang mirip. Dengan kata lain, wilayah dengan karhutla tinggi cenderung berada di dekat wilayah lain yang juga memiliki karhutla tinggi, sedangkan wilayah dengan karhutla rendah cenderung berdekatan dengan wilayah yang juga rendah.

Nilai Moran’s I sebesar 0,140 menunjukkan bahwa hubungan spasial tersebut ada, tetapi tidak terlalu kuat. Hal ini mengindikasikan bahwa selain faktor lokasi geografis, masih terdapat faktor-faktor lain yang turut memengaruhi variasi karhutla antarprovinsi, seperti kondisi iklim, tutupan lahan, dan aktivitas manusia. Temuan kuantitatif ini mengonfirmasi pola yang sudah terlihat secara visual pada peta tematik Langkah 5, dan sejalan dengan

Hukum Pertama Tobler

7 Kesimpulan

1. Karakteristik statistik deskriptif. Secara nasional, luas kebakaran hutan dan lahan (karhutla) per provinsi per tahun memiliki rata-rata 12.129 ha, sedangkan mediannya hanya 2.405 ha (sd = 28.933 ha). Perbedaan yang besar antara rata-rata dan median menunjukkan bahwa distribusi data sangat menceng ke kanan (right-skewed), dengan sebagian kecil observasi bernilai sangat tinggi yang menarik rata-rata ke atas. Secara spasial, rerata karhutla tertinggi ditemukan di beberapa provinsi Kalimantan, Sumatera Selatan, Nusa Tenggara Timur, dan Papua Selatan, sedangkan provinsi di Pulau Jawa umumnya memiliki nilai yang lebih rendah.

2. Pola distribusi, tren, dan sebaran spasial. Histogram dan boxplot menunjukkan bahwa sebagian besar observasi memiliki luas karhutla relatif rendah, sementara hanya sedikit observasi yang memiliki nilai sangat besar. Analisis tren tahunan menunjukkan bahwa tahun 2023 merupakan tahun dengan rerata karhutla tertinggi secara nasional, lebih dari tiga kali lipat dibandingkan tahun lainnya. Selain itu, terdapat indikasi hubungan negatif antara curah hujan dan luas karhutla, di mana tahun dengan curah hujan lebih rendah cenderung memiliki luas karhutla yang lebih besar. Peta tematik juga memperlihatkan adanya kecenderungan pengelompokan wilayah dengan tingkat karhutla yang serupa.

3. Provinsi/tahun dengan karhutla outlier. Deteksi outlier menggunakan metode IQR dan Z-score menunjukkan adanya sejumlah observasi ekstrem, terutama pada tahun 2023. Kedua metode secara konsisten mengidentifikasi provinsi-provinsi seperti Kalimantan Selatan, Kalimantan Tengah, Papua Selatan, Sumatera Selatan, Nusa Tenggara Timur, dan Kalimantan Barat sebagai lokasi dengan kejadian karhutla yang sangat tinggi dibandingkan pola umum data. Observasi ekstrem ini merepresentasikan fenomena nyata sehingga lebih tepat diperlakukan sebagai informasi penting daripada dihapus dari analisis.

4. Autokorelasi temporal dan spasial Analisis ACF pada curah hujan dekadal di Nusa Tenggara Timur menunjukkan adanya pola musiman tahunan yang kuat dan konsisten, ditandai oleh pengulangan pola curah hujan dari tahun ke tahun. Sementara itu, hasil uji Moran’s I (I = 0,140; p = 0,033) menunjukkan adanya autokorelasi spasial positif yang signifikan pada rerata karhutla antarprovinsi. Artinya, provinsi yang berdekatan secara geografis cenderung memiliki tingkat karhutla yang mirip. Meskipun demikian, kekuatan autokorelasi spasial tersebut relatif lemah hingga sedang, sehingga faktor-faktor lain seperti iklim, tutupan lahan, dan aktivitas manusia juga berperan dalam membentuk variasi karhutla di Indonesia.