Lewati ke isi
Profil penulisSeri Buku GIS Kehutanan dan Pertanian/ M2
Tampilkan bagian untuk:

BAB 5: Menebak di Antara Titik: Interpolasi

#Studi kasus: "Hujan di petak yang tidak punya stasiun"

KPH Contoh punya 30 stasiun hujan, tetapi 25 petaknya tidak punya satu pun. Untuk menilai kesesuaian lahan, Kepala Seksi butuh peta curah hujan yang menutup seluruh wilayah. Anda harus menebak nilai di tempat yang tidak diukur. Pertanyaannya: metode tebakan mana yang paling bisa dipercaya?

#Konsep: interpolasi dalam tiga kalimat

Interpolasi menebak nilai di tempat yang tidak diukur dari nilai di tempat yang diukur. Bayangkan menebak suhu kamar dari tiga termometer di sudut-sudutnya. Tebakan yang baik memberi bobot lebih besar kepada termometer yang dekat.

Istilah baru:

  • IDW (inverse distance weighting): bobot berbanding terbalik dengan jarak, dipangkatkan p. Makin besar p, makin dekat titik makin dominan.
  • TIN: jaring segitiga yang menghubungkan titik, nilai di dalam segitiga dihitung lurus dari ketiga sudutnya.
  • Kriging: tebakan berbobot yang bobotnya dihitung dari pola kemiripan antartitik (semivariogram).
  • Titik latih dan titik uji: sebagian titik dipakai menebak, sebagian disimpan untuk menilai tebakan.
  • RMSE: rata-rata besar selisih tebakan dan kenyataan. Makin kecil makin baik.

Untuk latihan ini, 30 stasiun dibagi: 24 titik latih dan 6 titik uji. Nilai 6 stasiun uji tidak boleh dipakai saat menebak. Anda hanya mengintip mereka saat menilai.

#Bagian A: QGIS

#Bagian A: Interpolasi di QGIS

Semua metode di bawah dijalankan oleh skrip m2_05a_interpolasi.py. Keluarannya ada di folder hasil.

#A1. IDW

  1. Buka Processing ► Toolbox, cari IDW interpolation.
  2. Pada Input layer(s), klik tombol +, pilih layer stasiun latih, atribut ch_tahunan_mm, dan jenis Points.
  3. Distance coefficient P: 2. Extent: samakan dengan DEM_10m (Calculate from layer). Output raster size: 10.
  4. Klik Run.

Ulangi dengan P = 1 dan 3. Perhatikan: kolom pertama di daftar atribut yang dilihat QGIS bisa berupa fid, sehingga nomor kolom bergeser satu. Pastikan atribut yang dipilih memang ch_tahunan_mm. (Saat membuat skrip, penulis sempat memilih kolom keliru dan hasilnya memuat angka elevasi.)

#A2. TIN

Alat TIN interpolation. Pilih metode Linear atau Clough-Toucher (cubic). TIN hanya menghitung di dalam jaring segitiga. Titik di luar tepi jaring tidak mendapat nilai.

#A3. Spline

QGIS memuat alat GRASS v.surf.bspline (B-spline) dan v.surf.rst (spline dengan tegangan). Nilai awal skrip: v.surf.bspline dengan langkah 500 m dan metode bikubik, v.surf.rst dengan tension 40 dan smooth 0,5. Parameter ini hasil pilihan penulis, bukan angka baku.

#A4. Kriging

QGIS 4.0.2 yang diuji penulis tidak punya alat kriging bawaan di toolbox. Skrip penulis menyediakan kriging biasa dengan NumPy:

  1. Hitung selisih nilai antarpasang titik, kelompokkan menurut jarak (semivariogram empiris).
  2. Cocokkan model eksponensial ke kelompok itu.
  3. Selesaikan sistem persamaan kriging untuk tiap sel.

Hasil penyesuaian untuk data latih: nugget 0, sill parsial sekitar 66.870 mm², jangkauan 5.000 m. Angka jangkauan menempel pada batas atas yang penulis izinkan dalam penyesuaian. Artinya pola kemiripan berlangsung sangat jauh, wajar karena hujan di data ini naik mengikuti elevasi (suatu tren). Di luar skrip ini, kriging dapat memakai plugin atau alat pihak ketiga. [CEK]

#A5. Regresi dengan elevasi lalu sisa dengan IDW

Hujan sering naik bersama ketinggian. Jadi Anda menyesuaikan garis lurus hujan terhadap elevasi dan koordinat X, lalu menginterpolasi sisanya (residu) dengan IDW p = 2 dan menjumlahkannya kembali. Hasil penyesuaian dari 24 titik latih:

KODE
CH = 1944,2 + 3,04 x elevasi + 0,1371 x (X - X0)      R² = 0,912

Artinya: tiap naik 1 m, hujan tahunan naik sekitar 3 mm. Elevasi tiap sel diambil dari DEM. Metode ini dikerjakan lewat skrip, bukan lewat tombol QGIS.

#A6. Menilai: siapa yang paling dekat dengan kenyataan?

Skrip mengambil nilai tiap raster hasil di enam stasiun uji dan menghitung selisih terhadap pengamatan. Hasil uji:

MetodeRMSE (mm/th)MAEBiasTitik uji terjangkau
IDW p = 1183,4175,2+150,16
IDW p = 2126,5119,7+97,46
IDW p = 394,683,8+61,06
TIN linear55,940,5+38,62
TIN Clough-Tocher98,374,6-45,63
B-spline bikubik100,784,1+9,56
RST (spline tegangan)89,976,2+41,06
Kriging biasa83,768,4+33,56
Regresi elevasi + residu IDW42,233,3+9,66
Ilustrasi 5.1: Galat metode interpolasi
Diagram batang RMSE tujuh metode pada enam stasiun uji, metode regresi elevasi plus residu IDW paling kecil

