MODEL ADITIF TERAMPAT VEKTOR DENGAN KOMPONEN UTAMA UNTUK PENDUGAAN CURAH HUJAN EKSTRIM (STUDI KASUS: INDRAMAYU)
EKA PUTRI NUR UTAMI
SEKOLAH PASCASARJANA INSTITUT PERTANIAN BOGOR BOGOR 2016
PERNYATAAN MENGENAI TESIS DAN SUMBER INFORMASI SERTA PELIMPAHAN HAK CIPTA* Dengan ini saya menyatakan bahwa tesis berjudul “Model Aditif Terampat Vektor dengan Komponen Utama untuk Pendugaan Curah Hujan Ekstrim (Studi Kasus: Indramayu)” adalah benar karya saya dengan arahan dari komisi pembimbing dan belum diajukan dalam bentuk apapun kepada perguruan tinggi manapun. Sumber informasi yang berasal atau dikutip dari karya yang diterbitkan maupun tidak diterbitkan dari penulis lain telah disebutkan dalam teks dan dicantumkan dalam Daftar Pustaka di bagian akhir tesis ini. Dengan ini saya melimpahkan hak cipta dari karya tulis saya kepada Institut Pertanian Bogor. Bogor, Januari 2016
Eka Putri Nur Utami NIM G151130091
* Pelimpahan hak cipta atas karya tulis dari penelitian kerjasama dengan pihak luar IPB harus didasarkan pada perjanjian kerjasama yang terkait
RINGKASAN EKA PUTRI NUR UTAMI. Model Aditif Terampat Vektor dengan Komponen Utama untuk Pendugaan Curah Hujan Ekstrim (Studi Kasus: Indramayu). Dibimbing oleh AJI HAMIM WIGENA dan ANIK DJURAIDAH. Perubahan iklim yang terjadi beberapa tahun terakhir ini mengakibatkan terjadinya curah hujan ekstrim. Curah hujan ekstrim ini berpotensi menimbulkan banjir yang akhirnya dapat mengakibatkan gagal panen. Oleh karena itu pemodelan curah hujan diperlukan untuk meminimumkan dampak yang terjadi. Statistical downscaling (SD) merupakan model statistik yang digunakan untuk menduga curah hujan (berskala lokal) dengan memanfaatkan informasi dari data luaran global circulation model (GCM). Data GCM adalah data hasil simulasi komputer yang memanfaatkan kaidah fisika, kondisi lautan, dan perubahan iklim pada atmosfer bumi yang dapat digunakan untuk menduga unsurunsur iklim saat ini ataupun di masa yang akan datang. Permasalahan yang muncul dalam metode SD adalah penentuan teknik statistika yang tepat untuk memodelkan curah hujan ekstrim. Sebaran nilai ekstrim yang biasa digunakan dalam teori nilai ekstrim adalah generalized extreme value (GEV) dan generalized pareto distribution (GPD). Penggunaan sebaran GEV memiliki kelemahan yaitu akan menghilangkan banyak data amatan. Oleh karena itu diusulkan untuk menggunakan GPD. Pemodelan yang dapat digunakan dalam SD antara lain model aditif terampat (generalized additive model/GAM) yang dapat mengatasi pengaruh nonlinier masing-masing peubah penjelas terhadap peubah respon dengan teknik pemulusan. Kekurangan dari metode ini adalah GAM hanya terbatas pada sebaran-sebaran yang termasuk ke dalam sebaran keluarga eksponensial, sedangkan GPD tidak termasuk ke dalam sebaran keluarga eksponensial. Selain itu, metode GAM hanya dapat membuat satu fungsi hubung dari parameter sebarannya. Permasalahan tersebut, diatasi dengan vector generalized additive model (VGAM). Teknik ini merupakan perluasan dari model aditif terampat karena VGAM yang tidak terbatas pada sebaran keluarga eksponensial. Tujuan dari penelitian ini adalah menduga curah hujan ekstrim menggunakan VGAM berdasarkan sebaran GPD di Kabupaten Indramayu. Data curah hujan di Indramayu digunakan sebagai peubah respon dan data presipitasi GCM sebagai peubah penjelas. Data penelitian dibagi menjadi dua yaitu data tahun 1979-2007 untuk pemodelan dan data tahun 2008 sebagai validasi model. Ukuran grid data GCM yang digunakan yaitu ukuran 6 × 6 dan 8 × 8. Sifat dari data GCM yang berdimensi tinggi dan terdapat multikolinearitas diatasi dengan menggunakan analisis komponen utama. Pemodelan VGAM dilakukan terhadap data yang lebih besar dari nilai ambang dengan menggunakan metode backfitting vector. Pemilihan nilai ambang ini dilakukan dengan menggunakan mean residual life plot (MRLP). Pemulusan data menggunakan pemulusan spline dan pemilihan derajat bebas optimum dilakukan dengan memilih nilai Akaike Information Criterion (AIC) terkecil. Pemodelan dilakukan dengan memodelkan keseluruhan bulan, dan membagi bulan dalam empat kelompok. Kelompok tersebut adalah bulan hujan, peralihan
hujan-kering, kering, dan peralihan kering-hujan. Model yang dibuat adalah model keseluruhan dan model setiap kelompok. Hasil analisis komponen utama menghasilkan dua komponen utama untuk data grid 6 × 6 dan empat komponen utama untuk data grid 8 × 8 . Proporsi keragaman kumulatif yang dihasilkan adalah 96% untuk data ukuran grid 6 × 6 dan 95% untuk grid 8 × 8. Berdasarkan hasil MRLP diperoleh nilai ambang secara keseluruhan adalah 145. Nilai ambang untuk kelompok bulan hujan, peralihan hujan-kering, kering, dan peralihan kering-hujan berturut-turut adalah 145, 100,10 dan 45. Hasil pemodelan VGAM pada tahap awal menunjukkan data grid 6 × 6 menghasilkan nilai root mean square error prediction (RMSEP) yang lebih rendah dibandingkan dengan data grid 8 × 8 . Oleh karena itu pada tahap selanjutnya pemodelan dilakukan dengan menggunakan data grid 6 × 6 . Secara keseluruhan hasil pendugaan curah hujan dengan metode VGAM yang berbasis GPD mampu menghasilkan pola dugaan yang mirip dengan data aktual. Penggunaaan kuantil yang berbeda untuk setiap kelompok bulan menghasilkan dugaan yang lebih baik. Kelompok bulan kering memiliki nilai dugaan yang mendekati aktual pada kuantil rendah yaitu kuantil 5 dan 10, kelompok bulan hujan pada kuantil 75, kelompok bulan peralihan hujan-kering pada kuantil 50, sedangkan pada kelompok bulan peralihan kering-hujan dugaan terbaik yaitu pada kuntil 25. Kebaikan model menunjukkan bahwa model VGAM dengan GPD menghasilkan dugaan yang konsisten meskipun panjang data dugaan ditambah, nilai RMSEP yang dihasilkan tidak jauh berbeda jika dibandingkan dengan melakukan dugaan pada satu tahun.
Kata kunci : generalized pareto distributions, model aditif terampat vektor, statistical downscaling, teori nilai ekstrim
SUMMARY EKA PUTRI NUR UTAMI. Vector Generalized Additive Model with Principal Components to Estimate Extreme Rainfall (Case Study: Indramayu). Supervised by AJI HAMIM WIGENA and ANIK DJURAIDAH. In the last few years, climate change happened and caused extreme rainfall. A high extreme rainfall potentially cause flooded and then impact on agricultural. Therefore modeling rainfall data needed to minimize the impact. Statistical downscaling (SD) is a statistical model to predict rainfall data (local-scale) by using the information of global circulation model (GCM). The GCM data is a computer simulation result of a large number interaction of physics, chemistry, and dynamics of the earth’s atmosphere. GCM can be used to predict climate elements in the future. The problems for SD are to choose the methods for modeling the extreme events. There are two distributions in extreme value theory. They are generalized extreme value (GEV) and generalized pareto distributions (GPD). Modeling based on GEV will elliminate some data, than proposed modeling with generalized pareto distributions (GPD). One of the models used in SD is generalized additive models (GAM). The advantage of GAM is it can be solved the non-linear effects of predictor variables to responses variable by using a smoothing technique. Unfortunately, GAM are restricted to one-parameter distributions within the classical exponential family, whereas GPD is not included of exponential family. So, there are a new method to solve this problems. The methods known as vector generalized models (VGAM). The method of VGAM is claimed to be a larger framework of distribution. VGAM was an extension from generalized additive models which is not restricted for exponential family. The purpose of this research is to predict the rainfall data in Indramayu based on extreme value distribution (GPD) with VGAM technique. This research used monthly rainfall data in Indramayu from1979 to 2008 as independent variable and precipitation from GCM as dependent variables. The data divided into two parts, first is 1979-2007 as data modeling and then 2008 as data validation. There are two size of grid for GCM data, that is 6 × 6 and 8 × 8. Principal compenent analysis used to reduce the dimention of the data and solved the multicollinearity. Data for modeling are the rainfall data that is greater than the threshold value. The threshold values choosed by using the mean residual life plot (MRLP). The method of VGAM are estimated by vector backfitting algorithm with vector smoothing. The optimum degree of freedom from smoothing is choosed by the minimum Akaike Information Criterion (AIC). For the first step, modeling the data for all month and then modeling data for group of month. The groups are rainfall, rainfall-dry transition, dry and transition dry-rainfall. The result of principal component analysis show that there are two components choosed for grid data 6 × 6 and four components for 8 × 8. The cumulative proportion of data variety is 96% for grid data 6 × 6 and 95% for grid data 8 × 8. Based on MRLP, the threshold value for all months is 145 and
the threshold value for group rainy, rainy-dry transition, dry and dry-rainy transition are 145, 100, 10 and 45 respectively. For the model of all month, VGAM result shows that the model from grid data 6 × 6 have the lowest value of RMSEP than the model from grid data 8 × 8. So for the next step, modeling is use with the data from grid 6 × 6. Overall, the result show that predicted value have the same pattern with the actual data. By using different quantile for each group of month can give a better estimation than using a same quantile. Dry group have a good predicted value at quantile 5 and 10, but for the rainy group is at quantile 75. Goodness of fit model showed that VGAM based on GPD can give consistent result for different length data. Keywords : extreme value theory, generalized pareto distribution, statistical downscaling, vector generalized additive model.
© Hak Cipta Milik IPB, Tahun 2016 Hak Cipta Dilindungi Undang-Undang Dilarang mengutip sebagian atau seluruh karya tulis ini tanpa mencantumkan atau menyebutkan sumbernya. Pengutipan hanya untuk kepentingan pendidikan, penelitian, penulisan karya ilmiah, penyusunan laporan, penulisan kritik, atau tinjauan suatu masalah, dan pengutipan tersebut tidak merugikan kepentingan IPB Dilarang mengumumkan dan memperbanyak sebagian atau seluruh karya tulis ini dalam bentuk apapun tanpa izin IPB
MODEL ADITIF TERAMPAT VEKTOR DENGAN KOMPONEN UTAMA UNTUK PENDUGAAN CURAH HUJAN EKSTRIM (STUDI KASUS: INDRAMAYU)
EKA PUTRI NUR UTAMI
Tesis sebagai salah satu syarat untuk memperoleh gelar Magister Sains pada Program Studi Statistika
SEKOLAH PASCASARJANA INSTITUT PERTANIAN BOGOR BOGOR 2016
Penguji Luar Komisi pada Ujian Tesis : Dr Ir Kusman Sadik, MSi
Judul Tesis : Model Aditif Terampat Vektor dengan Komponen Utama untuk Pendugaan Curah Hujan Ekstrim (Studi Kasus: Indramayu) Nama : Eka Putri Nur Utami NIM : G151130091
Disetujui oleh Komisi Pembimbing
Dr Ir Aji Hamim Wigena, MSc Ketua
Dr Ir Anik Djuraidah, MS Anggota
Diketahui oleh
Ketua Program Studi Statistika
Dekan Sekolah Pascasarjana
Dr Ir Kusman Sadik, MSi
Dr Ir Dahrul Syah, MScAgr
Tanggal Ujian : 19 November 2015
Tanggal Lulus :
PRAKATA Puji syukur penulis panjatkan kepada Allah SWT atas limpahan rahmat dan karunia-Nya sehingga penulis dapat menyelesaikan tesis yang berjudul “Model Aditif Terampat Vektor dengan Komponen Utama untuk Pendugaan Curah Hujan Ekstrim (Studi Kasus: Indramayu)”. Sesungguhnya penyelesaian tugas akhir ini semuanya adalah berkat kemudahan dan pertolongan dari-Nya. Terima kasih penulis ucapkan kepada Bapak Dr.Ir.Aji Hamim Wigena, MSc dan Ibu Dr.Ir.Anik Djuraidah,MS selaku pembimbing, yang telah sabar membimbing penulis dalam penyusunan tesis ini. Penulis juga mengucapkan ucapan terima kasih kepada Bapak Dr.Ir.Kusman Sadik MSi sebagai dosen penguji dan Ibu Dr. Ir. Indahwati,MSi selaku dosen perwakilan prodi statistika yang telah meluangkan waktunya hadir pada sidang tesis saya. Penulis juga mengucapkan ucapan terimakasih yang sebesar-besarnya kepada seluruh dosen departemen Statistika IPB yang telah mendidik dan memberikan ilmunya kepada penulis selama di bangku kuliah hingga berhasil menyelesaikan studi, serta seluruh staf departemen Statistika IPB atas bantuan, pelayanan, dan kerjasamanya selama ini. Terimakasih juga kepada rekan-rekan satu grup riset downscaling, Bapak Agus M Soleh serta Shynde Limar Kinanti dan Dewi Santri atas kerjasamanya selama ini dalam penelitian downscaling sehingga penyusunan tesis ini dapat selesai dengan baik. Ucapan terima kasih yang tulus dan penghargaan yang tak terhingga juga penulis ucapkan kepada kedua orangtua, Bapak Isman dan Ibu Nuriyah yang telah memberikan kasih sayangnya kepada penulis dan mendidik penulis dengan penuh kasih sayang. Selain itu kepada kedua adikku tersayang Lena Dwi Prestiyani dan Widiya Tri Astuti atas doa dan semangatnya. Terakhir tak lupa penulis juga mengucapkan terima kasih kepada seluruh mahasiswa Pascasarjana departemen Statistika atas kerjasama dan kebersamaannya selama menghadapi masa-masa perkuliahan baik itu senang ataupun susah. Terimakasih karena telah menciptakan kenangan yang takkan terlupakan. Semoga tesis ini dapat bermanfaat bagi semua pihak yang membutuhkan.
Bogor, Januari 2016
Eka Putri Nur Utami
DAFTAR ISI DAFTAR TABEL
vi
DAFTAR GAMBAR
vi
DAFTAR LAMPIRAN
vi
1 PENDAHULUAN Latar Belakang Tujuan Penelitian
1 1 2
2 TINJAUAN PUSTAKA Statistical Downscaling (SD) Teori Nilai Ekstrim Model Aditif Terampat Vektor (VGAM) Model Aditif Terampat Vektor untuk Nilai Ekstrim
3 3 4 6 10
3 METODE PENELITIAN Data Analisis Data
11 11 11
4 HASIL DAN PEMBAHASAN Eksplorasi Data Penentuan Nilai Ambang Pemodelan GCM
13 13 13 15
5 SIMPULAN DAN SARAN Simpulan Saran
25 25 25
DAFTAR PUSTAKA
26
LAMPIRAN
27
RIWAYAT HIDUP
38
DAFTAR TABEL 1. Nilai akar ciri dan proporsi keragaman dari kumulatif analisis komponen utama 2. Nilai AIC dari model VGAM pada beberapa derajat bebas 3. Nilai RMSEP untuk model awal 4. Hasil uji sebaran GPD dengan uji Kolmogorov-Smirnov 5. Nilai AIC dari model VGAM pada beberapa derajat bebas pemulus setiap komponen utama 6. Penduga parameter VGAM masing-masing kelompok bulan 7. Nilai RMSEP data dugaan pada berbagai kuantil 8. Nilai RMSEP pada panjang data dugaan yang berbeda 9. Nilai korelasi data dugaan dengan data aktual pada pada panjang data dugaan yang berbeda
15 16 18 19 19 21 23 23 24
DAFTAR GAMBAR 1 . Ilustrasi proses statistical downscaling 2 . Pola curah hujan Kabupaten Indramayu tahun 1979-2008 3. Grafik MRLP data curah hujan 4. Grafik hubungan pendugaan parameter GPD dengan nilai ambang data curah hujan rata-rata 5. Pengepasan model VGAM pada data berukuran 6 × 6 6. Pengepasan model VGAM pada data berukuran 8 × 8 7. MRLP curah hujan masing-masing kelompok bulan 8. Pengepasan model VGAM pada model kelompok curah hujan 9. Grafik hasil validasi silang model terhadap data curah hujan tahun 2008 pada kuantil atas 10. Grafik hasil validasi silang model terhadap data curah hujan tahun 2008 pada kuantil bawah sampai atas
3 13 14 14 17 17 18 20 22 22
DAFTAR LAMPIRAN 1. 2. 3. 4.
Plot diagnostik kecocokan data curah hujan dengan GPD Plot curah hujan ekstrim dengan skor komponen utama terpilih Nilai dugaan curah hujan tahun 2008 pada model awal Grafik hubungan pendugaan parameter dengan nilai ambang pada data curah hujan setiap kelompok bulan 5. Plot diagnostik untuk nilai ambang terpilih setiap kelomok bulan 6. Plot pemulusan pada data curah hujan rata-rata berdasarkan kelompok bulan 7. Nilai aktual dan dugaan curah hujan pada model setiap kelompok bulan pada tahun 2008
27 28 29 32 33 35 37
1 PENDAHULUAN Latar Belakang Indonesia merupakan negara kepulauan yang beriklim tropis. Perubahan iklim yang terjadi beberapa tahun terakhir ini mengakibatkan terjadinya curah hujan ekstrim. Curah hujan ekstrim merupakan kondisi curah hujan yang sangat tinggi atau sangat rendah. Curah hujan ekstrim tinggi berpotensi menimbulkan banjir yang akhirnya mengakibatkan gagal panen. Oleh karena itu diperlukan pendugaan melalui suatu pemodelan untuk mengantisipasi dampak yang akan terjadi. Penelitian mengenai pemodelan curah hujan menjadi salah satu hal penting untuk memprediksi curah hujan. Salah satu metode pemodelan curah hujan adalah dengan memanfaatan data Global Circulation Model (GCM). Menurut Wigena (2006) GCM merupakan sumber informasi primer dan dapat digunakan untuk memprediksi iklim dan cuaca secara numerik. GCM merupakan model peramalan curah hujan pada skala lokal dengan mempertimbangkan berbagai informasi atmosfer global yang disimulasikan. Teknik statistika yang menyatakan hubungan fungsional antara peubahpeubah luaran GCM dengan peubah respon lokal dikenal dengan statistical downscaling (SD). Peubah respon lokal yang digunakan untuk pemodelan curah hujan adalah curah hujan daerah setempat dan data GCM sebagai peubah skala global. Pemodelan SD dengan menggunakan data GCM menyertakan peubah berdimensi besar yang saling berkorelasi, sehingga untuk mengatasi hal ini pereduksian dimensi umumnya dilakukan dengan menggunakan analisis komponen utama. Permasalahan yang muncul dalam metode SD adalah menentukan teknik statistika yang tepat. Berbagai penelitian mengenai metode SD sudah dikaji dan menunjukkan gambaran perkembangan pemodelan SD banyak mengarah ke teknik statistika berbasis regresi parametrik, semi parametrik atau statistika nonparametrik. Penelitian mengenai pemodelan semi parametrik dilakukan oleh Rizki (2014) dan memperoleh kesimpulan bahwa pola hubungan fungsional curah hujan dengan komponen utama memiliki hubungan yang berupa gabungan parametrik dan nonparametrik yaitu semiparametrik kubik. Penelitian mengenai curah hujan yang berbasis nilai ekstrim telah dilakukan sebelumnya. Secara umum data curah hujan tidak mengikuti sebaran normal, sehingga pemodelan baku tidak dapat dilakukan. Beberapa penelitian menyatakan bahwa curah hujan memiliki sebaran gamma. Sedangkan untuk curah hujan ekstrim, banyak pendekatan dilakukan dari teori nilai ekstrim. Penelitian curah hujan ekstrim sudah dilakukan antara lain oleh Mondiana (2012) yang menggunakan regresi kuantil, dan Handayani (2014) yang melakukan penelitian mengenai SD untuk curah hujan berdasarkan sebaran nilai ekstrim Generalized Extreme Value (GEV). Menurut Mallor (2009) pemilihan nilai ekstrim dengan metode block maxima yang mengikuti sebaran GEV akan mengakibatkan hilangnya banyak data, sehingga baik digunakan jika data yang dimiliki sangat banyak. Alternatif lain dalam pemilihan nilai ekstrim adalah dengan metode peaks over threshold (POT). Nilai ekstrim yang diperoleh dari metode ini akan mengikuti sebaran generalized pareto distribution (GPD).
2 Contoh penelitian curah hujan yang berdasarkan sebaran GPD adalah Li et al. (2005) yang mengidentifikasi curah hujan ekstrim di wilayah Australia. Penelitian tersebut tidak menggunakan peubah penjelas dalam pemodelan curah hujan ekstrim, hanya menduga parameter-parameter GPD. Oleh karena itu penelitian kali ini akan melakukan pemodelan curah hujan ekstrim yang mengikuti sebaran GPD dengan mengikutsertakan pengaruh peubah luaran GCM. Model aditif terampat (generalized additive model/GAM) adalah suatu metode nonparametrik yang dapat digunakan untuk memodelkan data GCM. GAM merupakan generalisasi dari model linier terampat (generalized linear model /GLM). Metode GAM mengganti bentuk linier hubungan Y dan X sebagai suatu hubungan fungsional yang diharapkan mampu mengatasi pengaruh nonlinier terhadap peubah y (Hastie dan Tibshirani, 1986) Kelemahan dari metode GLM dan GAM adalah peubah respon terbatas pada sebaran keluarga eksponensial, sedangkan GPD tidak termasuk kedalam keluarga eksponensial. Yee dan Wild (1996) memperkenalkan teknik model aditif terampat vektor (vector generalized aditif model/VGAM) yang merupakan perluasan dari teknik GAM, teknik ini memperluas teknik GAM ataupun GLM dalam tiga cara utama: 1. Sebaran y tidak terbatas hanya untuk keluarga eksponensial 2. Mampu mengatasi multi respon y dengan prediktor linear atau prediktor aditif 3. ηj tidak harus merupakan fungsi dari mean μ, akan tetapidapat sebagai fungsi dari parameter lainnya, yaitu ηj = g j (θj ) untuk setiap parameter θj Sehingga VGAM dapat mengatasi banyak sebaran dan model yang memungkinkan. Yee dan Wild (1996) telah mengembangkan metode VGLM/VGAM yang menjadi dasar berkembangnya metode tersebut. Salah satunya dilakukan oleh Yee dan Stephenson (2007) mengenai VGLM/VGAM untuk nilai ekstrim. Penelitian tersebut menggunakan sebaran dari teori nilai ekstrim yaitu GEV dan GPD sebagai dasar dalam pemodelanVGAM. Sehingga dalam penelitian kali ini metode tersebut akan dikaji dan diikuti namun dengan menggunakan data presipitasi dari data luaran GCM sebagai bagian dari statistical downscaling. Teknik VGAM diharapkan mampu menduga curah hujan ekstrim yang terjadi di wilayah Indramayu. Tujuan Penelitian Tujuan penelitian ini adalah menduga curah hujan ekstrim menggunakan model aditif terampat vektor berdasarkan sebaran GPD di Kabupaten Indramayu.
3
2 TINJAUAN PUSTAKA Statistical Downscalling Data GCM adalah data hasil simulasi komputer yang memanfaatkan kaidah fisika, kondisi lautan, dan perubahan iklim pada atmosfer bumi yang dapat digunakan untuk menduga unsur-unsur iklim saat ini ataupun di masa yang akan datang. Data yang dihasilkan GCM berupa unsur-unsur iklim dalam bentuk grid (persegi) . Setiap grid berukuran 2.5° × 2.5° dan memiliki nilai unsur-unsur iklim masing-masing. Namun, data GCM masih berskala global dan tidak cocok digunakan untuk menduga unsur-unsur iklim pada skala yang lebih kecil. Resolusi data GCM yang rendah tidak cocok digunakan sebagai penduga data berskala lokal, tetapi masih memungkinkan sebagai penduga jika menggunakan teknik downscaling. Downscaling adalah suatu teknik transformasi hasil simulasi GCM pada skala besar ke dalam skala yang lebih kecil. Statistical downscaling (SD) adalah proses downscaling yang memanfaatkan ilmu statistika untuk menduga data berskala lokal. SD didasarkan pada pandangan bahwa iklim regional dikondisikan oleh dua faktor yaitu iklim skala global dan kondisi fisiografi lokal yaitu curah hujan di daratan (Wilby et al. 2004). Gambar 1 menunjukkan ilustrasi dari SD. Output dari data GCM adalah memberikan masukan untuk pemodelan statistik terhadap data daerah lokal yang berhubungan. SD dapat dimulai denganmenghubungkan variabel iklim skala global (X) dengan variabel-variabel lokal (Y). Luaran skala besar dari simulasi GCM dimasukkan sebagai input bagi model statistik tersebut untuk memperkirakan karakteristik iklim lokal atau regional yang bersangkutan. Bentuk umum dari model SD adalah: 𝐲 = f(𝐗) dengan , 𝐲 : peubah iklim lokal (curah hujan) 𝐗 : peubah GCM (presipitasi) t : banyaknya waktu (bulanan) g : banyaknya grid domain GCM
Gambar 1 Ilustrasi proses statistical downscaling Sumber:Wilby RL & Dawson (2007)
4 Teori Nilai Ekstrim Teori nilai ekstrim (extreme value theory, EVT) merupakan salah satu cara yang dapat digunakan untuk memodelkan kejadian-kejadian ekstrim. Hal yang membedakan antara EVT dengan teknik statisktika lain adalah EVT bertujuan untuk menghitung pola stokastik dari amatan yang sangat besar atau sangat kecil dibandingkan dengan amatan lain. Analisis nilai ekstrim biasanya memerlukan pendugaan dari sebaran amatan yang lebih ekstrim dari yang sudah diamati (Coles 2001). Model EVT didasarkan pada kejadian-kejadian ekstrim. Misalkan peubah acak Y1 , Y2 … , Yn yang bebas stokastik identik dengan sebaran F maka EVT didasarkan pada nilai Mn yaitu Mn = maks{Y1 , … , Yn } Sebaran Mn pada teori dapat diturunkan dengan mudah jika sebaran Yi diketahui. Pemilihan kejadian ekstrim terdiri dari dua cara, pertama memilih kejadian ekstrim dalam suatu blok atau yang dikenal dengan block maxima sehingga data mengikuti sebaran generalized extreme value distribution (GEV). Fungsi sebaran GEV adalah sebagai berikut : 1 y − μ −ξ F(y) = exp {− [1 + ξ ( )] } σ y−µ
dengan x ∶ 1 + ξ ( σ ) > 0, -∞ < µ < ∞, σ > 0 dan -∞ < ξ < ∞, dalam model ini ξ, σ, dan µ berturut-turut adalah parameter bentuk, skala, dan lokasi. Kelemahan dari metode blok maxima adalah keharusan membagi data ke dalam ukuran blok yang sama. Selain itu pemilihan nilai maksimum dari setiap blok akan mengakibatkan hilangnya banyak data amatan. Oleh karena itu metode kedua yang digunakan untuk menentukan kejadian ekstrim adalah dengan menentukan nilai ambang u dan dikenal dengan metode pemilihan nilai ambang (peaks over threshold/POT) . Metode POT ini akan mengakibatkan data mengikuti sebaran generalized pareto distribution (GPD). Kejadian-kejadian ekstrim didefinisikan sebagai data Y yang lebih besar dari u. Dalam terminologi kejadian ekstrim Y − u disebut excesses (Yee 2007). Fungsi sebaran (Y − u) untuk Y > u didekati dengan fungsi berikut yang dikenal dengan sebaran pareto terampat (generalized pareto distribution/GPD). 1
G(y) = {
1 − (1 1
σ
ξy − ξ + σ) −y
− exp ( σ )
,ξ ≠ 0 ,ξ = 0
dan fungsi kepekatan GPD: 1
1
g (y) = {σ
1
σ
(1 +
ξy − ξ −1 ) σ −y
exp ( σ )
,ξ ≠ 0 ,ξ = 0
(1)
5 dalam model ini ξ adalah parameter bentuk, dan σ adalah parameter skala. Sebaran GPD akan menjadi bentuk pareto jika nilai ξ > 0 . Hubungan antara sebaran GPD dengan pareto dapat diturunkan sebagai berikut: 𝜎 jika ξ > 0 maka 𝜇 = 𝜉 g (y) = = = =
1
(1 +
𝜎 1 𝜎 1 𝜎 1
(1 + (1 + 𝜉
𝜎
(𝜎)
1 𝜉
𝜎 𝜉 (𝜉 )
𝜉(𝑋−𝜇) 𝜎 𝜎 𝜉
−1 −1 𝜉
)
−1 −1 𝜉
𝜉(𝑋− ) 𝜎 𝜉𝑋 𝜎
)
− 1)
−1 −1 𝜉
𝑋
−1 −1 𝜉
−1 −1 𝜉
1
=
𝑋
1 𝑝
( +1) 1
𝜎 𝜉 ( ) 𝜉 𝜉 1
=
1
( +1)
𝑋𝑝 Sehingga fungsi kepekatan pareto adalah: 𝛼𝛽 𝛼 h (y) = (𝛼+1) 𝑋 1 𝜎 dengan 𝜉 = 𝛼 dan 𝜉 merupakan parameter dari sebaran pareto. Penentuan nilai ambang u pada GPD dapat ditentukan dengan berbagai cara salah satunya adalah dengan menggunakan mean residual life plot (MRLP). Pemilihan nilai ambang dengan menggunakan MRLP dianalogikan seperti menentukan ukuran blok yang akan dipilih sehingga ada keseimbangan antara bias dan keragaman. Nilai ambang yang terlalu rendah mengakibatkan pendekatan dari model yang gagal atau parameter yang diduga bias meskipun ragamnya kecil sedangkan nilai ambang yang terlalu tinggi akan mengakibatkan jumlah nilai excesses yang diamati sedikit sehingga mengakibatkan keragaman yang tinggi. Prosedur untuk menentukan nilai ambang ini dapat dilakukan dengan cara membentuk MRLP. Titik plot (X, Y) dalam MRLP ditentukan berdasarkan: {(u ,
1 nu
u ∑ni=1 (y(i) − u))
∶ u < ymaks }
dengan y(1) , … . , y(nu) adalah observasi sebanyak nu yang melebihi nilai u dan ymaks adalah nilai terbesar dari Yi. Interpretasi dari MRLP tidak selalu mudah dan sedikit subjektif. Nilai ambang batas yang dipilih adalah nilai yang berada diatas nilai dimana plot mendekati linear terhadap u ( Mallor et al. 2009). Pendugaan parameter dari sebaran GPD dapat dilakukan dengan menggunakan metode kemungkinan maksimum. Jika terdapat y1 , … , yk adalah k
6 banyaknya pengamatan yang melebihi nilai ambang u. Untuk ξ ≠ 0 fungsi logkemungkinan didapatkan dari persamaan (1) menjadi 1 ξy ℓ ( σ, ξ ) = −k log σ − ( 1 + ξ ) ∑ki=1 log( 1 + σi ) 𝜉𝑦
dengan syarat ( 1 + 𝜎 𝑖 ) > 0 untuk i= 1, ...,k . Maksimisasi dari fungsi logkemungkinan tersebut dibutuhkan teknik analisis numerik, dan salah satu analisis numerik yang digunakan adalah metode Newton Raphson. Model Aditif Terampat Vektor Regresi linear klasik memerlukan asumsi bahwa peubah respon mengikuti sebaran Normal. Kenyataannya, banyak ditemukan peubah respon yang tidak mengikuti sebaran Normal. Oleh karena itu dikembangkan suatu metode yang dapat mengatasi peubah respon yang tidak menyebar normal yaitu model linier terampat (generalized linear model/GLM). Regresi linear, regresi Poisson, regresi logistik, analisis probit dan beberapa model lainnya dapat diperlakukan sebagai kasus dari GLM. Peningkatan yang utama dari GLM dibandingkan regresi linear adalah kemampuan untuk mengatasi kelas sebaran peubah respon yang lebih besar tetapi masih termasuk ke dalam keluarga eksponensial. Terdapat tiga komponen utama dalam GLM yaitu (McCullagh dan Nelder 1983): 1. Komponen acak: komponen acak adalah peubah respon Y. Dalam model linier terampat peubah respon diasumsikan mengikuti sebaran keluarga eksponensial. 2. Komponen sistematik: komponen sistematik adalah kombinasi linier dari x1 , x2 , … , xp . Sehingga dapat dituliskan sebagai berikut: p η = ∑i=1 βi Xi η disebut juga sebagai penduga linier atau linear predictor dan βi adalah konstanta. 3. Fungsi hubung: yaitu fungsi yang menghubungkan antara komponen acak dengan komponen sistematik. Misalkan E (y) = μ selanjutnya dapat dibuat hubungan p g(μ) = η = ∑i=1 βi Xi Fungsi hubung η yang berpola linier memiliki keterbatasan yaitu tidak mampu mengatasi hubungan nonlinier antara peubah respon dengan peubah bebas. Oleh karena itu dikembangkan model aditif terampat atau generalized additive models (GAM) sebagai model semi-parametrik yang dapat mengatasi pengaruh tidak linier peubah respon. GAM merupakan bentuk aditif dari model GLM dengan mengubah bentuk linier dari ∑(βi xi ) menjadi penjumlahan hubungan fungsional fj (. ) antara peubah bebas dengan peubah respon yang dimodelkan secara nonparametrik dengan suatu fungsi pemulus (Hastie dan Tibshirani 1987). Bentuk umum model GAM adalah sebagai berikut: η = β 0 + ∑m j=1 fj (X k ) Batasan dalam metode GLM/GAM yaitu metode ini hanya dapat memodelkan satu parameter dari sebarannya. Contohnya pada model linier, peubah respon Y diasumsikan mengikuti sebaran Normal (𝜇, 𝜎 2 ) dengan fungs
7 hubung η1 = μ = 𝛃Ti 𝑥. Teori standar GLM/GAM memperlakukan parameter σ sebagai skala. Yee dan Wild (1996) mengusulkan suatu metode yang lebih luas yang mampu memodelkan lebih dari satu parameter dan tidak terbatas hanya pada keluarga eksponensial. Metode ini dikenal dengan vector generalized linear model (VGLM) dan vector generalized additive model (VGAM). VGLM/VGAM merupakan perluasan dari model GLM/GAM dengan mengikutsertakan lebih dari satu prediktor aditif ( η ) (Mackenzie dan Yee 2002). Banykanya prediktor aditif ( η ) yang dapat diikutseratkan maksimal sebanyak jumlah parameter dari sebaran yang digunakan. Sehingga VGLM/VGAM diklaim berpotensial untuk menangani berbagai model dari bermacam-macam sebaran. Bentuk umum dari VGLM adalah sebagai berikut: f(𝐲|𝐱; 𝐁) = h(𝐲, η1 , … , ηm )
(2)
m adalah banyaknya parameter dalam sebaran yang digunakan. T dengan 𝐁 = (𝛃1, … , 𝛃Tm, )T dan 𝛃Ti = ( β(1)1 … β(1)p ) dan p
η(j) = η(j) (x) = β0j + β1(j) x1 + ⋯ + βj(p) xp = β0j + ∑k=1 βk(j) xk (3) Persamaan (3) dalam bentuk umum VGAM berubah menjadi p ηj (𝐱) = β(j)0 + f(j)1 x1 + ⋯ + f(j)p xp = β(j)0 + ∑k=1 f(j)k (xk )
(4)
Model dalam persamaan (2) , dapat dibuat menjadi bentuk log-kemungkinan ℓ(𝛃) = ∑ni=1 𝓵i {η1 (𝐱1 ), … , ηm (𝐱 i )} dengan η(j) = η(j) (𝐱) = 𝛃Tj 𝐱 k . Nilai 𝐔(β) = ϑℓ⁄ϑβ menunjukkan skor vektor −ϑ2 ℓ(𝛃) ⁄ untuk model dan ℐ(𝛃) = adalah matriks informasi. Algoritma ϑ𝛃 ϑ𝛃T Newton-Raphson untuk memaksimalkan log-kemungkinan adalah 𝛃(a) = 𝛃(a−1) + ℐ(𝛃(a−1) )−1 𝐔 (𝛃(a−1) ) Untuk setiap iterasi a. dapat dituliskan dalam bentuk iterasi kuadrat terkecil terboboti yang pertama kali diperkenalkan oleh Green (1984): (𝑎−1)
𝛃(a) = (∑ni=1 𝐗 Ti 𝐖𝑖
−1
(𝑎−1)
𝐗 i ) (∑𝐧𝐢=𝟏 𝐗 𝐓𝐢 𝐖𝑖
𝐳i )
dengan 𝐖i adalah matriks berukuran m × m dengan element (j, k) ϑ2 ℓi (𝐖i )jk = − ϑηj ϑηk Persamaan (4) merupakan penjumlahan fungsi pemulusan dari peubah bebas yang diduga secara simultan. Proses pemulusan dalam model aditif terampat merupakah langkah penting dalam proses pendugaan. Pada metode GAM, pemulusan dilakukan dengan pemulusan spline. Oleh karena itu pada metode VGAM pemulusan dilakukan dengan menggunakan pemulusan vektor spline. Fungsi pemulusan vektor f(x) yang dituliskan (f1 (x), … , fM (x) )T dapat diduga dengan meminimumkan jumlah kuadrat terpenalti, yaitu:
8 m ′′ 2 ∑ni=1{𝐲i − 𝐟 (xi ) }T ∑−1 i {𝐲i − 𝐟 (𝐱 i )} + ∑j=1 λj ∫ fj (t) dt (5)
Pada persamaan (5), Σi adalah matriks ragam peragam galat yang bersifat definit positif. Setiap komponen dari fungsi fj memiliki parameter non-negatif pemulusan λj . Seperti halnya pada spline, persamaan (5) dapat diubah ke dalam suatu nilai dari vektor fungsi f pada observasi nilai x. Jika 𝑥1 < 𝑥2 < ⋯ < 𝑥𝑛 , T 𝐲 = (𝐲1T , … , 𝐲nT )T , 𝐟 = (f1 (x1 ), … , fm (x1 ), … , f1 (xn ), … , fm (xn )) dan 𝚺 = diag (𝚺1 , … , 𝚺n ). Maka kriteria jumlah kuadrat terpenalti pada persamaan (5) sama dengan (𝐲 − 𝐟)T 𝚺 −1 (𝐲 − 𝐟) + 𝐟 T 𝐊𝐟 dan 𝐊 = (𝐐𝐓 −1 𝐐T ) ⨂ diag (λ1 , … , λM ) . Elemen-elemen dari matriks T dan Q adalah sebagai berikut: (𝐐)𝑖𝑖 = h−1 i −1 (𝐐)𝑖+1,𝑖 = −(h−1 i + hi+1
(𝐐)i+2,i = −(h−1 i+1 ) (𝐓)𝑖𝑖 = ( hi + hi+1 )⁄3 (𝐓)𝑖,𝑖−1 = (𝐓)𝑖,𝑖+1 = hi ⁄6 hi = xi+1 − xi untuk i= 1,…., n-1 Pendugaan VGLM dan VGAM Pendugaan VGLM dilakukan dengan menggunakan metode iteratively reweighted least square dengan algoritma sebagai berikut: 1. Bentuk matriks 𝐗 VLM dari 𝐗 LM dan 𝐇1 , … , 𝐇p 𝐗 LM ≡ 𝐗 = (𝒙𝟏 , … , 𝒙𝒏 )𝑻 𝐗 VLM = 𝐗 LM ⊗ 𝐈M Sedangkan 𝐇p adalah matriks kendala dengan elemen 0 dan 1. Penentuan nilai ini berdasarkan fungsi hubung yang ditetapkan. (0)
2. Inisialisasi : a = 1 dan tentukan 𝛈i dari 𝛃(0) atau 𝛍(0) (0)
3. Jika diperlukan hitung 𝛍i , ℓ(0) , definisikan 𝛃(0) ℓ(𝛃) = ∑ni=1 𝓵i {η1 (𝐱1 ), … , ηM (𝐱 i )} (0) (0) 4. Hitung turunan pertama dan kedua 𝒖𝑖 dan 𝐖i (𝐮𝐢 )𝐣 = ϑℓ⁄ϑη j ϑ 2 ℓi (𝐖i )jk = − ϑηj ϑηk (0)T
5. Hitung 𝐳 (0) = (𝐳1 (0)
𝐳i
(0)
= 𝛈i
(0)T T
, … , 𝐳n
(0)
−1
) dengan fungsi:
(0)
+ (𝐖i ) 𝐮i
6. Setelah itu (a). Regresikan 𝐳 (a−1) terhadap 𝐗 VLM dengan bobot
9 (a−1)
W (a−1) = diag (w1 W1
(a−1)
, … , wn Wn
)
(a)
untuk menghasilkan 𝛃 ∗∗(a−1) 𝐳 ∗∗(a−1) = 𝐗 VLM 𝛃(a) + ε∗∗(a−1) 𝐗 VLM ≡ 𝐗 LM ⊗ 𝐈M = (𝐗1T , … , 𝐗 Tn )T (b) Tentukan η(a) = 𝐗 VLM 𝛃(a) (a) (c). Hitung 𝛍i , 𝓵(a) dari 𝛈(a) (d). Uji kekonvergenan, konvergen jika 𝓵(a) − 𝓵(a−1) < ε (e). Jika tidak konvergen dan a < iterasi maksimum , maka: (a−1) (a−1) a = a + 1 kemudian hitung 𝐮i , 𝐖i dan (a−1)
𝐳i
(a−1)
= 𝛈i
(a−1)
+ (𝐖i
(a−1)
) −1 𝐮i
7. Simpan a, 𝛈(a) , 𝛍(a) , 𝛃(a) sebagai penduga akhir Seperti halnya GAM, metode VGAM dilakukan dengan menggunakan algoritma backfitting. Pada prosedur lebih dari satu kovariate, algoritma backfitting dapat diaplikasikan terhadap peubah bebas terkoreksi zi dengan menggunakan vektor pemulus. Istilah ini disebut sebagai algoritma backfitting vektor (Yee & Wild 1996). Metode pendugaan VGAM menggunakan algoritma vektor backfitting yang diaplikasikan terhadap zi dapat dituliskan sebagai berikut: Inisialisasi: (0) 𝛃0 = rata − rata {𝐳1 , . . . , 𝐳n } dan 𝐟k = 0 , k= 1, . . .,p 𝐳i adalah peubah bebas terkoreksi, 𝐳i = 𝐗 i 𝛃(a+1) + 𝐖i−1 𝐝i di adalah elemen ke j dari (𝐝i )j = ϑℓi ⁄ϑηj dan (𝐖i )jk = −
ϑ 2 ℓi
ϑηj ϑηk
Iterasi untuk β=1,2… (a) Iterasikan untuk k=1,….p (i) Hitunglah fungsi vektor 𝐟k∗ sebagai fungsi vektor bobot penghalus (k) dari observasi (xik , 𝐳i ) dan (bobot Wi), i=1, . . ., n ′′ 2 ∑ni=1{𝐳i − 𝐟k (xik )}T 𝐖i {𝐳i − 𝐟k (xik )} + ∑m j=1 λ(j)k ∫ f(j)k (t) dt T
dengan 𝐟k (xik ) = (𝐟(1)k (xik ) … … . , 𝐟M (xik )) . (ii) Duga intersep β0 = β0 + rata − rata {𝐟k∗ (𝐱ik ), … . 𝐟k∗ (𝐱 nk )} (b) (iii) Duga 𝐟k (xik ) = 𝐟k∗ (xik ) − rata − rata {𝐟k∗ (xik ), … . 𝐟k∗ (xnk )} (b) Berhenti jika perbedaan antara 𝐟 (b−1) dengan 𝐟 (b) kecil. (k)
Pada iterasi diatas 𝐳i adalah penyesuaian dari zi sehingga efek dari semua peragam kecuali ke-k dihilangkan. Penyesuaian β0 dibutuhkan untuk memudahkan identifikasi. (k)
𝐳i
(b−1) (𝐱 il ) = 𝐳i − β0 − ∑l≠k 𝐟l
10 Proses pemulusan dalam model aditif terampat merupakah langkah penting dalam proses pendugaan. Derajat bebas dari fungsi pemulus mengukur banyaknya pemulusan yang dilakukan. Berdasarkan Yee (2015), derajat bebas didefinisikan sebagai teras dari matriks A, df = teras (𝐀) dan A didefinisikan sebagai berikut: 𝐀(λ) = (𝐈nM + ∑𝐊)−1 Matriks K adalah matriks penalti yang sudah didefinisikan sebelumnya dan ∑ merupakan matriks ragam peragam galat. Besarnya pemulusan yang dapat dilakukan adalah 2 sampai n. Semakin besar derajat bebas artinya semakin banyak pemulusan yang dilakukan.
VGAM untuk Nilai Ekstrim Analisis VGAM dengan sebaran GPD akan memiliki dua penduga sesuai dengan jumlah parameter dalam fungsi sebaran. Model VGAM GPD dapat dituliskan sebagai berikut (Yee & Stephenson, 2007) log σ ( x) = η1 = β(1) + f1 (x1 ) + f2 (x2 ) + ⋯ + fk (xk ) (6) log (ξ +
1 2
) = η2 = β2
(7)
Nilai dugaan dari VGAM GPD pada persentil-p adalah sebagai berikut: σ yp = μ + ξ [(1 − p)−ξ − 1] dengan nilai σ dan ξ yang diperoleh dari persamaan (6) dan (7), sedangkan p adalah besarnya persentil.
11
3 METODE PENELITIAN Data Data yang digunakan adalah data curah hujan bulanan Kabupaten Indramayu dari tahun 1979 sampai dengan 2008 sebagai peubah respon dan data GCM presipitasi sebagai peubah bebas. Data luaran GCM yang digunakan berasal dari Climate Model Intercomparison Project (CIMP5) dan diperoleh dari http://climexp.knmi.nl. Keterbatasan penelitian ini adalah data curah hujan lokal yang digunakan hanya sampai tahun 2008. Hal ini dikarenakan adanya program dari BMKG pada tahun tersebut sehingga data sampai tahun 2008 dapat dipertanggungjawabkan. Selain itu alasan lainnya adalah dikarenakan agar dapat membandingkan dengan penelitian-penelitian yang sudah dilakukan sebelumnya. Berdasarkan Wigena (2006) ukuran grid yang digunakan adalah 8 × 8. Penelitian kali ini akan mencoba ukuran lain yaitu 6 × 6 , dengan alasan korelasi antar grid pada ukuran grid tersebut lebih tinggi. Posisi wilayah untuk data grid 8 × 8 adalah 1.25 – 18.75 LS sampai 98.75 – 118.75 BT. Sedangkan posisi wilayah untuk data grid 6 × 6 adalah 1.25 – 11.25 LS sampai 101.25 – 113.75. Data yang dihasilkan GCM berupa grid-grid berdasarkan garis lintang dan garis bujur dengan ukuran 2.5° (±300 km). Analisis Data 1. Eksplorasi Data a. Deskripsi data curah hujan dengan membuat diagram kotak garis data curah hujan berdasarkan bulan amatan untuk melihat pola curah hujan yang terjadi. b. Mengidentifikasi data curah hujan ekstrim dengan menentukan nilai ambang melalui MRLP. c. Melakukan uji kesesuaian sebaran dengan menggunakan uji Kolmogorov-Smirnov. Pengujian ini dilakukan dengan menyesuaikan fungsi sebaran empirik curah hujan lokal yang diatas nilai ambang Fn(y) dengan sebaran teoritis tertentu F0(y) dengan hipotesis: H0 : Fn(y) = F0(y) ( data mengikuti sebaran GPD ) H1 : Fn(y) ≠ F0(y) ( data tidak mengikuti sebaran GPD ) Fn(y) adalah fungsi sebaran kumulatif tertentu. Statistik uji Kolmogorov-Smirnov ini adalah D-hitung = sup | Fn(y) - F0(y)| Untuk mendapatkan kesimpulan maka D-hitung dibandingkan dengan tabel ( D1-α ) pada tabel Kolmogorov Smirnov. 2. Data GCM a. Mereduksi data GCM yang bersesuaian dengan analisis komponen utama (AKU), dengan langkah-langkah sebagai berikut: i. Sekumpulan data GCM berukuran p × p, dengan p = 6 dan 8 ii. Menghitung matriks ragam peragam atau matriks korelasi A iii. Menghitung nilai akar ciri dengan menggunakan persamaan |𝐀 − λ𝐈| iv. Menghitung skor komponen utama (KU) dari setiap akar ciri, yaitu 𝐲i = a′i 𝐗 , dengan a′i adalah vektor ciri. v. Menentukan jumlah KU terpilih berdasarkan akar ciri lebih dari 1
12 b. Melihat pola sebaran skor KU sebagai peubah penjelas (x) dengan peubah respon curah hujan (y) 3. Pemodelan Pemodelan dilakukan dengan membagi dua data menjadi data untuk pemodelan sampai tahun 2007 dan data tahun 2008 untuk validasi model. Selanjutnya menghitung rata-rata curah hujan bulanan dari 15 stasiun. Nilai rata-rata inilah yang menjadi peubah respon dalam pemodelan. Tahap pemodelan adalah sebagai berikut: a. Menentukan derajat bebas pemulus pada peubah yang mempunyai hubungan nonlinear dengan pemulusan vektor spline dan nilai akaike information criterion (AIC) AIC = 2k − 2 ln 𝐿 b. Menduga fungsi fj pada model ηj dengan menggunakan algoritma backfitting vector. c. Melakukan prediksi dengan menggunakan data validasi dan model yang sudah dibentuk. Nilai dugaan curah hujan dapat dirumuskan: σ yp = μ + [(1 − p)−ξ − 1] ξ dengan p adalah besarnya kuantil. Penelitian kali ini menggunakan kuantil yang berbeda untuk setiap kelompok bulan, yaitu; i. Bulan hujan: 75, 90,95 ii. Bulan peralihan hujan-kering: 25, 70, 75 iii. Bulan kering: 5, 10, 25 iv. Bulan peralihan kering-hujan: 10, 25, 50 4. Validasi dan Kebaikan model Validasi dan kebaikan model dilakukan untuk memilih model terbaik yang sudah dibuat. Validasi ini dilakukan dengan cara menghitung korelasi dan nilai root mean square error of prediction (RMSEP). Korelasi antara nilai dugaan dan nilai aktul digunakan untuk melihat kesamaan pola dugaan. Semakin mendekati satu artinya data dugaan dan data aktual memiliki hubungan yang kuat dan memiliki pola yang sama. Rumus korelasi adalah sebagai berikut: n n n n yi yˆ i yi yˆ i i 1 i 1 i 1 r n 2 n 2 n 2 n 2 n yi yi n yˆ i yˆ i i 1 i 1 i 1 i 1
Selain itu nilai RMSEP digunakan untuk melihat rata-rata eror antara data dugaan dengan data aktual. Semakin kecil nilai RMSEP artinya semakin kecil selisih antara nilai dugaan dengan nilai aktual artinya, model yang dibangun semakin baik. Rumus dari RMSEP adalah: 1
2 RMSEP = √ ∑ni=1(yi − ŷi ) n
13
4 HASIL DAN PEMBAHASAN Eksplorasi Data Eksplorasi data dilakukan sebagai informasi awal untuk mengetahui karakteristik curah hujan di Indramayu. Deskriptif statistik menunjukkan bahwa rata-rata curah hujan di Kabupaten Indramayu pada tahun 1979 sampai 2008 adalah 123 mm dengan curah hujan minimum yaitu 0 dan curah hujan maksimum mencapai 583 mm. Nilai simpangan baku sebesar 110 mm menggambarkan keragaman data curah hujan yang cukup besar. Gambar 2 menunjukkan gambaran curah hujan di Kabupaten Indramayu yang berpola monsun. Pola monsun adalah pola curah hujan yang membentuk huruf U atau dapat dikatakan memiliki satu puncak musim hujan. Pola curah hujan pada Kabupaten Indramayu memiliki puncak curah hujan pada bulan Januari dan berada pada nilai terendah sekitar bulan Agustus. Rata-rata curah hujan bulan Januari yaitu sebesar 308 mm, sedangkan rata-rata curah hujan bulan Agustus adalah 16 mm. Berdasarkan Gambar 2 dapat dilihat bahwa terjadi kejadian ekstrim pada beberapa bulan yaitu Februari, Juli, Agustus, September, November dan Desember.
Gambar 2 Pola Curah hujan Kabupaten Indramayu tahun 1979 – 2008 Penentuan Nilai Ambang Identifikasi curah hujan ekstrim dilakukan dengan menentukan nilai ambang (u) dengan menggunakan MRLP. Gambar 3 adalah hasil MRLP untuk curah hujan rata-rata di Kabupaten Indramayu. Nilai ambang dipilih dengan cara melihat pola grafik yang membentuk garis lurus atau linier terhadap u. Ketika pola grafik sudah tidak beraturan maka nilai ambang dipilih pada titik awal terjadinya perubahan pola tersebut. Pada Gambar 3 terlihat pola mulai konstan di titik 100 dan mulai terjadi perubahan di titik antara 100 dan 200 seperti yang ditunjukkan pada garis putus-putus. Selanjutnya untuk menentukan titik nilai ambang digunakan grafik range pendugaan parameter dengan nilai ambang pada Gambar
14
Rata-rata nilai kejadian ekstrim
4. Berdasarkan grafik pada Gambar 4, nilai dugaan kedua parameter konstan di nilai 145 maka nilai ambang yang dipilih pada kasus ini adalah 145. Plot diagnostik kecocokan data dengan sebaran dapat dilihat pada Lampiran 1. Grafik fungsi kepekatan peluang dengan data menunjukkan bahwa data ekstrim yang diambil mengikuti sebaran GPD. Sementara itu hasil plot peluang dan plot kuantil-kuantilnya pada Lampiran 1 juga menunjukkan data mengikuti sebaran teoritisnya yaitu GPD. Hal ini dikarenakan kedua plot membentuk garis lurus mengikuti sebaran teoritisnya.
Parameter skala
Gambar 3 Grafik MRLP data curah hujan
Parameter bentuk
Nilai ambang
Nilai ambang
Gambar 4 Grafik hubungan pendugaan parameter GPD dengan nilai ambang pada data curah hujan rata-rata
15 Pengujian formal perlu dilakukan untuk meningkatkan kepercayaan data ekstrim yang diambil sudah mengikuti sebaran GPD. Pada penelitian kali ini pengujian formal dilakukan dengan uji Kolmogorov-Smirnov. Hasil pengujian menghasilkan statistik hitung 0.047, nilai ini lebih kecil dari nilai kritisnya 0.11 pada α = 5% sehingga hipotesis nol tidak ditolak. Hal ini menunjukkan bahwa data curah hujan ekstrim yang terpilih mengikuti sebaran GPD. Pemodelan Data GCM Data GCM yang digunakan adalah data presipitasi dengan ukuran grid 6 × 6 dan 8 × 8. Peubah GCM yang berjumlah banyak ini dianggap terlalu besar dan saling berkorelasi, sehingga akan mengakibatkan kondisi buruk. Analisis komponen utama dilakukan untuk mengatasai masalah ini. Hasil analisis komponen utama dapat dilihat pada Tabel 1. Jumlah komponen yang memiliki akar ciri lebih besar dari satu (𝜆𝑖 > 1) pada data grid 6 × 6 adalah sebanyak dua komponen. Proporsi kumulatif yang dihasilkan KU ke-1 dan KU ke-2 adalah 96% artinya dua KU mampu menjelaskan 96% peubah asalnya. Sedangkan pada data grid 8 × 8, akar ciri yang memiliki nilai yang lebih besar dari satu adalah empat komponen dengan proporsi kumulatif 95%. Dengan demikian jumlah komponen utama yang digunakan pada data grid 6 × 6 adalah sebanyak dua komponen dan pada data grid 8 × 8 sebanyak empat komponen dimana keragaman yang mampu dijelaskan keduanya hampir sama. Tabel 1 Nilai akar ciri dan proporsi keragaman kumulatif dari hasil analisis komponen utama Proporsi Ukuran Grid KU keAkar ciri Proporsi keragaman kumulatif 1 5.65 0.89 0.89 2 1.59 0.07 0.96 6 ×6 3 0.86 0.02 0.98 4 0.53 0.01 0.98 5 0.40 0.00 0.99 1 6.51 0.66 0.66 2 3.53 0.19 0.86 8×8 3 2.14 0.07 0.93 4 1.19 0.02 0.95 5 0.85 0.01 0.96 Analisis selanjutnya adalah data curah hujan rata-rata (y) akan dimodelkan dengan skor komponen utama terpilih, dua komponen utama untuk data grid 6 × 6 dan empat komponen utama untuk data grid 8 × 8. Pemodelan analisis VGAM dilakukan dengan memodelkan curah hujan yang berada diatas nilai ambang sebagai kombinasi aditif dari skor komponen utama terpilih. Eksplorasi skor komponen utama terhadap curah hujan rata-rata pada Lampiran 2 mengindikasikan bahwa hubungan antara peubah respon dengan masing-masing
16 peubah penjelas ada yang tidak linear, hal ini terlihat dari sebaran data yang membentuk pola tidak linear. Penelitian kali ini terdiri dari beberapa alternatif model untuk memperoleh model yang paling baik. Model pertama adalah model dasar tanpa dummy, model kedua adalah model dengan menambahkan satu peubah dummy (1 untuk bulan hujan dan 0 untuk bulan kering), dan model ketiga adalah model dengan empat kategori (hujan, peralihan hujan-kering, kering dan peralihan kering-hujan). Model aditif terampat menggunakan pemulusan dan teknik pemulusan yang digunakan dalam model ini adalah pemulusan vektor spline. Banyaknya pemulusan yang dilakukan ditunjukkan dengan besarnya derajat bebas. Semakin besar derajat bebas maka pemulusan yang dilakukan semakin banyak. Penentuan derajat bebas dilakukan dengan membuat model semua kemungkinan derajat bebas pada setiap komponen dari 2 sampai banyaknya observasi. Derajat bebas terbaik adalah model yang memiliki nilai AIC paling kecil dibandingkan dengan derajat bebas yang lain. Nilai AIC setiap model dapat dilihat pada Tabel 2. Berdasarkan hasil tersebut pada data grid 6 × 6, derajat bebas yang dipilih untuk KU1 adalah 9 dan 7 untuk KU2, sedangkan untuk data grid 8 × 8 derajat bebas yang dipilih untuk KU1, KU2, KU3, dan KU4 berturutturut adalah 9, 3, 2 dan 4. Tabel 2 Nilai AIC dari model VGAM pada beberapa derajat bebas Grid 6 × 6 Grid 8 × 8 derajat bebas KU1 KU2 KU1 KU2 KU3 KU4 2 1410.91 1405.08 1416.89 1405.95 1398.11 1414.23 3 1410.40 1403.45 1417.17 1405.11 1398.73 1411.60 4 1408.92 1402.81 1416.67 1405.68 1398.81 1410.36 5 1407.50 1402.20 1416.14 1406.03 1398.80 1410.38 6 1406.65 1401.76 1415.78 1406.24 1398.86 1410.72 7 1406.26 1401.66 1415.44 1406.32 1399.10 1410.96 8 1406.16 1401.89 1415.18 1406.33 1399.56 1411.05 9 1406.15 1402.32 1415.12 1406.33 1400.19 1411.08 10 1406.16 1402.83 1415.33 1406.42 1400.92 1411.17 Pengepasan model terhadap data 2008 dilakukan dengan menggunakan kuantil ke 50. Hasil pengepasan untuk masing-masing grid 6 × 6 dan 8 × 8 ditampilkan pada Gambar 5 dan Gambar 6. Berdasarkan kedua gambar tersebut, terlihat bahwa pola dugaan curah hujan data grid 6 × 6 dan 8 × 8 hampir sama. Pada kedua data terlihat bahwa hasil nilai dugaan curah hujan belum membentuk pola monsun seperti pola curah hujan aktual. Berdasarkan Gambar 5 dan Gambar 6, kedua model menduga curah hujan pada musim kering jauh diatas nilai aktualnya dan menduga curah hujan ekstrim pada musim hujan dibawah nilai aktual. Terdapat beberapa bulan yang memiliki nilai dugaan mendekati nilai aktual yaitu bulan November dan Desember. Nilai dugaan curah hujan tahun 2008 masing-masing model secara lengkap dapat dilihat pada Lampiran 4.
17 Penambahan peubah dummy tidak memberikan hasil yang signifikan pada pendugaan curah hujan. Dugaan yang dihasilkan hampir sama antara model awal tanpa dummy dengan model yang sudah ditambahkan dummy. Walau demikian, terdapat sedikit penurunan dalam nilai RMSEP. 500 450 400 350 300 250 200 150 100 50 0
Model tanpa dummy
Model dengan 2 dummy
Model 4 dummy
Aktual
Gambar 5 Pengepasan model VGAM pada data berukuran 6x6
500 450 400 350 300 250 200 150 100 50 0
Model tanpa dummy
Model dengan 2 dummy
Model 4 dummy
Aktual
Gambar 6 Pengepasan model VGAM pada data berukuran 8x8 Nilai RMSEP pada Tabel 3 menunjukkan bahwa terjadi sedikit penurunan nilai RMSEP dari model sebelum diberi peubah dummy dengan model setelah pemberian peubah dummy. Selanjutnya untuk mendapatkan hasil yang lebih baik, pemodelan dilakukan kembali dengan membagi bulan kedalam empat kelompok (Sutikno 2013). Data yang digunakan adalah data grid 6 × 6 karena berdasarkan nilai RMSEP, data grid 6 × 6 menghasilkan nilai RMSEP yang lebih kecil jika dibandingkan dengan data grid 8 × 8.
18 Tabel 3 Nilai RMSEP untuk model awal Grid 6 × 6 Grid 8 × 8 (pada kuantil) (pada kuantil) Model 50 75 90 50 75 90 Model tanpa dummy 126.66 134.32 155.63 129.67 138.32 158.56 Model dengan 2 kategori 122.23 126.48 144.95 125.80 132.76 152.39 Model dengan 4 kategori 121.07 122.97 138.76 126.19 132.61 151.36 Pada tahap selanjutnya, bulan dalam satu tahun dibagi menjadi empat kelompok yaitu bulan hujan, bulan peralihan hujan-kering, bulan kering dan bulan peralihan kering-hujan. Identifikasi nilai ekstrim dilakukan untuk setiap kelompok bulan. Gambar 7 adalah hasil MRLP untuk setiap kelompok bulan. Pemilihan nilai ambang untuk kelompok bulan hujan dapat dilihat pada Gambar 7(a), terlihat perubahan mulai terjadi di titik antara 100 dan 200 dan berdasarkan grafik hubungan pendugaan parameter GPD dengan nilai ambang pada Lampiran 4(a) nilai dugaan parameter konstan di nilai 145 maka nilai ambang yang dipilih pada kelompok bulan hujan adalah 145. Hal yang sam dilakukan untuk kelompok bulan lainnya dan menghasilkan nilai ambang untuk bulan peralihan hujan-kering, bulan kering dan bulan peralihan kering-hujan masing-masing adalah 100,10 dan 45.
(a)
(b)
(c)
(d)
Gambar 7 MRLP plot curah hujan masing-masing kelompok bulan. (a) bulan hujan, (b) bulan peralihan hujan-keirng , (c) bulan kering, (d) bulan peralihan kering-hujan Kecocokan sebaran GPD data dilihat berdasarkan dua hal yaitu plot diagnostik sisaan dan uji formal. Plot diagnostik kecocokan dapat dilihat pada Lampiran 5. Berdasarkan Lampiran 5 dapat dilihat bahwa sebaran data mengikuti
19 sebaran GPD. Plot kuantil-kuantil yang dihasilkan menunjukkan bahwa sebaran data mengikuti sebaran empiriknya. Tabel 4 Hasil Uji sebaran GPD dengan uji Kolmogorov-Smirnov Titik daerah Statistik Kesimpulan Kelompok kritis hitung Sebaran (α=5%) Bulan Hujan 0.05 0.15 menyebar GPD Bulan peralihan hujan-kering 0.07 0.17 menyebar GPD Bulan kering 0.07 0.17 menyebar GPD Bulan peralihan kering-hujan 0.07 0.19 menyebar GPD Hasil pengujian uji Kolmogorov-Smirnov dapat dilihat pada Tabel 4. Berdasarkan Tabel 4 setiap kelompok bulan menghasilkan statistik hitung yang lebih kecil dibandingkan dengan nilai kritisnya pada α = 5% sehingga tidak tolak hipotesis nol. Hal ini menunjukkan bahwa data curah hujan ekstrim setiap kelompok bulan mengikuti sebaran GPD. Tabel 5 Nilai AIC dari model VGAM pada beberapa derajat bebas pemulus setiap komponen utama GCM Kelompok Bulan
Hujan
Peralihan hujankering
db
KU1 2 3 4 5 6 7 2 3 4 5 2 3
Bulan kering
Peralihan keringhujan
4 5 2 3 4 5
893.47 893.78 894.69 895.50 896.35 897.40 603.91 605.38 606.89 608.57 573.11 573.22 573.24 572.89 573.41 497.13 501.22 Tidak Konvergen
KU2 891.34 889.32 889.82 890.37 890.83 891.32 603.69 605.02 606.64 608.47 587.32 588.34 589.12 589.63 590.11 505.75 506.40 506.76 507.65
Tabel 5 menunjukkan nilai AIC yang mungkin untuk setiap komponen utama pada masing-masing kelompok bulan. Tabel 5 memperlihatkan bahwa pada
20
S (KU2, df = 5)
S (KU1, df = 2)
kelompok bulan hujan derajat bebas yang terpilih untuk KU1 sampai KU4 masing-masing adalah 2, 5. Hal ini berarti pada bulan hujan, KU2 memiliki pengaruh nonlinier. Sedangkan untuk kelompok bulan peralihan hujan-kering, bulan kering dan peralihan kering-hujan derajat bebas yang terpilih untuk semua komponen adalah 2 artinya bentuk hubungan cruah hujan dengan komponen utama yang paling optimal adalah bentuk linier.
KU1
KU2
S (KU2, df = 2)
S (KU1, df = 2)
(a) Bulan hujan
KU1
KU2
S (KU2 df = 2)
S (KU1, df = 2)
(b) Bulan peralihan hujan-kering
KU2
KU1
S (KU2 df = 2)
S (KU1, df = 2)
(c) Bulan peralihan kering
KU1
KU2
(d) Bulan peralihan kering-hujan Gambar 8 Pengepasan model VGAM pada model kelompok curah hujan
21 Pada analisis VGAM, parameter bentuk (ξ) hanya dimodelkan menjadi intersep. Gambar 8 memperlihatkan hasil pengepasan dengan ±2 simpangan baku. Berdasarkan Gambar 8, skor KU memiliki hubungan yang linear dan nonlinear terhadap curah hujan bulanan. Misalnya pada bulan hujan, skor komponen 1 memiliki hubungan yang linear negatif sedangkan skor komponen 2 memiliki hubungan nonlinear. Hasil pemodelan dengan VGAM dapat dilihat pada Tabel 6. Berdasarkan Tabel 6 hasil dugaan untuk parameter bentuk (ξ) pada kelompok bulan hujan adalah -0.41, dan pada kelompok peralihan hujan-kering, kering serta peralihan kering-hujan memiliki dugaan parameter ξ masing-masing -0.47, -0.39, dan -0.49. Hasil ini sejalan dengan hasil yang ditunjukkan pada MRLP bahwa bentuk yang linear menurun mengindikasikan bahwa nilai dari parameter bentuk adalah negatif. Tabel 6. Penduga parameter VGAM masing-masing kelompok bulan Derajat Kelompok Komponen Bebas Koefisien Pemulus Intersep 1 0 3.84 Intersep 2 0 -2.32 Bulan Hujan f(KU1) 2 -0.13 f(KU2) 5 0.23 Intersep 1 0 3.82 Intersep 2 0 -3.85 Bulan peralihan hujan-kering f(KU1) 2 -0.05 f(KU2) 2 -0.22 Intersep 1 0 4.49 Intersep 2 0 -2.21 Bulan kering f(KU1) 2 -0.18 f(KU2) 2 -0.29 Intersep 1 0 7.44 Intersep 2 0 -10.91 Bulan peralihan kering - hujan f(KU1) 2 0.44 f(KU2) 2 -0.39 Setelah didapatkan model maka langkah selanjutnya adalah melakukan validasi dari model yang dibuat terhadap data curah hujan tahun 2008 . Visualisai prediksi kejadian ekstrim untuk data tahun 2008 dapat dilihat pada Gambar 9.. Terlihat pada Gambar 9 bahwa curah hujan pada bulan Februari tahun 2008 merupakan curah hujan maksimum untuk tahun tersebut dan nilai dugaan yang dihasilkan untuk bulan Februari juga merupakan curah hujan maksimum yang kemudian diikuti oleh bulan Januari. Bulan Agustus adalah bulan yang memiliki curah hujan minimum. Pola dugaan sudah mengikuti pola data aktual, terlihat pola data dugaan berbentuk pola monsun .
22
Curah Hujan (mm/bulan)
700 600 500
400 q95
300
q90
200
q75
100
aktual
0
Bulan Gambar 9. Grafik hasil validasi silang antara nilai dugaan dengan nilai aktual terhadap data curah hujan tahun 2008 pada quantil atas Bulan Februari yang memang sering terjadi curah hujan ekstrim memiliki dugaan sebesar 411 pada kuantil 75. Dugaan ini paling mendekati nilai aktualnya, sedangkan nilai dugaan yang paling mendekati bulan Januari adalah 379 terjadi pada kuantil 50. Pada Gambar 9 dapat dilihat bahwa posisi dugaan curah hujan untuk bulan kering ( Juli, Agustus dan September ) tidak berada pada kuantil atas. Menggunakan dugaan pada kuantil atas akan mengakibatkan jarak nilai dugaan curah hujan yang jauh dari nilai aktualnya. Oleh karena itu dugaan curah hujan untuk bulan-bulan kering dapat dilihat pada hasil dugaan kuantil bawah pada Gambar 10. Curah Hujan (mm/bulan)
700 600 q95
500
q90
400
q75
300
q50
200
q25
100
q10
0
q5 Aktual
Bulan Gambar 10. Grafik hasil validasi silang antara nilai dugaan dengan nilai aktual terhadap data curah hujan tahun 2008 pada kuantil bawah sampai kuantil atas.
23 Gambar 10 menyajikan dugaan curah hujan dari kuantil bawah sampai dengan kuantil atas. Berdasarkan gambar tersebut, dpaat dilihat bahwa untuk bulan-bulan kering, dugaan curah hujan yang tepat terjadi pada kuantil ke-5. Dugaan curah hujan pada kuantil ke-5 untuk bulan Juli, Agustus dan September masing-masing adalah 11.8, 11.6 dan 48.6. Nilai ini sudah mendekati nilai aktualnya untuk bulan Juli dan Agustus tetapi tidak untuk bulan September. Hal ini kemungkinan dapat disebabkan karena bulan September berbatasan dengan kelompok bulan peralihan. Dibandingkan dengan Gambar 5 dan Gambar 6, pola dugaan curah hujan bulan kering pada Gambar 10 jauh lebih mendekati aktualnya. Pembagian dugaan pada kuantil yang berbeda menghasilkan RMSEP yang lebih kecil dibandingkan degan menggunakan kuantil yang sama untuk setiap kelompok bulan. Hal ini dapat dilihat pada Tabel 7. Berdasarkan Tabel 7 terlihat bahwa nilai RMSEP untuk bulan kering pada kuantil 5 dan 10 kecil yaitu 9.2 dan 9.3. Hal serupa juga terjadi untuk bulan hujan dan bulan peralihan kering hujan, yaitu memiliki nilai RMSEP yang kecil pada kuantil 75 untuk bulan hujan , kuantil 25 untuk bulan peralihan kering-hujan. Hasil nilai RMSEP ini dapat dibandingkan dengan hasil RMSEP keseluruhan pada Tabel 8, yaitu rata-rata nilai RMSEP yang dihasilkan lebih besar. Tabel 7. Nilai RMSEP data dugaan pada berbagai kuantil Kategori Kuantil RMSEP 75 49.84 Bulan Hujan 90 117.62 95 159.85 25 98.55 Bulan peralihan hujan-kering 50 93.76 75 97.26 5 9.28 Bulan kering 10 9.32 25 13.02 10 51.85 Bulan peralihan kering-hujan 25 44.41 50 55.89
Tabel 8. Nilai RMSEP pada panjang data dugaan yang berbeda RMSEP pada kuantil Data Historis Data Dugaan 50 75 90 1979-2007 2008 74.17 79.17 113.19 1979-2006 2007-2008 72.04 80.44 113.29 1979-2005 2006-2008 68.14 82.60 117.92 1979-2004 2005-2008 63.14 85.45 123.81 1979-2003 2004-2008 73.17 92.25 127.14 Simpangan baku 4.53 5.19 6.26
95 136.44 136.27 141.03 147.11 148.92 5.89
Validasi dilakukan untuk mengetahui seberapa valid model yang dibangun. Hasil ini dapat dilihat pada Tabel 8. Nilai RMSEP yang dihasilkan pada
24 kuantil 50 lebih kecil dibandingkan dengan RMSEP pada kuantil yang lain. Berdasarkan Tabel 8, RMSEP yang dihasilkan tidak terlalu berbeda jauh meskipun panjang data dugaan bertambah. Hal ini terlihat dari nilai simpangan baku pada setiap kuantil. Secara keseluruhan kuantil 75 menghasilkan simpangan baku RMSEP 5.1, kuantil 90 menghasilkan 6.2 dan kuantil 95 menghasilkan 5.8. Sehingga dapat disimpulkan bahwa model yang dibangun dengan VGAM GPD sudah valid. Nilai korelasi dihitung untuk mengetahui kesamaan pola dugaan dengan aktual. Semakin tinggi nilai korelasi berarti dugaan memiliki pola yang mirip dengan curah hujan aktual. Hasil korelasi pada beberapa tahun dugaan dapat dilihat pada Tabel 9. Berdasarkan Tabel 9, nilai korelasi yang dihasilkan semuanya diatas 0.8 dengan simpangan baku 0.03. Hal ini berarti ada korelasi yang cukup kuat antara dugaan dengan data aktual. Korelasi yang positif menunjukkan pola yang sama antara dugaan dengan data aktual. Simpangan baku yang kecil untuk RMSEP dan korelasi menunjukkan bahwa meskipun panjang data dugaan ditambah sampai ke lima tahun ke depan, nilai dugaan tidak memili eror yang semakin besar dan tetap memiliki pola yang sama. Tabel 9. Nilai korelasi data dugaan dengan data aktual pada panjang data dugaan yang berbeda Korelasi pada kuantil Data Historis Data Dugaan 50 75 90 95 1979-2007 2008 0.91 0.92 0.93 0.93 1979-2006 2007-2008 0.86 0.85 0.85 0.84 1979-2005 2006-2008 0.87 0.88 0.87 0.87 1979-2004 2005-2008 0.86 0.86 0.85 0.85 1979-2003 2004-2008 0.82 0.82 0.82 0.81 Simpangan baku 0.03 0.04 0.04 0.04
25
5 SIMPULAN DAN SARAN Simpulan Berdasarkan penelitian ini dapat disimpulkan bahwa penggunaan ukuran grid 6 × 6 menghasilkan nilai dugaan yang lebih baik dibandingkan dengan menggunakan ukuran 8 × 8. Pemodelan curah hujan dengan menggunakan VGAM yang mengikuti sebaran teori nilai ekstrim GPD di Kabupaten Indramayu memberikan pola dugaan yang mendekati data aktual. Penggunaaan kuantil yang berbeda untuk setiap kelompok bulan menghasilkan dugaan yang lebih baik. Kelompok bulan kering memiliki nilai dugaan yang mendekati aktual pada kuantil rendah yaitu kuantil 5 dan 10, kelompok bulan hujan pada kuantil 75, kelompok bulan peralihan hujan-kering pada kuantil 50, sedangkan pada kelompok bulan peralihan kering-hujan dugaan terbaik yaitu pada kuntil 25. Kebaikan model menunjukkan bahwa model VGAM dengan GPD menghasilkan dugaan yang konsisten meskipun panjang data dugaan ditambah, nilai RMSEP yang dihasilkan tidak jauh berbeda jika dibandingkan dengan melakukan dugaan pada satu tahun. Saran Pemilihan nilai ambang yang menggunakan MRLP sedikit subjektif, oleh karena itu perlu dilakukan kajian mengenai penentuan nilai MRLP secara teoritik dan tidak subjektif. Selain itu, erlu dilakukan kajian lebih lanjut juga mengenai pemilihan derajat pemulusan yang optimum untuk metode VGAM yang mengikuti sebaran GPD sehingga menghasilkan model yang konvergen.
26
DAFTAR PUSTAKA Coles S. 2001. An Introduction to Statistical Modeling of Extreme Values. London: Springer Green P J. 1984. Iteratively Reweighted Least Square for Maximum Likelihood Estimation, and some Robust and Resistant Alternatives. Journal of the Royal Statistical Society Series B (Methodological). Vol 46, No 2 : 149192 Handayani L. 2014. Statistical Downscaling dengan Model Aditif Terampat untuk Pendugaan Curah Hujan EKstrim [Tesis]. Bogor: Institut Pertanian Bogor. Hastie T, Tibshirani R.1987. Generalized Additive Models: Some Applications. Journal of the American Statistical Association. Vol 82 No 398. Li Y, Cai W, Campbell. 2005. Statistical Modeling of Extreme Rainfall in Southwest Western Australia. Journal of Climate. Vol 18: 852-863. Mackenzie M, Yee TW. 2002. Vector Generalized Aditif Models in Plant Ecology. Science Direct- Ecological Modeling 157 : 141-156. Mallor F, Nualart E, Omey E. 2009. An Introduction to Statistical Modeling ef Extreme Values (Application to Calculate Extreme Wind Speeds). Hub Research Paper Economic & Management November/36. McCullagh P, Nelder JP.1983. Generalized Linear Model. London-New York: Chapman and Hall. Mondiana YQ. 2012. Pemodelan Statistical Downscaling dengan Regresi Kuantil untuk Pendugaan Curah Hujan Ekstrim. [tesis] Bogor : Sekolah Pascasarjana, Institut Pertanian Bogor. Rizki A. 2014. Pemodelan Semiparametrik Statistical Downscaling untuk Prediksi Curah Hujan di Kabupaten Indramayu [Tesis].Bogor: Institut Pertanian Bogor. Sutikno. 2013. Estimasi Parameter Generalized Pareto Distribution Pada Kasus Identifikasi Perubahan IKlim di Sentra Produksi Padi Jawa Timur. Jurnal Sains dan Seni POMITS Vol 2 No 2:141-146. Wigena AH. 2006. Pemodelan Statistical Downscaling dengan Regresi Projection Pursuit untuk Peramalan Curah Hujan Bulanan (Kasus Curah Hujan Bulanan di Indramayu). [disertasi] Bogor : Sekolah Pascasarjana, Institut Pertanian Bogor. Wilby RL, Charles SP, Zorita E, Timbal B, Whetton P, Mearns LO. 2004. Guideline for Use of Climate Scenarios Developed from Statistical Downscalling Methods. Task Group on Data and Scenario Support for Impact and Climate Analysis (TGICA). Wilby RL, Dawson CW. 2007. A Decision Support Tool for the Assesment of Regional Climate Change Impacts. UK. Yee TW, Wild CJ. 1996. Vector Generalized Aditif Models. Journal of Royal Statistics.58:481-493. Yee TW, Stephenson AG. 2007. “Vector Generalized Linear and Additive Extreme Value Models” in Article in Extremes Springer. Yee TW. 2007. VGAM Family Function for Extreme Value Data. University of Auckland. New Zealand. Yee TW. 2015. Vector Generalized Linear and Additive Models (with Implemantation in R). London: Springer
27 Lampiran 1 Plot diagnostik kecocokan data curah hujan dengan GPD
28 Lampiran 2 Plot curah hujan ekstrim dengan skor komponen utama terpilih 2a. Menggunakan data ukuran grid 6 × 6
2b. Menggunakan data ukuran grid 8 × 8
29 Lampiran 3 Nilai dugaan curah hujan tahun 2008 pada model awal 3a. Model tanpa dummy data ukuran grid 6 × 6 Bulan
Aktual
Jan Feb Mar Apr Mei Jun Jul Agustus Sept Okt Nov Des
351.47 439.33 260.87 97.13 19.73 23.47 0.00 7.33 2.20 68.00 136.47 198.40
Kuantil 50 271.60 273.13 193.34 178.46 179.97 170.23 177.02 169.86 165.07 167.46 210.86 209.50
75 365.90 368.56 229.34 203.38 206.02 189.02 200.87 188.38 180.02 184.18 259.91 257.54
90 454.69 458.42 263.24 226.85 230.55 206.71 223.33 205.82 194.10 199.93 306.10 302.78
3b. Model dengan 2 kategori data ukuran grid 6 × 6 Kuantil Bulan Aktual 50 75 90 Jan 351.47 272.01 365.56 452.54 Feb 439.33 272.88 367.08 454.65 Mar 260.87 204.82 248.88 289.84 Apr 97.13 178.40 203.00 225.87 Mei 19.73 185.14 214.70 242.19 Jun 23.47 164.20 178.33 191.48 Jul 0.00 164.51 178.87 192.23 Agustus 7.33 159.68 170.49 180.54 Sept 2.20 156.10 164.27 171.87 Okt 68.00 169.60 187.71 204.55 Nov 136.47 210.95 259.53 304.69 Des 198.40 212.08 261.49 307.43
95 502.27 506.58 281.41 239.43 243.69 216.19 235.36 215.16 201.65 208.37 330.85 327.02
95 498.54 500.97 311.51 237.96 256.72 198.43 199.29 185.85 175.89 213.46 328.58 331.72
30 3c. Model dengan 4 kategori data ukuran grid 6 × 6 Kuantil Bulan Aktual 50 75 90 95 Jan Feb Mar Apr Mei Jun Jul Agustus Sept Okt Nov Des
351.47 439.33 260.87 97.13 19.73 23.47 0.00 7.33 2.20 68.00 136.47 198.40
271.38 261.63 191.38 178.19 178.07 155.72 150.75 148.65 162.01 165.77 214.10 209.15
365.55 348.53 225.94 202.91 202.71 163.72 155.04 151.36 174.68 181.24 265.58 256.94
454.24 430.38 258.49 226.20 225.91 171.24 159.08 153.92 186.62 195.81 314.08 301.95
3d. Model tanpa dummy data ukuran grid 8 × 8 Kuantil Bulan Aktual 50 75 90 Jan 351.47 258.13 341.99 420.52 Feb 439.33 255.38 337.21 413.83 Mar 260.87 192.64 227.95 261.01 Apr 97.13 183.13 211.40 237.87 Mei 19.73 175.44 198.00 219.13 Jun 23.47 157.74 167.19 176.03 Jul 0.00 160.13 171.35 181.86 Agustus 7.33 172.36 192.64 211.64 Sept 2.20 168.48 185.89 202.19 Okt 68.00 201.59 243.54 282.83 Nov 136.47 210.43 258.94 304.36 Des 198.40 230.90 294.57 354.20
501.79 474.26 275.94 238.69 238.35 175.28 161.25 155.29 193.02 203.63 340.08 326.09
95 462.37 454.66 278.63 251.97 230.39 180.75 187.46 221.76 210.88 303.76 328.56 385.97
31 3e. Model dengan 2 kategori data ukuran grid 8 × 8 Kuantil Bulan Aktual 50 75 90 95 Jan 351.47 269.03 360.40 445.35 490.29 Feb 439.33 258.77 342.58 420.50 461.72 Mar 260.87 196.48 234.41 269.67 288.33 Apr 97.13 177.57 201.56 223.87 235.67 Mei 19.73 175.72 198.35 219.40 230.53 Jun 23.47 156.56 165.08 173.00 177.19 Jul 0.00 158.80 168.97 178.43 183.43 Agustus 7.33 164.49 178.85 192.19 199.26 Sept 2.20 160.44 171.81 182.39 187.98 Okt 68.00 194.66 231.24 265.25 283.24 Nov 136.47 212.33 261.93 308.05 332.44 Des 198.40 231.49 295.20 354.44 385.77
3f. Model dengan 4 kategori data ukuran grid 8 × 8 Bulan
Aktual
Jan Feb Mar Apr Mei Jun Jul Agustus Sept Okt Nov Des
351.467 439.333 260.867 97.133 19.733 23.467 0 7.333 2.2 68 136.467 198.4
50 266.65 253.57 190.78 188.01 174.79 157.10 151.53 154.72 166.89 193.56 211.83 230.99
Kuantil 75 90 356.34 439.83 333.62 408.13 224.53 255.95 219.72 249.24 196.76 217.21 166.02 174.32 156.34 160.82 161.89 168.56 183.04 198.06 229.37 262.70 261.11 306.98 294.40 353.42
95 484.05 447.59 272.59 264.87 228.04 178.72 163.19 172.10 206.02 280.35 331.27 384.67
32 Lampiran 4 Grafik hubungan pendugaan parameter dengan nilai ambang pada data curah hujan setiap kelompok bulan 4a. Bulan hujan
4b. Bulan peralihan hujan-kering
4c. Bulan kering
4d. Bulan peralihan kering-hujan
33 Lampiran 5 Plot diagnostik untuk nilai ambang terpilih setiap kelompok bulan 5a. Bulan Hujan
5b. Bulan peralihan hujan-kering
34 5c. Bulan Kering
5d. Bulan peralihan kering-hujan
35 Lampiran 6 Plot pemulusan pada data curah hujan rata-rata berdasarkan kelompok bulan 6a. Bulan Hujan
6b. Bulan peralihan hujan-kering
36 6c. Bulan Kering
6d. Bulan Peralihan Hujan-kering
37 Lampiran 7 Nilai aktual dan dugaan curah hujan model setiap kelompok bulan pada tahun 2008 Bulan
Aktual
Jan Feb Mar Apr Mei Jun Jul Agustus Sept Okt Nov Des
351.47 439.33 260.87 97.13 19.73 23.47 0.00 7.33 2.20 68.00 136.47 198.40
5 156.18 157.75 104.52 104.28 103.48 13.30 11.80 11.30 48.66 50.71 53.83 150.56
10 167.73 170.91 109.17 108.67 107.06 16.74 13.68 12.66 52.42 56.58 62.90 156.30
25 204.85 213.23 123.99 122.68 118.47 27.99 19.83 17.10 64.37 75.24 91.74 174.77
Kuantil 50 285.77 308.57 151.82 136.28 133.53 79.53 34.18 26.64 88.90 109.91 146.32 221.03
75 387.32 426.58 188.69 162.10 157.38 130.46 51.89 38.82 119.93 155.81 217.96 275.88
90 479.52 533.72 221.63 185.17 178.69 177.48 68.24 50.08 147.48 196.54 281.53 325.68
95 527.11 589.01 238.36 196.88 189.51 202.19 76.83 55.99 161.36 217.07 313.57 351.38
38
RIWAYAT HIDUP Penulis dilahirkan di Bogor pada tanggal 21 Januari 1989, sebagai anak pertama dari pasangan Isman dan Nuriyah. Pendidikan sekolah menengah ditempuh di SMA Negeri 1 Bogor Program IPA, lulus pada tahun 2007. Pada tahun yang sama penulis diterima di program studi Statistika Institut Pertanian Bogor dan menyelesaikannya pada tahun 2011. Setelah lulus, penulis sempat bekerja di Marketing Riset Indonesia (MRI) sebagai research executive dan data analyst di telkomsel masing-masing selama satu tahu. Kesempatan untuk melanjutkan program master (S2) pada program studi Statistika, Sekolah Pascasarjana IPB, diperoleh pada tahun 2013 dengan program Beasiswa Pendidikan Pascasarjana Dalam Negeri (BPPDN) dari Direktorat Jendral Pendidikan Tinggi (Dikti). Selama perkuliahan penulis juga aktif sebagai asisten dosen mata kuliah Analisis Statistika (STK 511), rancangan percobaan dan model linier pada Sekolah Pascasarjana IPB. Karya ilmiah ini telah dipresentasikan dalam SEAMS UGM International Conference on Mathematics and Its Applications di Universitas Gadjah Mada, Yogyakarta pada tanggal 18-21 Agustus 2015.