Bacalah tabel ini dengan hati-hati:

  • TIN tidak adil dibandingkan. Ia hanya mencakup 2 atau 3 dari 6 titik uji. RMSE-nya dihitung dari titik yang sedikit itu, jadi tidak ikut dalam peringkat.
  • IDW punya bias positif (menebak terlalu tinggi). Penyebabnya belum diperiksa. Dugaan penulis: IDW tidak bisa menebak di luar rentang nilai titik latih, padahal hujan naik mengikuti elevasi. [kemungkinan]
  • Metode regresi menang karena datanya memang dibuat naik bersama elevasi. Ini data sintetis. Pada data nyata, hubungan itu harus dibuktikan dulu.

Validasi silang tinggalkan-satu pada 24 titik latih untuk IDW (RMSE menurut pangkat p):

p0,51234
RMSE169,0149,4115,396,090,1

Pada data ini, p yang lebih besar memberi galat lebih kecil. Pangkat 2 adalah nilai awal umum, bukan jaminan terbaik.

Ilustrasi 5.2: Peta curah hujan hasil interpolasi
Peta sintetis curah hujan tahunan dari 2200 sampai sekitar 3000 milimeter dengan stasiun latih berupa lingkaran putih dan stasiun uji berupa belah ketupat merah

Peta pemenang disimpan sebagai hasil/CH_tahunan.tif (200 x 200 sel, 10 m) dan dipakai di Bab 6.

#Bagian B: ArcGIS Pro

#Bagian B: ArcGIS Pro

Kotak peralatan Spatial Analyst Tools memuat kelompok Interpolation: IDW, Kriging, Spline, Natural Neighbor, dan Trend. Pada IDW, dokumentasi Esri menyebut parameter Power dengan nilai awal 2, serta pilihan radius pencarian (variabel atau tetap) dan barrier berupa garis patahan. [Kriging, Spline, dan Trend: nama alat dari pengetahuan umum, belum dibaca di dokumentasi; UJI]

  1. Siapkan dua layer dari stasiun: latih dan uji. Pisahkan lewat Select By Attributes (peran = 'latih') dan ekspor. [CEK]
  2. Jalankan IDW pada stasiun latih, ch_tahunan_mm sebagai Z value field. Isi Power. [CEK]
  3. Jalankan Kriging (pilih model Spherical atau Exponential) dan Spline. [CEK]
  4. Ambil nilai di titik uji dengan Extract Values to Points, lalu hitung selisihnya di tabel. [CEK]
  5. Ekstensi Geostatistical Analyst menyediakan wizard dan validasi silang yang lebih lengkap, termasuk kriging bayes empiris. [CEK: ekstensi dan nama alat]

#Bagian C: ArcMap 10.8

#Bagian C: ArcMap 10.8

  1. Spatial Analyst Tools ► Interpolation: IDW, Kriging, Spline, Natural Neighbor, Trend. [CEK]
  2. Extract Values to Points untuk membaca nilai di titik uji. [CEK]
  3. Geostatistical Analyst ► Geostatistical Wizard untuk kriging dengan validasi silang. [CEK]

#Cek paham

  1. Mengapa nilai enam stasiun uji tidak boleh ikut dipakai saat menebak?
  2. TIN memberi nilai hanya di 2 dari 6 titik uji. Apa penyebabnya?
  3. Mengapa metode regresi menang pada latihan ini, dan mengapa itu belum tentu berlaku pada data nyata?

Jawaban:

  1. Bila ikut dipakai, tebakannya "mencontek". Nilai galat jadi terlalu bagus dan tidak menunjukkan kemampuan menebak tempat baru.
  2. TIN hanya menghitung di dalam segitiga yang dibentuk titik latih. Titik uji di luar tepi jaring tidak punya nilai.
  3. Data dibuat sehingga hujan naik bersama elevasi, jadi regresi elevasi cocok. Pada data nyata, hubungan itu harus diperiksa dulu.

#Kesalahan umum

  • Menilai metode dengan titik yang dipakai untuk melatihnya. Hasilnya selalu tampak sempurna (interpolator eksak).
  • Membandingkan RMSE dari jumlah titik yang berbeda. Seperti pada TIN. Laporkan juga jumlah titik.
  • Memilih kolom atribut yang keliru. Periksa nilai terkecil dan terbesar hasil. Hujan 120 sampai 275 jelas bukan milimeter, itu elevasi.
  • Memakai IDW untuk menebak jauh di luar sebaran titik. Hasil tertarik ke rata-rata dan tidak pernah melampaui rentang data.

#Ringkasan dan latihan

Ringkasan: interpolasi menebak nilai di tempat tak terukur. Sisihkan titik uji, bandingkan metode dengan RMSE, laporkan jumlah titik, dan curigai hasil yang terlalu bagus.

Latihan: jalankan IDW p = 2 dan p = 3 pada data latih, lalu hitung sendiri RMSE di enam stasiun uji dengan Sample raster values. Cocokkan dengan tabel A6.

#Tabel perbandingan: interpolasi

MetodeQGISArcGIS ProArcMap 10.8
IDWIDW interpolationIDW [CEK]IDW [CEK]
TINTIN interpolationCreate TIN, TIN to Raster [CEK]sama [CEK]
Splinev.surf.bspline, v.surf.rst (GRASS)Spline [CEK]Spline [CEK]
Krigingtidak ada alat bawaan yang diuji; skrip NumPyKriging, Geostatistical Analyst [CEK]Kriging, Geostatistical Analyst [CEK]
Nilai di titik ujiSample raster valuesExtract Values to Points [CEK]sama [CEK]