MODEL MATEMATIK DEMAM BERDARAH DENGUE DENGAN NYAMUK Aedes albopictus SEBAGAI VEKTOR
JAMES U. L. MANGOBI
SEKOLAH PASCASARJANA INSTITUT PERTANIAN BOGOR BOGOR 2011
i
PERNYATAAN MENGENAI TESIS DAN SUMBER INFORMASI
Dengan ini saya menyatakan bahwa tesis dengan judul Model Matematik Demam Berdarah Dengue dengan Nyamuk Aedes albopictus sebagai Vektor adalah karya saya dengan arahan 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 penulis lain telah disebutkan dalam teks dan dicantumkan dalam Daftar Pustaka di bagian akhir tesis ini.
Bogor, Agustus 2011 James U.L. Mangobi NIM G551090351
ii
iii
ABSTRACT
JAMES URIEL LIVINGSTONE MANGOBI. Mathematical Model of Dengue Hemorrhagic Fever with Aedes albopictus Mosquitos as Vector. Supervised by PAIAN SIANTURI and N.K. KUTHA ARDANA. Dengue Hemorrhagic Fever (DHF) is an acute febrile illness caused by a dengue virus. This virus has four serotypes, i.e. Dengue I - IV. The dengue virus is transmitted by various species of Aedes mosquitoes. Mathematical model can be used to study the spread of the disease. The mathematical model discussed in this paper is SEIR model. The main vector of the disease is mosquito of the Aedes albopictus type. In the SEIR model, an analysis is performed to assess the stability of the equilibrium point and numerical simulations. There are two equilibrium points obtained. The first equilibrium point is a disease-free equilibrium (DFE), which is stable, given the basic reproductive number ℜ < 1. The second equilibrium point is called an endemic point, which stability is guaranteed if the value ℜ > 1. The numerical simulations show that increasing mosquitoes mortality rate makes the number of exposed susceptible humans decrease. Furthermore, increase in the average bite of infected mosquito will increase the number of exposed susceptible humans. For the mosquito population, increasing mosquitoes mortality rate will decrease the number of exposed susceptible mosquitoes. Finally, increase in the average bite of infected mosquito will increase the number of exposed susceptible mosquitoes. Keywords: dengue hemorrhagic fever, Aedes albopictus, SEIR model, equilibrium point, stability analysis
iv
v
RINGKASAN
JAMES U.L. MANGOBI. Model Matematik Demam Berdarah Dengue dengan Nyamuk Aedes albopictus sebagai Vektor. Dibimbing oleh PAIAN SIANTURI dan N.K. KUTHA ARDANA. Demam Berdarah Dengue (DBD) merupakan penyakit demam akut yang disebabkan oleh virus Dengue. Virus ini memiliki empat serotype virus, yaitu Dengue I – IV (Gubler 1998). Virus ini ditularkan oleh berbagai nyamuk spesies Aedes. Nyamuk ini merupakan vektor yang sangat efisien, sehingga penyakit ini menjadi wabah (epidemi). Berbagai program pengendalian epidemi DBD menjadi prioritas utama WHO dan departemen kesehatan di banyak negara selama ini. Di Indonesia, upaya ini terbilang belum berhasil karena adanya berbagai kendala baik secara teknis maupun non-teknis. Sehubungan dengan banyaknya kendala tersebut, perlu adanya suatu penelitian dan pemikiran yang dilakukan. Pemodelan Matematika dapat membantu memahami dan mengidentifikasi hubungan penyebaran penyakit DBD dengan berbagai parameter epidemiologi. Model matematik yang dimaksud diantaranya ialah model Susceptible, Infected, Recovered (SIR) dan model Susceptible, Exposed, Infected, Recovered (SEIR). Dalam tesis ini dibahas model SEIR yang dikenalkan oleh Erickson et al. (2010). Vektor utama dalam model ini adalah nyamuk Aedes albopictus, sehubungan dengan banyaknya kasus DBD yang disebabkan oleh nyamuk ini dan juga, mempunyai cakupan yang lebih besar dan lebih sulit dikendalikan. Pada model SEIR dilakukan analisis kestabilan dan simulasi numerik dengan pemrograman berbasis fungsional menggunakan software Mathematica 8.0 (Wolfram Research, Inc, Champaign, IL). Dalam proses analisis kestabilan, ditentukan titik-titik tetap, bilangan reproduksi dasar dan kestabilan dari titik-titik tetap tersebut. Simulasi dilakukan untuk melihat pengaruh perubahan laju kematian nyamuk dan rata-rata gigitan nyamuk terinfeksi terhadap populasi manusia dan nyamuk. Bilangan reproduksi dasar (ℜ ) merupakan suatu ukuran potensi penyebaran penyakit dalam suatu populasi. Bilangan reproduksi dasar didefinisikan sebagai nilai harapan banyaknya populasi rentan yang menjadi terinfeksi selama masa infeksi berlangsung. Jika ℜ < 1, maka satu nyamuk terinfeksi akan menginfeksi kurang dari satu manusia rentan atau satu manusia terinfeksi akan menginfeksi kurang dari satu nyamuk rentan, sehingga lambat laun penyakit DBD akan hilang dari populasi. Jika ℜ > 1, maka satu nyamuk terinfeksi akan menginfeksi lebih dari satu manusia rentan atau satu manusia terinfeksi akan menginfeksi lebih dari satu nyamuk rentan, sehingga penyakit DBD akan bertahan di dalam populasi. Nilai ℜ ini dapat ditentukan antara lain dari the next generation matrix, yaitu matriks yang dibentuk hanya pada sub-subpopulasi yang menyebabkan infeksi saja ( ܧ , ܧ௩ , ܫ , dan ܫ௩ ).
vi
Analisis kestabilan menghasilkan dua titik tetap. Pertama, titik tetap tanpa penyakit, yang selalu ada dan merupakan titik yang stabil jika ℜ < 1. Kedua, titik tetap endemik. Kestabilan titik ini dijamin apabila nilai ℜ > 1. Simulasi menunjukkan bahwa jumlah tiap subpopulasi pada populasi manusia dan nyamuk akan stabil ke titik tetap tanpa penyakit ܶଵ ketika ℜ < 1 dan stabil ke titik tetap endemik ܶଶ ketika ℜ > 1. Pada populasi manusia, semakin besar laju kematian nyamuk maka jumlah manusia rentan yang menjadi terpapar semakin sedikit. Sedangkan semakin besar rata-rata gigitan nyamuk terinfeksi maka jumlah manusia rentan yang menjadi terpapar semakin banyak. Pada populasi nyamuk, semakin besar laju kematian nyamuk maka jumlah nyamuk rentan yang menjadi terpapar semakin sedikit. Sedangkan semakin besar rata-rata gigitan nyamuk terinfeksi maka jumlah nyamuk rentan yang menjadi terpapar semakin banyak. Kata kunci: demam berdarah dengue, Aedes albopictus, model SEIR, titik tetap, analisis kestabilan
vii
© Hak Cipta milik IPB, tahun 2011
Hak Cipta dilindungi Undang-Undang 1. Dilarang mengutip sebagian atau seluruh karya tulis ini tanpa mencantumkan atau menyebutkan sumbernya. a. Pengutipan hanya untuk kepentingan pendidikan, penelitian penulisan karya ilmiah, penyusunan laporan, penulisan kritik, atau tinjauan suatu masalah. b. Pengutipan tersebut tidak merugikan kepentingan yang wajar IPB 2. Dilarang mengumumkan dan memperbanyak sebagian atau seluruh Karya tulis dalam bentuk apa pun tanpa izin IPB
viii
ix
MODEL MATEMATIK DEMAM BERDARAH DENGUE DENGAN NYAMUK Aedes albopictus SEBAGAI VEKTOR
JAMES U. L. MANGOBI
Tesis sebagai salah satu syarat untuk memperoleh gelar Magister Sains pada Program Studi Matematika Terapan
SEKOLAH PASCASARJANA INSTITUT PERTANIAN BOGOR BOGOR 2011
x
Penguji Luar Komisi pada Ujian Tesis: Dr. Toni Bakhtiar, M.Sc.
xi
Judul Tesis : Model Matematik Demam Berdarah Dengue dengan Nyamuk Aedes albopictus sebagai Vektor Nama : James U.L. Mangobi NIM : G551090351
Disetujui Komisi Pembimbing
Dr. Paian Sianturi Ketua
Ir. N.K. Kutha Ardana, M.Sc. Anggota
Diketahui
Ketua Program Studi Matematika Terapan
Dekan Sekolah Pascasarjana
Dr. Ir. Endar H. Nugrahani, M.S.
Dr. Ir. Dahrul Syah, M.Sc.Agr.
Tanggal Ujian: 22-08-2011
Tanggal Lulus: ............................
xii
xiii
“I can do all things through Christ Who strengthens me” (Philippians 4 : 13)
ku persembahkan tesis ini untuk orang tuaku, isteri tercinta dan keluargaku
xiv
xv
PRAKATA
Puji dan syukur penulis panjatkan kepada Tuhan Yang Maha Esa, karena atas kasih dan pertolongan-Nya sehingga penulis dapat menyelesaikan penelitian serta penulisan tesis ini dengan judul Model Matematik Demam Berdarah Dengue dengan Nyamuk Aedes albopictus sebagai Vektor. Penulisan tesis ini merupakan salah satu syarat memperoleh gelar Magister Sains pada program studi Matematika Terapan Sekolah Pascasarjana Institut Pertanian Bogor. Penulis menyadari bahwa bantuan-bantuan dan arahan-arahan dari kedua pembimbing sangat membantu dalam menyelesaikan karya tulis ini. Penulis sangat berterima kasih kepada Dr. Paian Sianturi selaku pembimbing I dan Ir. N.K. Kutha Ardana, M.Sc. selaku pembimbing II. Penulis menyampaikan terima kasih juga kepada: 1. Prof. Dr. Ph.E.A. Tuerah, M.Si. DEA selaku Rektor Universitas Negeri Manado. 2. Prof.Dr.Ir. Herry Suhardiyanto, M.Sc. selaku Rektor Institut Pertanian Bogor. 3. Dr. Ir. Dahrul Syah, M.Sc.Agr. selaku Dekan Sekolah Pascasarjana Institut Pertanian Bogor. 4. Dr. S.M. Salajang, M.Si. selaku Direktur Eksekutif Proyek I-MHERE Batch IV Universitas Negeri Manado yang telah memberikan beasiswa kepada penulis. 5. Dr. Ir. Endar H. Nugrahani, M.S. selaku Ketua Program Studi Matematika Terapan merangkap penguji dari Departemen Matematika. 6. Dr. Toni Bakhtiar, M.Sc. selaku penguji luar komisi pembimbing. 7. Seluruh dosen dan staf pegawai tata usaha Departemen Matematika. 8. Papa (alm.) dan Mama yang dengan tabah mendidik, membesarkan dan memberikan doa restu. 9. Isteriku tercinta dan seluruh keluargaku yang selalu memberikan dorongan dan mendoakan keberhasilan studiku. 10. Teman-teman penghuni Asrama Mahasiswa Sulawesi Utara “Sam Ratulangi” di Sempur, Bogor Baru I dan Bogor Baru II. 11. Seluruh mahasiswa Departemen Matematika khususnya teman-teman angkatan tahun 2009 di program studi Matematika Terapan. 12. Sahabat-sahabatku yang tak dapat saya sebutkan satu persatu yang telah banyak membantu penulis dalam penyelesaian tesis ini. Dengan harapan disertai dengan keyakinan kiranya Tuhan Yang Maha Esa akan membalas segala kebaikan dan bantuan dari bapak, ibu, saudara sekalian. Semoga penulisan tesis ini dapat memperkaya pengalaman belajar serta wawasan kita semua. Bogor, Agustus 2011 James U.L. Mangobi NIM G551090351
xvi
xvii
RIWAYAT HIDUP
Penulis dilahirkan di Desa Dagho Kecamatan Tamako Kabupaten Kepulauan Sangihe Provinsi Sulawesi Utara pada tanggal 15 Juli 1977 sebagai anak bungsu dari pasangan Bapak Franssiscus Mangobi (alm) dan Ibu Martha Pereman. Penulis mempunyai isteri bernama Meyni E. Lepa. Penulis menamatkan SD, SMP dan SMA di Tamako. Setelah lulus dari SMA Negeri Tamako, penulis melanjutkan studi S1 di Jurusan Matematika FMIPA Universitas Negeri Manado dan lulus pada tahun 2002. Setelah memperoleh gelar sarjana, penulis menjadi Dosen Luar Biasa di Jurusan Matematika dan Statistika FMIPA Universitas Sam Ratulangi Manado hingga akhir tahun 2004. Tahun 2005, penulis diangkat sebagai Pegawai Negeri Sipil dan menjadi staf pengajar di Jurusan Matematika FMIPA Universitas Negeri Manado hingga sekarang. Pada tahun 2009 penulis lulus seleksi masuk Program Magister pada Program Studi Matematika Terapan Institut Pertanian Bogor.
xviii
xix
DAFTAR ISI
Halaman DAFTAR TABEL .......................................................................................
xxi
DAFTAR GAMBAR .................................................................................. xxiii DAFTAR LAMPIRAN ............................................................................... xxv I.
PENDAHULUAN ............................................................................. 1.1 Latar Belakang ......................................................................... 1.2 Tujuan Penelitian ......................................................................
1 1 2
II.
LANDASAN TEORI ......................................................................... 2.1 Sistem Persamaan Diferensial (SPD) ........................................ 2.2 Titik Tetap ................................................................................ 2.3 Nilai Eigen dan Vektor Eigen ................................................... 2.4 Pelinearan ................................................................................. 2.5 Kestabilan Titik Tetap .............................................................. 2.6 Kriteria Routh-Hurwitz ............................................................ 2.7 Bilangan Reproduksi Dasar .......................................................
3 3 4 4 4 5 5 6
III.
MODEL MATEMATIK PENYEBARAN PENYAKIT DBD ............. 3.1 Model SIR ................................................................................. 3.2 Model SEIR ..............................................................................
8 8 11
IV.
HASIL DAN PEMBAHASAN ........................................................... 4.1 Penentuan Titik Tetap ................................................................ 4.2 Penentuan Bilangan Reproduksi Dasar ....................................... 4.3 Analisis Kestabilan Titik Tetap .................................................. 4.3.1 Kestabilan Titik Tetap Tanpa Penyakit ........................... 4.3.2 Kestabilan Titik Tetap Endemik ..................................... 4.4 Simulasi Dinamika Populasi Penularan Virus Dengue ............... 4.4.1 Nilai Parameter .............................................................. 4.4.2 Dinamika Populasi untuk Kondisi ℜ < 1 ...................... 4.4.3 Dinamika Populasi untuk Kondisi ℜ > 1 ......................
15 15 16 16 16 18 19 19 20 25
V.
SIMPULAN DAN SARAN ................................................................. . 5.1 Simpulan ................................................................................... 5.2 Saran .........................................................................................
31 31 32
DAFTAR PUSTAKA .................................................................................
33
LAMPIRAN ................................................................................................
35
xx
xxi
DAFTAR TABEL
Halaman 1.
Kondisi Kestabilan Titik Tetap ...........................................................
19
2.
Definisi dan Nilai Parameter Model SEIR dalam Simulasi Numerik ...
20
3.
Simulasi untuk Kondisi ℜ < 1 ..........................................................
21
4.
Simulasi untuk Kondisi ℜ > 1 ..........................................................
26
xxii
xxiii
DAFTAR GAMBAR
Halaman 1.
Skema penyebaran penyakit DBD model SIR .....................................
8
2.
Skema penyebaran penyakit DBD model SEIR ..................................
12
3.
Dinamika populasi manusia (a) dan populasi nyamuk (b) terhadap waktu ݐuntuk kondisi ℜ < 1 ...........................................................
20
Dinamika populasi manusia (a) rentan ܵ , (b) terpapar ܧ , (c) terinfeksi ܫ , dan (d) sembuh ܴ , serta populasi nyamuk (e) rentan ܵ ௩ , (f) terpapar ܧ௩ dan (g) terinfeksi ܫ௩ terhadap waktu ݐpada kondisi ℜ < 1 dan nilai parameter ߤ௩ diubah ...............................................
22
Dinamika populasi manusia (a) rentan ܵ , (b) terpapar ܧ , (c) terinfeksi ܫ , dan (d) sembuh ܴ , serta populasi nyamuk (e) rentan ܵ ௩ , (f) terpapar ܧ௩ dan (g) terinfeksi ܫ௩ terhadap waktu ݐpada kondisi ℜ < 1 dan nilai parameter ܾ diubah ................................................
24
Dinamika populasi manusia (a) dan populasi nyamuk (b) terhadap waktu ݐuntuk kondisi ℜ > 1 ...........................................................
25
Dinamika populasi manusia (a) rentan ܵ , (b) terpapar ܧ , (c) terinfeksi ܫ , dan (d) sembuh ܴ , serta populasi nyamuk (e) rentan ܵ ௩ , (f) terpapar ܧ௩ dan (g) terinfeksi ܫ௩ terhadap waktu ݐpada kondisi ℜ > 1 dan nilai parameter ߤ௩ diubah ...............................................
27
Dinamika populasi manusia (a) rentan ܵ , (b) terpapar ܧ , (c) terinfeksi ܫ , dan (d) sembuh ܴ , serta populasi nyamuk (e) rentan ܵ ௩ , (f) terpapar ܧ௩ dan (g) terinfeksi ܫ௩ terhadap waktu ݐpada kondisi ℜ > 1 dan nilai parameter ܾ diubah ................................................
29
4.
5.
6.
7.
8.
xxiv
xxv
DAFTAR LAMPIRAN
Halaman 1.
Penentuan Titik Tetap ........................................................................
35
2.
Penentuan Bilangan Reproduksi Dasar (ℜ ) .......................................
36
3.
Analisis Kestabilan Titik Tetap Tanpa Penyakit .................................
37
4.
Analisis Kestabilan Titik Tetap Endemik ...........................................
39
5.
Simulasi untuk Kondisi ℜ < 1 .........................................................
44
6.
Simulasi untuk Kondisi ℜ > 1 .........................................................
46
1
I. PENDAHULUAN
1.1
Latar Belakang Demam Berdarah Dengue (DBD) merupakan penyakit demam akut yang
disebabkan oleh virus Dengue. Virus ini memiliki empat serotype virus, yaitu Dengue I – IV (Gubler 1998). Virus Dengue ditularkan oleh berbagai nyamuk spesies Aedes. Nyamuk ini merupakan vektor yang sangat efisien karena adanya asosiasi nyamuk dengan kehidupan manusia. Juga, perilaku menggigit dan menghisap darah pada beberapa orang oleh satu nyamuk betina dewasa. Dengan demikian, begitu mudah penyakit ini menjadi wabah (epidemi) di dalam populasi manusia. Penyakit DBD ini banyak ditemukan di Indonesia. Tercatat telah empat kali Kejadian Luar Biasa (KLB) yakni tahun 1988, 1998, 2004 dan 2006. WHO memperkirakan sekitar 2,5 miliar penduduk dunia menghadapi risiko penyakit DBD (Anonym 2009). Dengan fakta tersebut, program pengendalian epidemi DBD menjadi prioritas utama WHO dan Departemen Kesehatan RI. Sejak tahun 1962, pencegahan epidemi DBD telah difokuskan pada pemberantasan nyamuk pembawa virus Dengue. Namun demikian kita pahami, upaya penanggulangan epidemi DBD di Indonesia masih jauh dari memuaskan. Berbagai kendala seperti sedikitnya anggaran pemerintah untuk penanggulangan epidemi, keterbatasan infrastruktur, dan kurangnya data dan informasi menjadi penyebab utama keterbelakangan kita dalam pencegahan dan penanggulangan epidemi ini. Pemodelan matematika dapat membantu memahami dan mengidentifikasi hubungan penyebaran penyakit DBD dengan berbagai parameter epidemiologi. Model matematik di antaranya ialah model Susceptible-Infected-Recovered (SIR) dan model Susceptible-Exposed-Infected-Recovered (SEIR). Dalam tesis ini dibahas model SEIR yang mengacu pada kajian Erickson et al. (2010). Vektor utama dalam model ini adalah nyamuk Aedes albopictus, sehubungan dengan banyaknya kasus DBD yang disebabkan oleh nyamuk
2
ini (Gratz 2004). Juga, nyamuk Aedes albopictus mempunyai cakupan yang lebih besar dan lebih sulit dikendalikan (Estrada-Franco and Craig 1995). Pada model SEIR dilakukan analisis kestabilan dan simulasi numerik dengan pemrograman berbasis fungsional menggunakan software Mathematica 8.0 (Wolfram Research, Inc, Champaign, IL). 1.2
Tujuan Penelitian Penelitian ini bertujuan untuk:
1.
Mengkaji model-model matematik penyebaran penyakit DBD.
2.
Melakukan analisis kestabilan titik tetap model SEIR.
3.
Melakukan simulasi numerik terhadap model SEIR untuk melihat pengaruh perubahan laju kematian nyamuk dan rata-rata gigitan nyamuk terinfeksi terhadap populasi manusia dan nyamuk.
3
II. LANDASAN TEORI
2.1
Sistem Persamaan Diferensial Biasa
Definisi 1 (Sistem Persamaan Diferensial Biasa Linear) Misalkan suatu sistem persamaan diferensial biasa dinyatakan sebagai = + ; = , ∈ ℝ
(1)
dengan adalah matriks koefisien konstan berukuran × dan adalah vektor Linear Orde Satu dengan kondisi awal = . Jika = , maka sistem konstan. Sistem persamaan (1) dinamakan Sistem Persamaan Diferensial Biasa
dikatakan homogen dan jika ≠ , maka sistem dikatakan takhomogen.
(Tu 1994)
Definisi 2 (Sistem Persamaan Diferensial Biasa Taklinear) Misalkan suatu sistem persamaan diferensial biasa dinyatakan sebagai
dengan
= ,
(2)
, , , … , , = dan , = , , … , ⋮ ⋮ , , , … ,
adalah fungsi taklinear dalam , , … , . Sistem persamaan (2) disebut Sistem Persamaan Diferensial Biasa Taklinear.
(Braun 1983)
Definisi 3 (Sistem Persamaan Diferensial Biasa Mandiri) Misalkan suatu Sistem Persamaan Diferensial Biasa dinyatakan sebagai = , ∈ ℝ
(3)
dengan merupakan fungsi kontinu bernilai real dari dan mempunyai turunan parsial kontinu. Sistem persamaan (3) disebut Sistem Persamaan Diferensial Biasa Mandiri (autonomous) karena tidak memuat t secara eksplisit di dalamnya. (Tu 1994)
4
2.2
Titik Tetap
disebut titik tetap, jika = . Titik tetap disebut juga titik sistem (3). Titik Misalkan diberikan sistem persamaan diferensial biasa sebagaimana pada
kritis atau titik kesetimbangan. Untuk selanjutnya digunakan istilah titik tetap.
(Tu 1994) 2.3
Nilai Eigen dan Vektor Eigen
Diberikan matriks koefisien konstan berukuran × dan sistem
persamaan diferensial biasa homogen = , = , ∈ ℝ . Suatu vektor
taknol di dalam ℝ disebut vektor eigen dari jika untuk suatu skalar berlaku:
= .
Nilai skalar dinamakan nilai eigen dari .
(4)
Untuk mencari nilai dari , maka sistem persamaan (4) dapat ditulis
dengan
− = .
(5)
! = | − | = 0.
(6)
adalah matriks identitas. Sistem persamaan (5) mempunyai solusi taknol
jika dan hanya jika
Persamaan (6) merupakan persamaan karakteristik matriks . 2.4
(Anton 1995)
Pelinearan Analisis kestabilan sistem persamaan diferensial taklinear dapat dilakukan
melalui pelinearan. Misalkan diberikan sistem persamaan diferensial biasa taklinear
= , ∈ ℝ
(7)
dengan ∈ ℝ adalah suatu fungsi bernilai vektor dalam (waktu) dan : % → ℝ adalah suatu fungsi mulus yang terdefinisi pada subhimpunan
% ⊂ ℝ .
5
, maka sistem Dengan menggunakan ekspansi Taylor di sekitar titik tetap
persamaan (7) dapat ditulis sebagai
≡ ) = *) + + ).
dengan * adalah matriks Jacobi
012
/032 . 014 , * = , = .032 . ⋮ .016 -032
012
034 014 034
⋮
⋯
⋯
⋱ 016 ⋯
034
(8) 012
036 : 014 9 036 9
⋮ 9
016 9
036 8
dan + ) adalah suku berorde tinggi yang bersifat lim)→> + ) = 0, dengan
)=− . *) pada sistem persamaan (8) disebut pelinearan sistem persamaan (7).
(Tu 1994)
2.5
Kestabilan Titik Tetap
Misalkan diberikan sistem persamaan diferensial biasa sebarang = ,
sebagai titik tetap. Kestabilan titik tetap dapat ditentukan ∈ ℝ dengan
dengan memperhatikan nilai-nilai eigen, yaitu ? , @ = 1, 2, … , , yang diperoleh
dari persamaan karakteristik. Secara umum, kestabilan titik tetap mempunyai perilaku: a. Re ? < 0, untuk setiap @, atau
1. Stabil, jika:
b. Terdapat ReFG H = 0 untuk sebarang I dan Re ? < 0 untuk setiap @ ≠ I.
2. Tidak stabil, jika terdapat paling sedikit satu @ sehingga Re ? > 0. 2.6
(Tu 1994)
Kriteria Routh-Hurwitz Kriteria Routh-Hurwitz ini digunakan ketika nilai eigen persamaan
karakteristik ! = K + L KM + ⋯ + LK = 0, maka didefinisikan N matriks karakteristik tidak dapat ditentukan dengan mudah. Jika diberikan persamaan
sebagai berikut:
O = PL Q, O = R
L LS
L 1 0 / L L L 1 . S T, ..., OG = . LU LV LS L ⋮ ⋮ . ⋮ -LGM LGM LGMS
L 1 0 / L L L . S OK = . LU LV LS ⋮ ⋮ . ⋮ 0 0 - 0
0 1 L ⋮ 0
⋯ 0 : ⋯ 0 9 ⋯ 09 ⋱ ⋮9 ⋯ LK 8
6
0 ⋯ 0 : 1 ⋯ 0 9 L ⋯ 0 9, ..., ⋮ ⋱ ⋮ 9 LGMV ⋯ LW 8
dengan syarat setiap unsur X, Y pada matriks OG adalah
L[M\ , untuk 0 < 2X − Y ≤ N ℎ[\ = ] 1 , untuk 2X = Y 0 , untuk 2X < Y atau 2X > N + Y.
stabil jika dan hanya jika det OG > 0, untuk Dengan demikian, titik tetap
setiap I = 1, 2, … , N.
Untuk N = 4 dan N = 5, kriteria Routh-Hurwitz diberikan berikut ini.
N = 4: L > 0, LS > 0, LV > 0, L L LS > LS + L LV .
N = 5: L? > 0; @ = 1, … , 5, L L LS > LS + L LV , dan L LV − LU L LV > L L LS − LS − LU L L − LS + L LU .
2.7
(Edelstein-Keshet 1988)
Bilangan Reproduksi Dasar
Bilangan reproduksi dasar, dinotasikan dengan ℜ> , merupakan suatu
ukuran potensi penyebaran penyakit dalam suatu populasi. Bilangan reproduksi dasar didefinisikan sebagai nilai harapan banyaknya populasi rentan yang menjadi terinfeksi selama masa infeksi berlangsung (van den Driessche dan Watmough 2008). Kondisi yang timbul adalah:
1. Jika ℜ> < 1, maka satu nyamuk terinfeksi akan menginfeksi kurang dari satu
manusia rentan atau satu manusia terinfeksi akan menginfeksi kurang dari satu nyamuk rentan, sehingga penyakit DBD akan hilang dari populasi.
7
2. Jika ℜ> > 1, maka satu nyamuk terinfeksi akan menginfeksi lebih dari satu manusia rentan atau satu manusia terinfeksi akan menginfeksi lebih dari satu nyamuk rentan, sehingga penyakit DBD akan bertahan di dalam populasi. ℜ> dalam penelitian ini ditentukan dari nilai eigen taknegatif dengan
modulus terbesar the next generation matrix. Matriks ini merupakan suatu matriks
model umum dengan i kompartemen penyakit dan j kompartemen tanpa yang dikonstruksi dari sub-subpopulasi yang menyebabkan infeksi saja. Untuk penyakit, nilai ℜ> dapat dihitung untuk setiap kompartemen.
Misalkan diberikan sistem persamaan diferensial taklinear = ,
∈ ℝ dan misalkan k ∈ ℝl dan ∈ ℝm adalah sub-subpopulasi pada setiap
kompartemen. Selanjutnya, dinotasikan no sebagai laju kenaikan infeksi pada
kompartemen penyakit ke-@ dan po sebagai laju pergerakan penyakit, kematian
dan penurunan kesembuhan dari kompartemen ke-@. Model kompartemen dapat ditulis sebagai
k o = no k, − po k, q = rq k,
, @ = 1,2, … , , I = 1,2, … , Y
maka sistem persamaan diferensial taklinear = , ∈ ℝ dapat ditulis
sebagai
k = s − tk
dengan s dan t adalah matriks-matriks berukuran × serta s=
0no 0uv
0, w> dan t =
0, w> adalah titik tetap tanpa penyakit.
0po
0uv
0, w> ;
The next generation matrix x untuk suatu sistem persamaan diferensial pada
titik tetap tanpa penyakit berbentuk
x = yz M .
Nilai eigen taknegatif dengan modulus terbesar matriks x, yaitu { yz M ,
yang nantinya digunakan sebagai nilai ℜ> , sehingga dapat ditulis { yz M = ℜ> .
(van den Driessche dan Watmough 2008)
8
III. MODEL MATEMATIK PENYEBARAN PENYAKIT DBD
3.1
Model SIR Model SIR pada uraian berikut mengacu pada kajian Derouich et al. (2003).
Asumsi yang digunakan adalah: 1. Total populasi nyamuk dan total populasi manusia adalah konstan. 2. Populasi manusia dan nyamuk adalah populasi yang tertutup.
Dari asumsi di atas, misalkan |} adalah populasi manusia dan |~ adalah
manusia rentan (susceptible) } , manusia terinfeksi (infected) } , dan manusia
populasi nyamuk. Populasi manusia dibagi menjadi tiga subpopulasi, yaitu sembuh (recovered) } . Populasi nyamuk dibagi menjadi dua subpopulasi, yaitu nyamuk rentan (susceptible) ~ dan nyamuk terinfeksi (infected) ~ .
Manusia rentan adalah manusia yang bukan imun dan belum tertular virus
dengue. Manusia terinfeksi adalah manusia yang telah tertular virus dan dapat menularkan virus tersebut. Manusia sembuh dianggap tidak dapat tertular lagi. Nyamuk rentan adalah nyamuk yang belum tertular virus. Nyamuk terinfeksi adalah nyamuk yang telah tertular virus dan dapat menularkan virus tersebut. Secara skematis, pola penyebaran penyakit DBD dapat digambarkan dalam diagram kompartemen berikut:
Keterangan :
Sv
/
Sh
Iv
/
Ih
Perpindahan Individu Pengaruh
Gambar 1 Skema penyebaran penyakit DBD model SIR.
Rh
9
Arti diagram kompartemen di atas adalah: 1. Laju pertumbuhan manusia rentan mempertimbangkan faktor kelahiran, kematian, fraksi acak manusia rentan yang terimunisasi dan proporsi perpindahan
dari
manusia
= } |} − } + i +
rentan
ke
manusia
terinfeksi,
} = } |} − } + i +
} ,
ditulis: dimana
diambil } = } . Proporsi perpindahan manusia rentan ke manusia terinfeksi dipengaruhi oleh peluang kontak antara nyamuk terinfeksi dengan manusia rentan (~} ). Nilai peluang ini ialah perkalian antara peluang transmisi virus
dari nyamuk terinfeksi ke manusia rentan (i~} ) dengan rata-rata gigitan
nyamuk terinfeksi (w? ). Jadi, ~} = i~} w? .
2. Laju pertumbuhan manusia terinfeksi mempertimbangkan faktor kematian, proporsi
perpindahan
manusia
rentan
ke
manusia
terinfeksi
dan
proporsi perpindahan manusia terinfeksi ke manusia sembuh, ditulis:
=
} − } + } } .
3. Laju pertumbuhan manusia sembuh mempertimbangkan faktor kematian, fraksi acak manusia rentan yang terimunisasi dan proporsi perpindahan manusia terinfeksi ke manusia sembuh, ditulis:
= i} + } } − } } .
4. Laju pertumbuhan nyamuk rentan mempertimbangkan faktor kelahiran, kematian dan proporsi perpindahan nyamuk rentan ke nyamuk terinfeksi,
ditulis:
= ~ |~ − ~ +
~ = ~ |~ − ~ +
~ ,
dimana
diambil ~ = ~ . Proporsi perpindahan nyamuk rentan ke nyamuk terinfeksi dipengaruhi oleh peluang kontak antara nyamuk rentan dengan manusia terinfeksi (}~ ). Nilai peluang ini ialah perkalian antara peluang transmisi
virus dari manusia terinfeksi ke nyamuk rentan (i}~ ) dengan rata-rata gigitan
nyamuk rentan (w ). Jadi, }~ = i}~ w .
5. Laju pertumbuhan nyamuk terinfeksi mempertimbangkan faktor kematian dan proporsi perpindahan nyamuk rentan ke nyamuk terinfeksi, ditulis:
=
~ − ~ ~ .
10
Berdasarkan uraian di atas, model SIR dinyatakan sebagai berikut: Populasi Manusia
Populasi Nyamuk
serta
= } − } + } } = i} + } } − } }
dengan kondisi
= } |} − } + i +
= ~ |~ − ~ +
=
~ − ~ ~
} (9)
~
(10)
} + } + } = |} dan ~ + ~ = |~
(11)
|} : total populasi manusia. |~ : total populasi nyamuk. } : laju kelahiran manusia. ~ : laju kelahiran nyamuk. } : laju kematian manusia. ~ : laju kematian nyamuk. i : fraksi acak manusia rentan yang terimunisasi. } : proporsi perpindahan manusia terinfeksi ke manusia sembuh. }~ : peluang terjadinya kontak antara nyamuk rentan dengan manusia terinfeksi. ~} : peluang terjadinya kontak antara nyamuk terinfeksi dengan manusia rentan.
Selanjutnya, sistem-sistem (9) dan (10) serta kondisi (11) dapat
disederhanakan dengan pemisalan } = , } = , } = , ~ = , dan
~ = , sehingga sistem tersebut dapat ditulis:
~ } = } − } + i + ~} = ~} ~ } − } + } } } ~ ~ = }~ 1 − − ~
(12)
} + } + } = 1 dan ~ + ~ = 1
(13)
dengan = , serta kondisi
11
Karena virus dengue membutuhkan masa inkubasi intrinsik dan ekstrinsik sebelum menyebar (Heymann 2008), maka model SIR ini dimodifikasi menjadi model SEIR. Modifikasi dilakukan dengan menambahkan tahap exposed. Pada tahap ini, manusia atau nyamuk rentan yang telah tertular virus menyelesaikan masa inkubasi intrinsik atau ekstrinsik sebelum terinfeksi. 3.2
Model SEIR
Pada model ini, populasi manusia |} dibagi menjadi empat subpopulasi,
yaitu manusia rentan (susceptible) } , manusia terpapar (exposed) } , manusia
terinfeksi (infected) } , dan manusia sembuh (recovered) } sedangkan populasi
nyamuk |~ dibagi menjadi tiga subpopulasi, yaitu nyamuk rentan (susceptible) ~ , nyamuk terpapar (exposed) ~ , dan nyamuk terinfeksi (infected) ~ . Asumsi yang digunakan ialah:
1. Total populasi nyamuk adalah konstan sedangkan total populasi manusia tidak konstan. 2. Populasi manusia dan nyamuk adalah populasi yang tertutup. Penularan virus dari nyamuk ke manusia terjadi melalui gigitan pada saat virus tersebut berada di kelenjar ludah nyamuk. Setelah itu, virus memerlukan 4-6 hari yang menunjukkan masa inkubasi intrinsik sebelum menimbulkan penyakit. Dalam masa inkubasi ini, manusia rentan dianggap telah terbuka untuk diinfeksi virus. Dengan demikian, manusia rentan tersebut selanjutnya dikelompokkan ke dalam subpopulasi manusia terpapar. Penularan virus dari manusia ke nyamuk hanya dapat terjadi jika nyamuk rentan menggigit manusia terinfeksi yang sedang mengalami viremia, yaitu suatu kondisi medis dimana virus Dengue berada di dalam darah manusia. Kondisi ini berlangsung selama 2 hari sebelum demam sampai 5 hari setelah demam. Selanjutnya, virus memerlukan 8-10 hari yang menunjukkan masa inkubasi ekstrinsik sebelum menimbulkan penyakit. Ketika masa inkubasi ini, nyamuk rentan dianggap telah terbuka untuk diinfeksi oleh virus. Nyamuk-nyamuk tersebut selanjutnya dikelompokkan ke dalam suatu subpopulasi nyamuk terpapar.
12
Secara skematis, pola penyebaran penyakit DBD dapat digambarkan dalam diagram kompartemen berikut: λ
λ
Sh
Sv
/
£¤
Ev
/
£¤
Eh
¥
Ih
Iv
£o
Rh
Keterangan : Perpindahan Individu Pengaruh
Gambar 2 Skema penyebaran penyakit DBD model SEIR. Arti diagram kompartemen di atas adalah: 1. Laju pertumbuhan manusia rentan mempertimbangkan faktor kelahiran, kematian dan proporsi perpindahan manusia rentan ke manusia terpapar, ditulis:
= } |} −
+ } } = |} −
+ } } ,
dimana
diambil = } . Proporsi perpindahan manusia rentan ke manusia terpapar dipengaruhi oleh peluang kontak antara nyamuk terinfeksi dengan manusia rentan (~} ). Nilai peluang ini ialah perkalian antara peluang transmisi virus
dari nyamuk terinfeksi ke manusia rentan (i~} ) dengan rata-rata gigitan
nyamuk terinfeksi (w? ). Jadi, ~} = i~} w? .
2. Laju pertumbuhan manusia terpapar mempertimbangkan faktor kematian, proporsi perpindahan manusia rentan ke manusia terpapar dan proporsi perpindahan
=
manusia
terpapar
} − ¡3} + } } .
ke
manusia
terinfeksi,
ditulis:
3. Laju pertumbuhan manusia terinfeksi mempertimbangkan faktor kematian baik kematian secara alami maupun kematian karena DBD, proporsi perpindahan manusia terpapar ke manusia terinfeksi dan proporsi perpindahan manusia terinfeksi ke manusia sembuh, ditulis:
= ¡3} } − ?} + ¢ + } } .
13
4. Laju pertumbuhan manusia sembuh mempertimbangkan faktor kematian dan proporsi perpindahan manusia terinfeksi ke manusia sembuh, ditulis:
= ?} } − } } .
5. Laju pertumbuhan nyamuk rentan mempertimbangkan faktor kelahiran, kematian dan proporsi perpindahan nyamuk rentan ke nyamuk terpapar, ditulis:
= ~ |~ −
+ ~ ~ = ~ |~ −
+ ~ ~ , dimana
diambil ~ = ~ . Proporsi perpindahan nyamuk rentan ke nyamuk terpapar dipengaruhi oleh peluang kontak antara nyamuk rentan dengan manusia terinfeksi (}~ ). Nilai peluang ini ialah perkalian antara peluang transmisi
virus dari manusia terinfeksi ke nyamuk rentan (i}~ ) dengan rata-rata gigitan
nyamuk rentan (w ). Jadi, }~ = i}~ w .
6. Laju pertumbuhan nyamuk terpapar mempertimbangkan faktor kematian, proporsi perpindahan adalah proporsi perpindahan nyamuk rentan ke nyamuk terpapar dan proporsi perpindahan nyamuk terpapar ke nyamuk nyamuk terinfeksi, ditulis:
=
~ − ¡3~ + ~ ~ .
7. Laju pertumbuhan nyamuk terinfeksi mempertimbangkan faktor kematian dan proporsi perpindahan nyamuk terpapar ke nyamuk nyamuk terinfeksi, ditulis:
= ¡3~ ~ − ~ ~
Berdasarkan uraian di atas, model SEIR dapat dinyatakan sebagai berikut: Populasi Manusia
Populasi Nyamuk
= |} − + } } = − } ¡3} + } } = ¡3} } − ?} + ¢ + } } = ?} } − } }
= ~ |~ − + ~ ~ = ~ − ¡3~ + ~ ~ = − ¡3~ ~ ~ ~
(14)
(15)
14
dengan kondisi
serta |} |~ ~ } ¢ ¡3} ¡3~ ?} }~ ~}
: : : : : : : : : : :
} + } + } + } = |} dan ~ + ~ + ~ = |~
(16)
total populasi manusia. total populasi nyamuk. laju kelahiran manusia laju kematian nyamuk. laju kematian manusia secara alami. laju kematian manusia karena DBD. proporsi perpindahan manusia terpapar ke manusia terinfeksi. proporsi perpindahan nyamuk terpapar ke nyamuk terinfeksi. proporsi perpindahan manusia terinfeksi ke manusia sembuh. peluang terjadinya kontak antara nyamuk rentan dengan manusia terinfeksi. peluang terjadinya kontak antara nyamuk terinfeksi dengan manusia rentan.
Selanjutnya, sistem-sistem (14) dan (15) serta kondisi (16) dapat
disederhanakan dengan pemisalan } = , } = , } = , } = , ~ = ,
~ = dan ~ = , dan juga dalam model ini dianggap bahwa nilai }~ =
~} = ¦, maka sistem tersebut dapat ditulis:
~ } = − ¦ + } = ¦~ } − ¡3} + } } = ¡3} } − ?} + ¢ + } } = ¦ } 1 − ~ − ~ − ¡3~ + ~ ~ ~ ~ = ¡3~ − ~
(17)
} + } + } + } = 1 dan ~ + ~ + ~ = 1
(18)
dengan =
serta kondisi
Sistem (17) dan kondisi (18) ini yang dibahas lebih lanjut pada bab berikut. Pembahasannya meliputi analisis kestabilan dan simulasi numerik untuk melihat dinamika populasinya.
15
IV. HASIL DAN PEMBAHASAN
4.1
Penentuan Titik Tetap
makna secara biologi, disebut Ω, dengan
Pada sub-bab ini dicari titik tetap sistem (17) pada daerah yang memiliki Ω = ¨ } , } , } , ~ , ~ ∈ ℝU© | } + } + } ≤ 1, ~ + ~ ≤ 1ª .
Titik tetap ini diperoleh dengan menyelesaikan sistem (17) tersebut. Solusinya merupakan suatu solusi yang diperoleh pada saat
=
= 0, sehingga sistem tersebut dapat ditulis
=
− ¦ ~ + } } = 0 ~ } } ¦ − ¡3} + } = 0 ¡3} } − ?} + ¢ + } } = 0 } ~ ~ ~ ¦ 1 − − − ¡3~ + ~ = 0 ¡3~ ~ − ~ ~ = 0
=
=
(19)
Sistem (19) di atas memiliki dua jenis titik tetap, yaitu titik tetap tanpa penyakit merupakan titik yang memuat nilai } = 0 dan ~ = 0, sedangkan titik penyakit (disease-free equilibrium/DFE) dan titik tetap endemik. Titik tetap tanpa tetap endemik merupakan titik yang memuat nilai } ≠ 0 atau ~ ≠ 0.
Dengan menggunakan software Mathematica, diperoleh titik tetap tanpa
penyakit
« } , } , } , ~ , ~ = «
dan titik tetap endemik
dengan
¬
, 0, 0, 0, 0
« ∗} , ∗} , ∗} , ∗~ , ∗~
(20)
(21)
∗} = ¨ ~ + ¡3~ P¦ ¡3} + ~ } + ¡3} ¢ + } + ?} Qª/¨¦ ¡3} P¦ ¡3~ + } ~ + ¡3~ Qª
∗} = −P−¦ ¡3} ¡3~ + }S ~ ~ + ¡3~ + ¢} ~ } + ¡3} ~ + ¡3~ + } ~ ¡3} ~ + ¡3~ ?} + } ~ ~ + ¡3~ ¡3} + ?} Q/¨¦ ¡3} } + ¡3} P¦ ¡3~ + } ~ + ¡3~ Qª
∗} = −} ~ ~ + ¡3~ + ¨ ¦ ¡3} ¡3~ ⁄P } + ¡3} ¢ + } + ?} Qª/ ¨¦P¦ ¡3~ + } ~ + ¡3~ Qª
16
∗~ = ¨~ P−} ¢ + } ~ } + ¡3} − P−¦ ¡3} + } ¢ + } ~ } + ¡3} Q ¡3~ − } ~ } + ¡3} ~ + ¡3~ ?} Qª/¨¦ ¡3~ ~ + ¡3~ P¦ ¡3} + ~ } + ¡3} ¢ + } + ?} Qª,
∗~ = 1/ −} /¦ + ¨ ¡3} P¦ ¡3~ + } ~ + ¡3~ Qª/¨ ~ + ¡3~ P¦ ¡3} + ~ } + ¡3} ¢ + } + ?} Qª. Penentuan titik-titik tetap di atas dapat dilihat pada Lampiran 1. 4.2
Penentuan Bilangan Reproduksi Dasar
Dengan menentukan the next generation matrix x untuk sistem (17) pada
titik tetap tanpa penyakit,
¦ ¡3~ / 0 } ~ ~ + ¡3~ . . ¦ + ~ ¡3} ~ ¡3} ¡3~ °=. 0 .~ } + ¡3} ~ + ¡3~ ¢ + } + ?} . 0 0 0 0
diperoleh bilangan reproduksi dasar ℜ> =
±√√¬³´µ¶ ³´µ¶
³ ©´µ¶ ©´µ¶ ·© ©´¸
0 ¦
¢ + } + ?} 0 0
.
¦ : } ~ 9 9 0 9 9 0 9 0 8
(22)
Selanjutnya, dari hasil (22) di atas diperoleh juga ¹ = ℜ> =
± 4 ¬´µ¶ ´µ¶
©´µ¶ ©´µ¶ ·© ©´¸
.
(23)
Penentuan ℜ> dan ¹ dapat dilihat pada Lampiran 2. 4.3
Analisis Kestabilan Titik Tetap
4.3.1 Kestabilan Titik Tetap Tanpa Penyakit Misalkan sistem (22) ditulis sebagai
} , } , } , ~ , ~ = − ¦ ~ + } } } } } ~ ~ ~ } } , , , , = ¦ − ¡3} + } S } , } , } , ~ , ~ = ¡3} } − ?} + ¢ + } } } } } ~ ~ } ~ ~ ~ V , , , , = ¦ 1 − − − ¡3~ + ~ U } , } , } , ~ , ~ = ¡3~ ~ − ~ ~ .
(24)
17
« } , } , } , ~ , ~ = « /} , 0, 0, 0, 0, digunakan pelinearan pada sistem Untuk
menentukan
kestabilan
(24) di sekitar « , diperoleh matriks Jacobi *º»
/ . =. . -
−} 0 0 0 0
0 −} − ¡3} ¡3} 0 0
titik
0 0 −¢ − } − ?} ¦ 0
tetap
0 0 0 −~ − ¡3~ ¡3~
tanpa
penyakit
−¦/} ¦/} : 9 0 9. 0 9 −~ 8
Dari matriks *º» di atas diperoleh lima nilai eigen. Nilai eigen yang pertama
adalah −} dan empat nilai eigen lainnya merupakan akar-akar persamaan karakteristik
¼ V + L ¼S + L ¼ + LS ¼ + LV = 0
dengan
L = ¢ + 2} + 2~ + ¡3} + ¡3~ + ?}
L = } + ~ + 2~ ¡3} + ~ ¡3~ + ¡3} ¡3~ + ¢ } + 2~ + ¡3} + ¡3~ + 2~ + ¡3} + ¡3~ ?} + } 4~ + ¡3} + 2 ¡3~ + ?}
LS = ¢~ 2} + ~ + 2 ¡3} + ¢ } + ~ + ¡3} ¡3~ + ~ ¡3} ~ + ¡3~ + } 2~ + ¡3~ + P~ ~ + 2 ¡3} + ~ + ¡3} ¡3~ Q ?} + } P2~ + ¡3~ ¡3} + ?} + 2~ ¡3} + ¡3~ + ?} Q
(25)
LV = ~ } + ¡3} ~ + ¡3~ ¢ + } + ?} 1 − ¹; ¹ ada pada persamaan (23).
berderajat 4, kondisi kestabilan sistem (17) pada titik tetap « adalah Berdasarkan
kriteria
Routh-Hurwitz
untuk
persamaan
L > 0, LS > 0, LV > 0 dan L L LS > LS + L LV .
karakteristik
(26)
Karena semua parameter bernilai positif, maka L , L dan LS pada (25)
bernilai positif. Koefisien LV dan L L LS − LS − L LV akan bernilai positif, nol atau negatif bergantung pada nilai ¹.
Jika ¹ < 1 maka LV > 0 dan L L LS − LS − L LV > 0. Jika ¹ = 1 maka
LV = 0 dan L L LS − LS − L LV = 0. Jika ¹ > 1 maka LV < 0 dan L L LS −
LS − L LV < 0. Jadi, kondisi (26) terpenuhi ketika ¹ < 1.
Dengan demikian, karena nilai eigen yang pertama −} < 0 dan kriteria
Routh-Hurwitz telah ditunjukkan terpenuhi, maka « stabil ketika ¹ < 1. Dari
18
persamaan (23), didapat hubungan ℜ> = ³¹, maka nilai ¹ dibatasi pada interval
P0, 1 sehingga, jika ℜ> = ³¹ < 1 untuk 0 ≤ ¹ < 1 maka titik tetap tanpa penyakit « adalah stabil. Sebaliknya, jika ℜ> = ³¹ > 1 maka titik tetap tanpa
penyakit « menjadi tidak stabil.
Dalam model ini, ℜ> = ³¹ adalah bilangan reproduksi dasar sebagaimana
yang diperlihatkan pada persamaan (22). Pelinearan, penentuan nilai eigen dan persamaan karakteristik serta bukti di atas dapat diperhatikan pada Lampiran 3. 4.3.2 Kestabilan Titik Tetap Endemik
Untuk menentukan kestabilan titik tetap endemik « ∗} , ∗} , ∗} , ∗~ , ∗~ ,
digunakan pelinearan pada sistem (24) di sekitar « , diperoleh matriks *º½
L /L . = .LS .LV -LU
L L LS LV LU
LS LS LSS LVS LUS
LV LV LSV LVV LUV
dengan L , L , … , LUU dapat dilihat pada Lampiran 4.
LU LU : LSU 99 LVU 9 LUU 8
Nilai eigen matriks *º½ merupakan akar-akar persamaan karakteristik ¼ U + L ¼ V + L ¼S + LS ¼ + LV ¼ + LU = 0
dengan L , L , LS , LV , dan LU dapat dilihat pada Lampiran 4.
derajat 5, kondisi kestabilan sistem (17) pada titik tetap « adalah
Berdasarkan kriteria Routh-Hurwitz untuk persamaan karakteristik berL > 0, L > 0, LS > 0, LV > 0, LU > 0, L L LS > LS + L LV , dan L LV − LU L L LS − LS − L LV > LU L L − LS + L LU .
(27)
Karena semua parameter bernilai positif, maka L bernilai positif. Koefisien
L , LS , LV , LU , L L LS − LS − L LV dan L LV − LU L L LS − LS − L LV − LU
L L − LS − L LU akan bernilai positif, nol atau negatif bergantung pada nilai ¹. Jika ¹ < 1 maka L < 0, LS < 0, LV < 0, LU < 0 dan L L LS LV + L LS LU +
2L LV LU − LS LV − L LV − L L LU − LU < 0. Jika ¹ = 1 maka L = 0, LS = 0,
19
LV = 0, LU = 0 dan L L LS LV + L LS LU + 2L LV LU − LS LV − L LV − L L LU −
LU = 0. Jika ¹ > 1 maka L > 0, LS > 0, LV > 0, LU > 0 dan L L LS LV + L LS LU
+2L LV LU − LS LV − L LV − L L LU − LU > 0. Jadi, kondisi (27) terpenuhi ketika ¹ > 1.
Dengan demikian, karena kriteria Routh-Hurwitz telah ditunjukkan
terpenuhi, maka « stabil ketika ℜ> = ³¹ > 1. Sebaliknya, jika ℜ> = ³¹ < 1
untuk 0 ≤ ¹ < 1 maka titik tetap endemik « menjadi tidak stabil.
Pelinearan, penentuan persamaan karakteristik dan bukti di atas dapat
diperhatikan pada Lampiran 4. Berikut ini adalah tabel kondisi kestabilan kedua titik tetap yang diperoleh. Tabel 1 Kondisi Kestabilan Titik Tetap Kondisi
ℜ> < 1 ℜ> > 1
Titik Tetap Tanpa Penyakit º» Stabil Tidak Stabil
Titik Tetap Endemik º½
Tidak Stabil Stabil
Dari Tabel 1 dapat dilihat bahwa kondisi kestabilan dari titik tetap yang diperoleh saling bertentangan. Ketika titik tetap pertama stabil, titik tetap kedua tidak stabil dan ketika titik tetap pertama tidak stabil, titik tetap kedua stabil. 4.4
Simulasi Dinamika Populasi Penularan Virus Dengue
4.4.1 Nilai Parameter
Dinamika populasi yang dianalisis adalah untuk kondisi ℜ> < 1 dan
ℜ> > 1. Dalam hal ini, ℜ> adalah bilangan reproduksi dasar (persamaan 22).
Untuk menganalisis dinamika populasi, dilakukan perubahan laju kematian nyamuk (~ ) dan rata-rata gigitan nyamuk terinfeksi (w? ). Dua parameter ini
dipilih karena dianggap berpengaruh dalam penanggulangan wabah.
Nilai ~ yang diambil berada pada [0.01, 0.09] dengan langkah 0.01,
sedangkan nilai w? yang diambil pada [0.25, 0.60] dengan langkah 0.01 (Hawley 1988; Vazeille et al. 2003 dan Richards et al. 2006). Nilai-nilai parameter lain dapat dilihat pada Tabel 2 berikut.
20
Tabel 2 Definisi dan Nilai Parameter Model SEIR dalam Simulasi Numerik Parameter
Notasi
Laju kelahiran manusia per hari
Nilai
}
0,003
2,244 × 10−5
i~}
Peluang transmisi virus dari nyamuk terinfeksi ke manusia per hari
¢
Laju kematian manusia karena DBD per hari Laju kematian manusia secara alami per hari
0,4 1/28 000
¡3}
Proporsi perpindahan manusia terpapar ke manusia terinfeksi per hari
1/10
¡3~
Proporsi perpindahan nyamuk terpapar ke nyamuk terinfeksi per hari
1/9
?}
Proporsi perpindahan manusia terinfeksi ke manusia sembuh per hari Sumber: Erickson et al. (2010) dan Derouich et al. (2003)
1/4
Nilai awal total populasi manusia yang seluruhnya rentan adalah 1. Nilai awal total populasi nyamuk adalah 1 dengan jumlah nyamuk terinfeksi adalah 20%. Berikut adalah simulasi untuk melihat dinamika populasi manusia dan
nyamuk yang dilakukan dengan mengubah nilai laju kematian nyamuk (~ ) dan
rata-rata gigitan nyamuk terinfeksi (w? ) pada kondisi ℜ> < 1 dan kondisi ℜ> > 1.
4.4.2 Dinamika Populasi untuk Kondisi ¾ < 1
populasi manusia maupun populasi nyamuk, untuk kondisi ℜ> < 1. Berdasarkan
Gambar 3 berikut menunjukkan kestabilan tiap subpopulasi, baik pada
nilai-nilai parameter yang ada pada Tabel 2 dan dengan mengambil nilai ~ dan w?
bawah ini untuk nilai ~ = 0,07 dan w? = 0,3 dengan nilai ℜ> = 0,56.
pada interval yang sudah ditetapkan, diperoleh gambar dinamika populasi di
1.0
1.0 0.8
0.8
Sh Ht L 0.6
Eh Ht L I h Ht L
0.4
Sv Ht L
0.6
Ev Ht L 0.4
I v Ht L
Rh Ht L 0.2 0.0
0.2
0
50
100
150
200
(a)
250
300
350
Hari t
0.0
0
50
100
150
200
250
300
350
Hari t
(b)
Gambar 3 Dinamika populasi manusia (a) dan populasi nyamuk (b) terhadap waktu untuk kondisi ℜ> < 1.
21
Gambar 3a menunjukkan bahwa jumlah subpopulasi manusia rentan ( } ) setelah tertular virus, dari awal simulasi mengalami penurunan hingga stabil ke } = 0,63. Lain halnya yang terjadi pada subpopulasi manusia terpapar ( } )
dan terinfeksi (} ), awalnya mengalami peningkatan kemudian menurun
hingga stabil ke } = 0 dan } = 0. Pada subpopulasi manusia sembuh ( } ), dari
awal simulasi mengalami peningkatan hingga stabil ke } = 1 − } + } +
} = 0,37.
Pada Gambar 3b, jumlah subpopulasi nyamuk terpapar ( ~ ), awalnya
mengalami peningkatan kemudian menurun hingga stabil ke ~ = 0. Lain halnya
yang terjadi pada subpopulasi nyamuk terinfeksi (~ ), dari awal simulasi mengalami penurunan hingga stabil ke ~ = 0. Pada subpopulasi nyamuk rentan ( ~ ), mengalami peningkatan hingga stabil ke ~ = 1 − ~ + ~ = 1.
titik tetap tanpa penyakit « } , } , } , ~ , ~ = « /} , 0, 0, 0, 0 dengan Dengan demikian, dapat dikatakan bahwa jumlah tiap subpopulasi stabil ke
/} = 0,63. Ini menunjukkan bahwa sub-subpopulasi manusia terpapar dan terinfeksi serta nyamuk terpapar dan terinfeksi menuju nol.
mengubah nilai parameter ~ dan w? . Pengambilan nilai kedua parameter ini Selanjutnya dilakukan simulasi pada populasi manusia dan nyamuk dengan
memenuhi kondisi ℜ> < 1, sehingga dapat disimulasikan untuk beberapa kondisi berbeda sebagaimana yang tertera pada Tabel 3 serta Gambar 4 dan 5 berikut. Tabel 3 Simulasi untuk Kondisi ℜ> < 1 Parameter Model
~ = 0,03
w? = 0,30
~ = 0,09
w? = 0,30
~ = 0,05 ~ = 0,07 ~ = 0,05 ~ = 0,05 ~ = 0,05 ~ = 0,05
w? = 0,30 w? = 0,30 w? = 0,25 w? = 0,30 w? = 0,35 w? = 0,40
Bilangan Reproduksi Dasar (¾ ) 0,97 0,70 0,56 0,47 0,59 0,70 0,82 0,94
22
Gambar 4 di bawah ini menunjukkan perubahan jumlah tiap subpopulasi
ketika nilai laju kematian nyamuk (~ ) diubah. S h Ht L 1.0
E h Ht L 0.14
0.9
0.12
0.8
mv = 0.03
0.10
mv = 0.03
0.7
mv = 0.05
0.08
mv = 0.05
mv = 0.07
0.06
mv = 0.07
mv = 0.09
0.04
mv = 0.09
0.6 0.5 0.4
0.02
0
10
20
30
40
Hari t
0.00
0
10
30
20
(a)
40
Hari t
(b)
I h Ht L 0.05
Rh Ht L 0.5
0.04
0.4 mv = 0.03
0.03
mv = 0.05 mv = 0.07
0.02
mv = 0.03
0.3
mv = 0.05 mv = 0.07
0.2
mv = 0.09
0.01 0.00
mv = 0.09
0.1
0
10
30
20
40
Hari t
0.0
0
10
20
(c)
40
Hari t
(d)
S v Ht L
E v Ht L
1.00
0.04
0.95
0.03 mv = 0.03
0.90
mv = 0.05
mv = 0.03 mv = 0.05
0.02
mv = 0.07
0.85
mv = 0.09
0.80 0.75
30
0
10
20
30
40
Hari t
mv = 0.07 mv = 0.09
0.01
0.00
0
10
(e)
20
30
40
Hari t
(f) I v Ht L 0.20
0.15 mv = 0.03 mv = 0.05
0.10
mv = 0.07 mv = 0.09
0.05
0.00
0
10
20
30
40
Hari t
(g)
Gambar 4
Dinamika populasi manusia (a) rentan } , (b) terpapar } , (c) terinfeksi } , dan (d) sembuh } , serta populasi nyamuk (e) rentan ~ , (f) terpapar ~ dan (g) terinfeksi ~ terhadap waktu pada kondisi ℜ> < 1 dan nilai parameter ~ diubah.
23
Pada populasi manusia sebagaimana ditunjukkan oleh Gambar 4a – 4d, jika laju kematian nyamuk naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi manusia rentan semakin bertambah sedangkan jumlah subpopulasi manusia lainnya semakin berkurang. Hal ini dikarenakan peningkatan laju kematian nyamuk menyebabkan penurunan pada jumlah nyamuk termasuk nyamuk terinfeksi. Akibatnya, proporsi perpindahan manusia rentan ke manusia terpapar semakin berkurang sehingga jumlah manusia rentan semakin bertambah. Pada populasi nyamuk sebagaimana ditunjukkan oleh Gambar 4e – 4g, jika laju kematian nyamuk naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi nyamuk rentan semakin bertambah sedangkan jumlah subpopulasi nyamuk lainnya semakin berkurang. Peningkatan laju kematian nyamuk ini menyebabkan penurunan pada jumlah nyamuk terinfeksi sehingga jumlah manusia terinfeksi pun semakin berkurang. Akibatnya, proporsi perpindahan nyamuk rentan ke nyamuk terpapar semakin berkurang sehingga jumlah nyamuk rentan semakin bertambah. Bertambah atau berkurangnya jumlah tiap subpopulasi cenderung tidak sama untuk setiap kenaikan laju kematian nyamuk, baik pada populasi manusia terjadi pada hari ke-15 dengan proporsi 12% dan laju kematian nyamuk sebesar maupun populasi nyamuk. Maksimum jumlah subpopulasi manusia terpapar
dengan proporsi 5% dan laju kematian nyamuk sebesar 0,03.
0,03. Pada subpopulasi manusia terinfeksi, maksimum terjadi pada hari ke-21
Gambar 5 berikut ini menunjukkan perubahan jumlah tiap subpopulasi
ketika nilai rata-rata gigitan nyamuk terinfeksi (w? ) diubah.
24 S h Ht L 1.0
E h Ht L 0.14
0.9
0.12
0.8
bi = 0.25
0.10
bi = 0.25
0.7
bi = 0.30
0.08
bi = 0.30
bi = 0.35
0.06
bi = 0.35
bi = 0.40
0.04
bi = 0.40
0.6 0.5 0.4
0.02
0
10
20
30
Hari t
40
0.00
0
10
30
20
(a)
40
Hari t
(b)
I h Ht L
Rh Ht L
0.06
0.4
0.05 0.3 0.04
bi = 0.25
0.03
bi = 0.30
bi = 0.25 bi = 0.30
0.2
bi = 0.35
0.02
bi = 0.40
bi = 0.35 bi = 0.40
0.1
0.01 0.00
0
10
30
20
40
Hari t
0.0
0
10
20
(c)
30
40
Hari t
(d)
S v Ht L
E v Ht L
1.00
0.05
0.95
0.04 bi = 0.25
0.90
bi = 0.25
0.03
bi = 0.30 bi = 0.35
0.85
bi = 0.30 bi = 0.35
0.02
bi = 0.40 0.80 0.75
bi = 0.40 0.01
0
10
20
30
40
Hari t
0.00
0
10
(e)
20
30
40
Hari t
(f) I v Ht L 0.20
bi = 0.25
0.15
bi = 0.30 bi = 0.35 0.10
0.05
bi = 0.40
10
20
30
40
Hari t
Dinamika populasi manusia (a) rentan } , (b) terpapar } , (c) terinfeksi } , dan (d) sembuh } , serta populasi nyamuk (e) rentan ~ , (f) terpapar ~ dan (g) terinfeksi ~ terhadap waktu pada kondisi ℜ> < 1 dan nilai parameter w? diubah. (g)
Gambar 5
Pada populasi manusia sebagaimana ditunjukkan oleh Gambar 5a – 5d, jika rata-rata gigitan nyamuk terinfeksi naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi manusia rentan semakin berkurang sedangkan jumlah
25
subpopulasi manusia lainnya semakin bertambah. Peningkatan rata-rata gigitan nyamuk terinfeksi dapat meningkatkan nilai peluang kontak antara nyamuk terinfeksi dengan manusia rentan. Akibatnya, proporsi perpindahan manusia rentan ke manusia terpapar semakin bertambah. Pada populasi nyamuk sebagaimana ditunjukkan oleh Gambar 5e – 5g, jika rata-rata gigitan nyamuk terinfeksi naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi nyamuk rentan semakin berkurang sedangkan jumlah subpopulasi nyamuk lainnya semakin bertambah. Hal ini disebabkan karena meningkatnya nilai peluang kontak antara nyamuk rentan dengan manusia terinfeksi sehingga proporsi perpindahan nyamuk rentan ke nyamuk terpapar semakin bertambah. Bertambah atau berkurangnya jumlah tiap subpopulasi cenderung tidak sama untuk setiap kenaikan rata-rata gigitan nyamuk terinfeksi, baik pada populasi terpapar terjadi pada hari ke-13 dengan proporsi 14% dan rata-rata gigitan nyamuk manusia maupun populasi nyamuk. Maksimum jumlah subpopulasi manusia
hari ke-18 dengan proporsi 5% dan laju kematian nyamuk sebesar 0,4.
terinfeksi sebesar 0,4. Pada subpopulasi manusia terinfeksi, maksimum terjadi pada 4.4.3 Dinamika Populasi untuk Kondisi ¾ > 1
populasi manusia maupun populasi nyamuk, untuk kondisi ℜ> > 1. Berdasarkan Gambar 6 berikut menunjukkan kestabilan tiap subpopulasi, baik pada
nilai-nilai parameter yang ada pada Tabel 2 dan dengan mengambil nilai ~ dan w?
ini untuk nilai ~ = 0,01 dan w? = 0,41 dengan ℜ> = 2,47.
pada interval yang sudah ditetapkan, diperoleh gambar dinamika populasi berikut
1.0
1.0
0.8
0.8
Sh Ht L 0.6
Eh Ht L I h Ht L
0.4
Sv Ht L
0.6
Ev Ht L 0.4
I v Ht L
Rh Ht L 0.2 0.0
0.2
0
50
100
150
200
(a)
250
300
350
Hari t
0.0
0
50
100
150
200
250
300
350
Hari t
(b)
Gambar 6 Dinamika populasi manusia (a) dan populasi nyamuk (b) terhadap waktu untuk kondisi ℜ> > 1.
26
Gambar 6a menunjukkan bahwa jumlah subpopulasi manusia rentan ( } ) setelah tertular virus, dari awal simulasi mengalami penurunan hingga stabil ke ∗} = 0,00248. Lain halnya yang terjadi pada subpopulasi manusia terpapar ( } ) dan terinfeksi ( } ), awalnya mengalami peningkatan kemudian menurun hingga
stabil ke ∗} = 0,00010 dan ∗} = 0,00004. Pada subpopulasi manusia sembuh ( } ),
∗} = 1 − ∗} + ∗} + ∗} = 0,99738.
di awal simulasi mengalami peningkatan kemudian menurun hingga stabil ke Pada Gambar 6b, jumlah subpopulasi nyamuk terpapar ( ~ ), awalnya
mengalami peningkatan kemudian menurun hingga stabil ke ∗~ = 0,00006. Lain
halnya yang terjadi pada subpopulasi nyamuk terinfeksi ( ~ ), dari awal simulasi mengalami penurunan hingga stabil ke ∗~ = 0,02357. Pada subpopulasi nyamuk rentan ( ~ ), mengalami peningkatan hingga stabil ke ∗~ = 1 − ∗~ + ∗~ = 0,97637.
Jadi, dapat dikatakan bahwa jumlah tiap subpopulasi stabil ke titik tetap
endemik « ∗} , ∗} , ∗} , ∗~ , ∗~ . Ini menunjukkan bahwa sub-subpopulasi manusia terpapar dan terinfeksi serta nyamuk terpapar dan terinfeksi menuju ke nilai yang tidak nol.
terhadap waktu dengan mengubah nilai parameter ~ dan w? . Pengambilan nilai Selanjutnya dilakukan simulasi pada tiap subpopulasi manusia dan nyamuk
kedua parameter ini disesuaikan dengan nilai ℜ> yang memenuhi kondisi ℜ> > 1,
sehingga dapat disimulasikan untuk tiga kondisi berbeda sebagaimana yang tertera pada Tabel 4 serta Gambar 7 dan 8 berikut. Tabel 4 Simulasi untuk Kondisi ℜ> > 1 Parameter Model
~ = 0,01
w? = 0,60
~ = 0,07
w? = 0,60
~ = 0,03 ~ = 0,05 ~ = 0,03 ~ = 0,03 ~ = 0,03 ~ = 0,03
w? = 0,60 w? = 0,60 w? = 0,40 w? = 0,45 w? = 0,50 w? = 0,55
Bilangan Reproduksi Dasar (¾ ) 3,62 1,94 1,40 1,12 1,29 1,45 1,61 1,78
27 S h Ht L 1.0
E h Ht L 0.20
0.8
0.15 mv = 0.01
0.6
mv = 0.03
mv = 0.03
mv = 0.05
0.4
mv = 0.07
0.2 0.0
mv = 0.01
0.10
0
10
20
30
40
Hari t
mv = 0.05 mv = 0.07
0.05
0.00
0
10
20
30
(a)
40
Hari t
(b) Rh Ht L
I h Ht L
0.7
0.08
0.6 0.06
0.5
mv = 0.01
0.04
0.02
mv = 0.01
0.4
mv = 0.03
mv = 0.05
0.3
mv = 0.05
mv = 0.07
0.2
mv = 0.07
mv = 0.03
0.1 0.00
0
10
20
30
40
0.0
Hari t
0
10
20
(c)
30
40
Hari t
(d)
S v Ht L
E v Ht L
1.0
0.07 0.06
0.9
0.05 mv = 0.01 mv = 0.03
0.8
0.7
mv = 0.01
0.04
mv = 0.03
mv = 0.05
0.03
mv = 0.05
mv = 0.07
0.02
mv = 0.07
0.01 0.6
0
10
20
30
40
0.00
Hari t
0
10
20
(e)
30
40
Hari t
(f) I v Ht L 0.35 0.30 0.25
mv = 0.01
0.20
mv = 0.03
0.15
mv = 0.05
0.10
mv = 0.07
0.05 0.00
0
10
20
30
40
Hari t
Dinamika populasi manusia (a) rentan } , (b) terpapar } , (c) terinfeksi } , dan (d) sembuh } , serta populasi nyamuk (e) rentan ~ , (f) terpapar ~ dan (g) terinfeksi ~ terhadap waktu pada kondisi ℜ> > 1 dan nilai parameter ~ diubah (g)
Gambar 7
Perubahan (bertambah atau berkurang) jumlah tiap subpopulasi, baik pada populasi manusia maupun populasi nyamuk, karena naiknya laju kematian nyamuk sebagaimana yang ditunjukkan dalam Gambar 7 di atas memiliki pola
28
yang sama dengan perubahan jumlah tiap subpopulasi untuk kondisi ℜ> < 1
sebagaimana yang ditunjukkan pada Gambar 4. Perbedaan kedua gambar hanya pada jumlah maksimum atau minimum tiap subpopulasi, tetapi dalam simulasi ini lebih difokuskan pada jumlah maksimum subpopulasi manusia terpapar dan manusia terinfeksi. Pada populasi manusia sebagaimana ditunjukkan oleh Gambar 7a – 7d, jika laju kematian nyamuk naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi manusia rentan semakin bertambah sedangkan jumlah subpopulasi manusia lainnya semakin berkurang. Hal ini dikarenakan peningkatan laju kematian nyamuk menyebabkan penurunan pada jumlah nyamuk termasuk nyamuk terinfeksi. Akibatnya, proporsi perpindahan manusia rentan ke manusia terpapar semakin berkurang. Pada populasi nyamuk sebagaimana ditunjukkan oleh Gambar 7e – 7g, jika laju kematian nyamuk naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi nyamuk rentan semakin bertambah sedangkan jumlah subpopulasi nyamuk lainnya semakin berkurang. Peningkatan laju kematian nyamuk ini menyebabkan penurunan pada jumlah nyamuk terinfeksi sehingga jumlah manusia terinfeksi pun semakin berkurang. Berkurangnya manusia terinfeksi menyebabkan proporsi perpindahan nyamuk rentan ke nyamuk terpapar semakin berkurang sehingga jumlah nyamuk rentan semakin bertambah. Bertambah atau berkurangnya jumlah tiap subpopulasi cenderung tidak sama untuk setiap kenaikan laju kematian nyamuk, baik pada populasi manusia terjadi pada hari ke-19 dengan proporsi 18% dan laju kematian nyamuk sebesar maupun populasi nyamuk. Maksimum jumlah subpopulasi manusia terpapar
dengan proporsi 7% dan laju kematian nyamuk sebesar 0,01.
0,01. Pada subpopulasi manusia terinfeksi, maksimum terjadi pada hari ke-24
Gambar 8 berikut ini menunjukkan perubahan jumlah tiap subpopulasi
ketika nilai rata-rata gigitan nyamuk terinfeksi (w? ) diubah.
29 S h Ht L 1.0
E h Ht L 0.25
0.8
0.20 bi = 0.40
0.6
bi = 0.45 bi = 0.50
0.4
bi = 0.40
0.15
bi = 0.45
bi = 0.50
0.10
bi = 0.55
0.2 0.0
bi = 0.55 0.05
0
10
20
30
40
Hari t
0.00
0
10
20
(a)
30
40
Hari t
(b)
I h Ht L
Rh Ht L
0.08
0.7 0.6
0.06
0.5
bi = 0.40 0.04
0.02
bi = 0.40
bi = 0.45
0.4
bi = 0.45
bi = 0.50
0.3
bi = 0.50
bi = 0.55
0.2
bi = 0.55
0.1 0.00
0
10
20
30
40
Hari t
0.0
0
10
20
(c)
30
40
Hari t
(d)
S v Ht L
E v Ht L
0.85 0.08 0.80
bi = 0.40
0.75
0.06
bi = 0.40
bi = 0.45
bi = 0.50
0.70
bi = 0.55 0.65 0.60
0
10
20
30
40
Hari t
bi = 0.45 0.04
bi = 0.50 bi = 0.55
0.02 0.00
0
10
20
(e)
30
40
Hari t
(f) I v Ht L 0.30
0.25 bi = 0.40 bi = 0.45
0.20
bi = 0.50 bi = 0.55
0.15
0.10
0
10
20
30
40
Hari t
Dinamika populasi manusia (a) rentan } , (b) terpapar } , (c) terinfeksi } , dan (d) sembuh } , serta populasi nyamuk (e) rentan ~ , (f) terpapar ~ dan (g) terinfeksi ~ terhadap waktu pada kondisi ℜ> > 1 dan nilai parameter w? diubah (g)
Gambar 8
Seperti uraian pada Gambar 7 sebelumnya, perubahan (bertambah atau berkurang) jumlah tiap subpopulasi, baik pada populasi manusia maupun populasi nyamuk, karena naiknya rata-rata gigitan nyamuk terinfeksi (Gambar 8), memiliki
30
pola yang sama dengan perubahan jumlah tiap subpopulasi untuk kondisi ℜ> < 1
(Gambar 5). Perbedaan kedua gambar hanya pada jumlah maksimum atau minimum tiap subpopulasi, tetapi dalam simulasi ini lebih difokuskan pada jumlah maksimum subpopulasi manusia terpapar dan manusia terinfeksi. Pada populasi manusia sebagaimana ditunjukkan oleh Gambar 8a – 8d, jika rata-rata gigitan nyamuk terinfeksi naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi manusia rentan semakin berkurang sedangkan jumlah subpopulasi manusia lainnya semakin bertambah. Peningkatan rata-rata gigitan nyamuk terinfeksi dapat meningkatkan nilai peluang kontak antara nyamuk terinfeksi dengan manusia rentan. Akibatnya, proporsi perpindahan manusia rentan ke manusia terpapar semakin bertambah. Pada populasi nyamuk sebagaimana ditunjukkan oleh Gambar 8e – 8g, jika rata-rata gigitan nyamuk terinfeksi naik dan nilai parameter lainnya tetap, maka jumlah subpopulasi nyamuk rentan semakin berkurang sedangkan jumlah subpopulasi nyamuk lainnya semakin bertambah. Hal ini disebabkan karena meningkatkan nilai peluang kontak antara nyamuk rentan dengan manusia terinfeksi sehingga proporsi perpindahan nyamuk rentan ke nyamuk terpapar semakin bertambah. Bertambah atau berkurangnya jumlah tiap subpopulasi cenderung tidak sama untuk setiap kenaikan rata-rata gigitan nyamuk terinfeksi, baik pada manusia terpapar terjadi pada hari ke-15 dengan proporsi 20% dan rata-rata
populasi manusia maupun populasi nyamuk. Maksimum jumlah subpopulasi
maksimum terjadi pada hari ke-21 dengan proporsi 8% dan laju kematian nyamuk
gigitan nyamuk terinfeksi sebesar 0,55. Pada subpopulasi manusia terinfeksi,
sebesar 0,55.
31
V. SIMPULAN DAN SARAN
5.1
Simpulan Secara umum model yang dihasilkan dapat menunjukkan adanya endemik di
suatu daerah untuk nilai parameter tertentu. Hal ini dapat dilihat dari perhitungan titik tetap model SEIR. Rincian hasil-hasil utama dalam tesis ini dapat disarikan pada uraian berikut. tetap tanpa penyakit « } , } , } , ~ , ~ = « /} , 0, 0, 0, 0 yang selalu ada Berdasarkan hasil dan pembahasan diperoleh 2 titik tetap. Pertama, titik
dan merupakan titik yang stabil apabila nilai bilangan reproduksi dasar ℜ> < 1.
Kedua, titik tetap endemik « ∗} , ∗} , ∗} , ∗~ , ∗~ dengan ∗} adalah jumlah subpopulasi manusia rentan, ∗} adalah jumlah subpopulasi manusia terpapar, ∗}
adalah jumlah subpopulasi manusia terinfeksi, ∗~ adalah jumlah subpopulasi nyamuk terpapar, ∗~ adalah jumlah subpopulasi nyamuk terinfeksi. Kestabilan titik tetap endemik ini dijamin apabila nilai ℜ> > 1.
Selanjutnya, melalui pengamatan yang dilakukan secara simulasi numerik
pemilihan nilai ℜ> . Dalam tesis ini, nilai ℜ> dipengaruhi oleh beberapa nilai
didapat hasil dinamik untuk masing-masing subpopulasi dipengaruhi oleh parameter, tetapi yang menjadi fokus simulasi adalah laju kematian nyamuk (~ )
dan rata-rata gigitan nyamuk terinfeksi (w? ).
Informasi tentang pengaruh ~ dan w? diuraikan sebagai berikut:
titik tetap tanpa penyakit « ketika ℜ> < 1 dan stabil ke titik tetap endemik «
1. Jumlah tiap subpopulasi pada populasi manusia dan nyamuk akan stabil ke ketika ℜ> > 1.
2. Pada populasi manusia: a. Semakin besar laju kematian nyamuk maka jumlah manusia rentan yang menjadi terpapar semakin sedikit. b. Semakin besar rata-rata gigitan nyamuk terinfeksi maka jumlah manusia rentan yang menjadi terpapar semakin banyak.
32
3. Pada populasi nyamuk: a. Semakin besar laju kematian nyamuk maka jumlah nyamuk rentan yang menjadi terpapar semakin sedikit. b. Semakin besar rata-rata gigitan nyamuk terinfeksi maka jumlah nyamuk rentan yang menjadi terpapar semakin banyak. 4. Bertambah atau berkurangnya jumlah tiap subpopulasi cenderung tidak sama untuk setiap kenaikan laju kematian nyamuk atau untuk setiap kenaikan ratarata gigitan nyamuk terinfeksi, baik pada populasi manusia maupun populasi nyamuk. 5.2
Saran Karena kestabilan sistem model SEIR terjadi setelah waktu yang relatif
cukup lama untuk beberapa nilai parameter, khususnya pada subpopulasi manusia rentan, maka perlu dilakukan pengembangan model SEIR ini.
33
DAFTAR PUSTAKA
Anonym. 2009. World Health Organization: Dengue and Dengue Haemorrhagic Fever. Geneva. Anton H. 1995. Aljabar Linear Elementer (Edisi ke-5). Terjemahan Pantur Silaban dan I Nyoman Susila. Jakarta. Erlangga Braun M. 1983. Differential Equations and Their Applications. New York. Springer-Verlag. Derouich M, Boutayeb A, Twizell EH. 2003. A Model of Dengue Fever. BioMedical Engineering OnLine 2 (4). Edelstein-Keshet L. 1988. Mathematical Models in Biology. New York. Random House Erickson RA, Presley SM, Allen LJS, Long KR, Cox SB. 2010. A Dengue Model with A Dynamic Aedes albopictus Vector Population. Ecological Modelling 221, 2899–2908. Estrada-Franco JGGB, Craig J. 1995. Biology, Disease Relationships, and Control of Aedes albopictus. Technical report. Pan American Health Organization. Gratz NG. 2004. Critical Review of The Vector Status of Aedes albopictus. Medical and Veterinary Entomology 18 (3), 215–227. Gubler DJ. 1998. Dengue and Dengue Hemorrhagic Fever. Clinical Microbiology Reviews 11 (3), 480–496. Hawley WA. 1988. The Biology of Aedes albopictus. Journal of the American Mosquito Control Association 4, 1–39. Heymann DL. 2008. Control of Communicable Diseases Manual, 18th edition. American Public Health Association, Washington, DC. Richards SL, Ponnusamy L, Unnasch TR, Hassan HK, Apperson CS. 2006. HostFeeding Patterns of Aedes albopictus (Diptera: Culicidae) in Relation to Availability of Human and Domestic Animals in Suburban Landscapes of Central North Carolina. Journal of Medical Entomology 43 (3), 543–551.
34
Tu PNV. 1994. Dinamical System, An Introduction with Applications in Economics and Biology. New York: Springer-Verlag van den Driessche P, Watmough J. 2008. Chapter 6: Further Notes on the Basic Reproduction Number. In: Brauer, F., van den Driessche, P., Wu, J. (Eds.), Mathematical Epidemiology, Vol. 1945. Lecture Notes in Mathematics, Springer, pp. 159–178. Vazeille M, Mousson L, Failloux AB, 2003. Low Oral Receptivity for Dengue Type 2 Viruses of Aedes ablbopictus from Southeast Asia Compared with That of Aedes aegypti. American Journal of Tropical Medicine and Hygiene 68 (2), 203–208. Watmough J. 2008. Computation of The Basic Reproduction Number. MITACSPIMS Summer School on Mathematical Modelling of Infectious Disease. University of Alberta
35
LAMPIRAN
Lampiran 1.
Penentuan Titik Tetap
Clear@λ, α, µh, µv, τexh, τih, τexv, c, bi, pvh, n, nh, nv, sh, eh, ih, ev, ivD H∗ Definisikan sistem persamaan diferensial 33 ∗L f1@sh_, eh_, ih_, ev_, iv_D := λ − Hn ∗ c ∗ iv + µhL ∗ sh; f2@sh_, eh_, ih_, ev_, iv_D := n ∗ c ∗ iv ∗ sh − Hτexh + µhL ∗ eh; f3@sh_, eh_, ih_, ev_, iv_D := τexh ∗ eh − Hτih + α + µhL ∗ ih; f4@sh_, eh_, ih_, ev_, iv_D := c ∗ ih ∗ H1 − ev − ivL − Hτexv + µvL ∗ ev; f5@sh_, eh_, ih_, ev_, iv_D := τexv ∗ ev − µv ∗ iv; sol = FullSimplify@Solve@ 8f1@sh, eh, ih, ev, ivD 0, f2@sh, eh, ih, ev, ivD 0, f3@sh, eh, ih, ev, ivD 0, f4@sh, eh, ih, ev, ivD 0, f5@sh, eh, ih, ev, ivD 0<, 8sh, eh, ih, ev, iv
:sh → HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê Hc τexh Hc n τexv + µh Hµv + τexvLLL, eh → − I− c2 n λ τexh τexv + µh3 µv Hµv + τexvL + α µh µv Hµh + τexhL Hµv + τexvL + µh µv τexh Hµv + τexvL τih + µh2 µv Hµv + τexvL Hτexh + τihLM ë Hc τexh Hµh + τexhL Hc n τexv + µh Hµv + τexvLLL, ih → I− µh µv Hµv + τexvL + Ic2 n λ τexh τexvM ë HHµh + τexhL Hα + µh + τihLLM ë
Hc Hc n τexv + µh Hµv + τexvLLL, ev → Iµv I− µh Hα + µhL µv2 Hµh + τexhL −
I− c2 n λ τexh + µh Hα + µhL µv Hµh + τexhLM τexv −
µh µv Hµh + τexhL Hµv + τexvL τihMM ë
Hc n τexv Hµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL, 1 iv → H− µh ê c + Hλ τexh Hc n τexv + µh Hµv + τexvLLL ê n HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLL>>
Lampiran 2.
Penentuan Bilangan Reproduksi Dasar H√0 )
Dengan mengambil sistem persamaan (22) dan menyusunnya kembali dalam urutan sub-sub popuasi yang menyebabkan infeksi saja, yaitu Eh , Ev , Ih , dan Iv diperoleh sistem persamaan sebagai berikut: dEh ë dt dEv ê dt dIh ë dt dIv ê dt
= = = =
ncIv Sh − Hτexh + µh L Eh cIh H1 − Ev − Iv L − Hτexv + µv L Ev τexh Eh − Hτih + α + µh L Ih τexv Ev − µv Iv
Dari sistem di atas, diperoleh matriks-matriks ncIv Sh cIh H1 − Ev − Iv L = 0 0
Hτexh + µh L Eh Hτexv + µv L Ev dan = − τexh Eh + Hτih + α + µh L Ih − τexv Ev + µv Iv
sehingga, 0 0 0
F= 0 0 c 0 0 0 0 0 0
ncλ µh
0 0 0
τexh + µh 0 0 0 0 τexv + µv 0 0 − τexh 0 τih + α + µh 0 0 − τexv 0 µv
dan V =
Selanjutnya, dihitung matriks K = FV-1 sebagai berikut:
K=
0
c n λ τexv µh µv Hµv +τexv L
0
cnλ µh µv
c Hµ2v τexh +µv τexh τexv L µv Hµh +τexh L Hµv +τexv L Hα+µh +τih L
0
c α+µh +τih
0
0 0
0 0
0 0
0 0
sehingga bilangan reproduksi dasar adalah √0 =
c
n
λ
τexh
τexv
µh µv Hµh +τexh L Hµv +τexv L Hα+µh +τih L
dan
x = √0 2 =
c2 n λ τexh τexv µh µv Hµh +τexh L Hµv +τexv L Hα+µh +τih L
Lampiran 3.
Analisis Kestabilan Titik Tetap Tanpa Penyakit
Clear@jD H∗ Definisikan Matriks Jacobi ∗L j = D@8f1@sh, eh, ih, ev, ivD, f2@sh, eh, ih, ev, ivD, f3@sh, eh, ih, ev, ivD, f4@sh, eh, ih, ev, ivD, f5@sh, eh, ih, ev, ivD<, 88sh, eh, ih, ev, iv<
1 a1 = FullSimplify@CoefficientList@perskar1, ψD@@4DDD α + 2 µh + 2 µv + τexh + τexv + τih
a2 = FullSimplify@CoefficientList@perskar1, ψD@@3DDD µh2 + µv2 + 2 µv τexh + µv τexv + τexh τexv + α Hµh + 2 µv + τexh + τexvL + H2 µv + τexh + τexvL τih + µh H4 µv + τexh + 2 τexv + τihL a3 = FullSimplify@CoefficientList@perskar1, ψD@@2DDD α µv H2 µh + µv + 2 τexhL + α Hµh + µv + τexhL τexv + µv τexh Hµv + τexvL + µh2 H2 µv + τexvL + Hµv Hµv + 2 τexhL + Hµv + τexhL τexvL τih + µh I2 µv2 + τexv Hτexh + τihL + 2 µv Hτexh + τexv + τihLM a4 = FullSimplify@CoefficientList@perskar1, ψD@@1DDD 1 ê µh I− c2 n λ τexh τexv + µh3 µv Hµv + τexvL + α µh µv Hµh + τexhL Hµv + τexvL + µh µv τexh Hµv + τexvL τih + µh2 µv Hµv + τexvL Hτexh + τihLM
Selanjutnya dilakukan penyederhanaan terhadap a4 sehingga diperoleh a4 = 1 ê µh I− c2 n λ τexh τexv + µh3 µv Hµv + τexvL + α µh µv Hµh + τexhL Hµv + τexvL + µh µv τexh Hµv + τexvL τih + µh2 µv Hµv + τexvL Hτexh + τihLM; = 1 ê µh I− c2 n λ τexh τexv + µh µv Hµh + τexhL Hµv + τexvL Hα + µh + τihLM; = µv Hµh + τexhL Hµv + τexvL Hα + µh + τihL H1 − ξL; ξ = Ic2 n λ τexh τexvM ë Hµh µv Hµh + τexhL Hµv + τexvL Hα + µh + τihLL
H∗ Masukan Nilai−nilai Parameter ∗L λ = 2.244 ∗ 10−5 ; α = 0.003; µh = 1 ê 28 000; µv = 1 ê 20; τexh = 1 ê 10; τih = 1 ê 4; τexv = 1 ê 9; bi = 0.3; pvh = 0.4; c = bi pvh; nh = 1; nv = 1; n = nh ê nv; H∗ Definisikan ξ ∗L Clear@ξD ξ = Ic2 n λ τexh τexvM ë Hµh µv Hµv + τexvL Hµh + τexhL Hµh + α + τihLL; H∗ Hitung a4 yang Belum Disederhanakan ∗L a4 H∗ Hitung a4 yang Sudah Disederhanakan ∗L Clear@a4D; a4 = µv Hµh + τexhL Hµv + τexvL Hα + µh + τihL H1 − ξL
0.000103376 0.000103376 Clear@ξ, a4D H∗ Ambil sebarang nilai ξ < 1 ∗L ξ = 0.5; H∗ Definisikan Koefisien a4 yang Sudah Disederhanakan ∗L a4 = µv Hµh + τexhL Hµv + τexvL Hα + µh + τihL H1 − ξL; H∗ Pembuktian Kriteria Routh−Hurwitz ∗L a1 > 0 a2 > 0 a3 > 0 a4 > 0 a1 a2 a3 > a32 + a12 a4 True True True True True
Lampiran 4.
Analisis Ketabilan Titik Tetap Endemik
Clear@jt2, brs, perskar2, λ, α, µh, µv, τexh, τih, τexv, c, bi, pvh, n, nh, nvD H∗ Matriks Jacobi untuk Titik Tetap Tanpa Penyakit ∗L jt2 = FullSimplify@j ê. sol@@2DDD; H∗ Unsur−unsur Baris Ke−p Matriks JT2 ∗L brs@p_ D := jt2@@pDD H∗ Unsur−unsur Baris Ke−1 Matriks JT2 ∗L brs@1D 8− Hc λ τexh Hc n τexv + µh Hµv + τexvLLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL, 0, 0, 0, − Hn Hµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê Hτexh Hc n τexv + µh Hµv + τexvLLL< H∗ Unsur−unsur Baris Ke−2 Matriks JT2 ∗L brs@2D 8c H− µh ê c + Hλ τexh Hc n τexv + µh Hµv + τexvLLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLL, − µh − τexh, 0, 0, Hn Hµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê Hτexh Hc n τexv + µh Hµv + τexvLLL< H∗ Unsur−unsur Baris Ke−3 Matriks JT2 ∗L brs@3D 80, τexh, − α − µh − τih, 0, 0< H∗ Unsur−unsur Baris Ke−4 Matriks JT2 ∗L brs@4D 90, 0, Hµv Hµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLL ê Hn τexv Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL, τexv H− 1 − Hc n Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê HHµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLLL, − I− µh µv Hµv + τexvL + Ic2 n λ τexh τexvM ë HHµh + τexhL Hα + µh + τihLLM ë Hc n τexv + µh Hµv + τexvLL= H∗ Unsur−unsur Baris Ke−5 Matriks JT2 ∗L brs@5D 80, 0, 0, τexv, − µv<
H∗ Ruas Kiri Persamaan Karakteristik ∗L perskar2 = Collect@HDet@jt2 − ψ IdentityMatrix@5DDL ∗ H− 1L, ψD; H∗ Koefisien Persamaan Karakteristik ∗L Clear@a0, a1, a2, a3, a4, a5D a0 = CoefficientList@perskar2, ψD@@6DD; a1 = CoefficientList@perskar2, ψD@@5DD; a2 = CoefficientList@perskar2, ψD@@4DD; a3 = CoefficientList@perskar2, ψD@@3DD; a4 = CoefficientList@perskar2, ψD@@2DD; a5 = CoefficientList@perskar2, ψD@@1DD;
Selanjutnya dilakukan penyederhanaan terhadap a2 , a3 , a4 , dan a5 sehingga diperoleh a2 = α µh + µh2 + α µv + 2 µh µv + α τexh + µh τexh + µv τexh + µh τih + µv τih + τexh τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLL Hα + 2 µh + µv + τexh + τihLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + τexv Hα + 2 µh + µv + τexh + τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLL H1 + Hc n Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê HHµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLLL + Hµh µv τexv Hµv + τexvLL ê Hc n τexv + µh Hµv + τexvLL H− 1 + ξL; a3 = α µh µv + µh2 µv + α µv τexh + µh µv τexh + µh µv τih + µv τexh τih + Ic λ τexh Hc n τexv + µh Hµv + τexvLL Iµh2 + µv τexh + α Hµh + µv + τexhL + Hµv + τexhL τih + µh H2 µv + τexh + τihLMM ë HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + τexv Iα µh + µh2 + α µv + 2 µh µv + α τexh + µh τexh + µv τexh + µh τih + µv τih + τexh τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLL τihL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLM H1 + Hc n Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê HHµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLLL + Hµh µv τexv Hµv + τexvL HHα + 2 µh + τexh + τihL ê Hc n τexv + µh Hµv + τexvLL + Hc λ τexhL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLLL H− 1 + ξL; a4 = Ic λ τexh Iµv Hµh + τexhL2 Hc n τexv + µh Hµv + τexvLL Hα + µh + τihL2 + τexv Ic2 n λ τexh + c n Hµh + τexhL Hµv + τexvL Hα + µh + τihL + µh Hµh + τexhL Hµv + τexvL Hα + µh + τihLM Iµh2 + µv τexh + α Hµh + µv + τexhL + Hµv + τexhL τih + µh H2 µv + τexh + τihLMMM ë HHµh + τexhL Hµv + τexvL Hα + µh + τihL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + IIµh µv Hµh + τexhL Hµv + τexvL2 Hα + µh + τihLM ë Hc n τexv + µh Hµv + τexvLL + Hc λ µh µv τexh τexv Hα + 2 µh + τexh + τihLL ê Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLM H− 1 + ξL;
a5 = µh µv Hµh + τexhL Hµv + τexvL Hα + µh + τihL H− 1 + ξL;
H∗ Masukan Nilai−nilai Parameter ∗L λ = 2.244 ∗ 10−5 ; α = 0.003; µh = 1 ê 28 000; µv = 1 ê 20; τexh = 1 ê 10; τih = 1 ê 4; τexv = 1 ê 9; bi = 0.3; pvh = 0.4; c = bi pvh; nh = 1; nv = 1; n = nh ê nv; H∗ Definisikan ξ ∗L Clear@ξD ξ= Ic2 n λ τexh τexvM ë Hµh µv Hµv + τexvL Hµh + τexhL Hµh + α + τihLL; H∗ Hitung a2 yang Belum Disederhanakan ∗L a2 H∗ Hitung a2 yang Sudah Disederhanakan ∗L Clear@a2D a2 = α µh + µh2 + α µv + 2 µh µv + α τexh + µh τexh + µv τexh + µh τih + µv τih + τexh τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLL Hα + 2 µh + µv + τexh + τihLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + τexv Hα + 2 µh + µv + τexh + τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLL H1 + Hc n Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê HHµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLLL + Hµh µv τexv Hµv + τexvLL ê Hc n τexv + µh Hµv + τexvLL H− 1 + ξL
0.10791 0.10791 H∗ Hitung a3 yang Belum Disederhanakan ∗L a3 H∗ Hitung a3 yang Sudah Disederhanakan ∗L Clear@a3D a3 = α µh µv + µh2 µv + α µv τexh + µh µv τexh + µh µv τih + µv τexh τih + Ic λ τexh
Hc n τexv + µh Hµv + τexvLL Iµh2 + µv τexh + α Hµh + µv + τexhL + Hµv + τexhL τih + µh H2 µv + τexh + τihLMM ë HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + τexv Iα µh + µh2 + α µv + 2 µh µv + α τexh + µh τexh + µv τexh + µh τih + µv τih + τexh τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLL τihL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLM H1 + Hc n Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê HHµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLLL + Hµh µv τexv Hµv + τexvL HHα + 2 µh + τexh + τihL ê Hc n τexv + µh Hµv + τexvLL + Hc λ τexhL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLLL H− 1 + ξL 0.00818896 0.00818896
H∗ Hitung a4 yang Belum Disederhanakan ∗L a4 H∗ Hitung a4 yang Sudah Disederhanakan ∗L Clear@a4D a4 = Ic λ τexh
Iµv Hµh + τexhL2 Hc n τexv + µh Hµv + τexvLL Hα + µh + τihL2 + τexv Ic2 n λ τexh + c n Hµh + τexhL Hµv + τexvL Hα + µh + τihL + µh Hµh + τexhL Hµv + τexvL Hα + µh + τihLM Iµh2 + µv τexh + α Hµh + µv + τexhL + Hµv + τexhL τih + µh H2 µv + τexh + τihLMMM ë HHµh + τexhL Hµv + τexvL Hα + µh + τihL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + IIµh µv Hµh + τexhL Hµv + τexvL2 Hα + µh + τihLM ë Hc n τexv + µh Hµv + τexvLL + Hc λ µh µv τexh τexv Hα + 2 µh + τexh + τihLL ê Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLM H− 1 + ξL
9.9597 × 10−8 9.9597 × 10−8 H∗ Hitung a5 yang Belum Disederhanakan ∗L a5 H∗ Hitung a5 yang Sudah Disederhanakan ∗L Clear@a5D a5 = µh µv Hµh + τexhL Hµv + τexvL Hα + µh + τihL H− 1 + ξL
− 3.692 × 10−9 − 3.692 × 10−9 Clear@ξ, a2, a3, a4, a5D H∗ Ambil sebarang nilai ξ > 1 ∗L ξ = 1.5;
H∗ Definisikan Koefisien−koefisien a2,a3, a4 dan a5 yang Sudah Disederhanakan ∗L
a2 = α µh + µh2 + α µv + 2 µh µv + α τexh + µh τexh + µv τexh + µh τih + µv τih + τexh τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLL Hα + 2 µh + µv + τexh + τihLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + τexv Hα + 2 µh + µv + τexh + τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLLL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLL H1 + Hc n Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê HHµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLLL + Hµh µv τexv Hµv + τexvLL ê Hc n τexv + µh Hµv + τexvLL H− 1 + ξL;
a3 = α µh µv + µh2 µv + α µv τexh + µh µv τexh + µh µv τih + µv τexh τih + Ic λ τexh Hc n τexv + µh Hµv + τexvLL Iµh2 + µv τexh + α Hµh + µv + τexhL + Hµv + τexhL τih + µh H2 µv + τexh + τihLMM ë HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + τexv Iα µh + µh2 + α µv + 2 µh µv + α τexh + µh τexh + µv τexh + µh τih + µv τih + τexh τih + Hc λ τexh Hc n τexv + µh Hµv + τexvLL τihL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLM H1 + Hc n Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL ê HHµh + τexhL Hc n τexv + µh Hµv + τexvLL Hα + µh + τihLLL + Hµh µv τexv Hµv + τexvL HHα + 2 µh + τexh + τihL ê Hc n τexv + µh Hµv + τexvLL + Hc λ τexhL ê HHµv + τexvL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLLLL H− 1 + ξL; a4 = Ic λ τexh Iµv Hµh + τexhL2 Hc n τexv + µh Hµv + τexvLL Hα + µh + τihL2 + τexv Ic2 n λ τexh + c n Hµh + τexhL Hµv + τexvL Hα + µh + τihL + µh Hµh + τexhL Hµv + τexvL Hα + µh + τihLM Iµh2 + µv τexh + α Hµh + µv + τexhL + Hµv + τexhL τih + µh H2 µv + τexh + τihLMMM ë HHµh + τexhL Hµv + τexvL Hα + µh + τihL Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLL + IIµh µv Hµh + τexhL Hµv + τexvL2 Hα + µh + τihLM ë Hc n τexv + µh Hµv + τexvLL + Hc λ µh µv τexh τexv Hα + 2 µh + τexh + τihLL ê Hc λ τexh + µv Hµh + τexhL Hα + µh + τihLLM H− 1 + ξL; a5 = µh µv Hµh + τexhL Hµv + τexvL Hα + µh + τihL H− 1 + ξL; H∗ Pembuktian Kriteria Routh−Hurwitz ∗L a1 > 0 a2 > 0 a3 > 0 a4 > 0 a5 > 0 a1 a2 a3 > a32 + a12 a4 Ha1 a4 − a5L Ia1 a2 a3 − a32 − a12 a4M > a5 Ha1 a2 − a3L2 + a1 a52
True True True True True True True
Lampiran 5.
Program Simulasi untuk Kondisi ¬0 < 1
H∗ Scipt Gambar 3 ∗L Clear@µv, c, bi, T, sh0, eh0, ih0, ev0, iv0D; Needs@"PlotLegends`"D; ManipulateAModuleA8sol1, sol2, sh, eh, ih, rh, sv, ev, iv, t<,
rh@t_D := nh − Hsh@tD + eh@tD + ih@tDL; sv@t_D := nv − Hev@tD + iv@tDL; sol1 = With@8c = bi ∗ pvh<, NDSolve@8sh '@tD λ − Hn ∗ c ∗ iv@tD + µhL ∗ sh@tD, eh '@tD n ∗ c ∗ iv@tD ∗ sh@tD − Hτexh + µhL ∗ eh@tD, ih '@tD τexh ∗ eh@tD − Hτih + α + µhL ∗ ih@tD, ev '@tD c ∗ ih@tD ∗ H1 − ev@tD − iv@tDL − Hτexv + µvL ∗ ev@tD, iv '@tD τexv ∗ ev@tD − µv ∗ iv@tD, sh@0D sh0, eh@0D eh0, ih@0D ih0, ev@0D ev0, iv@0D iv0<, 8sh, eh, ih, ev, iv<, 8t, 0, T
H∗ Scipt Gambar 4 dan 5 ∗L Clear@sol3, plot1, plot2, x, y, z, µv, c, bi, R0, sh, eh, ih, rh, sv, ev, ivD; rh@t_D := nh − Hsh@tD + eh@tD + ih@tDL; sv@t_D := nv − Hev@tD + iv@tDL; sol3@x_, µv_, bi_D := First@x ê. With@8c = bi ∗ pvh<, NDSolve@8sh '@tD λ − Hn ∗ c ∗ iv@tD + µhL ∗ sh@tD, eh '@tD n ∗ c ∗ iv@tD ∗ sh@tD − Hτexh + µhL ∗ eh@tD, ih '@tD τexh ∗ eh@tD − Hτih + α + µhL ∗ ih@tD, ev '@tD c ∗ ih@tD ∗ H1 − ev@tD − iv@tDL − Hτexv + µvL ∗ ev@tD, iv '@tD τexv ∗ ev@tD − µv ∗ iv@tD, sh@0D sh0, eh@0D eh0, ih@0D ih0, ev@0D ev0, iv@0D iv0<, 8sh, eh, ih, ev, iv<, 8t, 0, T
plot1@x_, y_, z_, xx_StringD := Plot@ Evaluate@Table@sol3@x, i, 0.30D, 8i, 0.03, 0.09, 0.02
Lampiran 6.
Program Simulasi untuk Kondisi ¬0 > 1
H∗ Scipt Gambar 6 ∗L Clear@µv, c, biD ManipulateA
ModuleA8sol4, sol5, sh, eh, ih, rh, sv, ev, iv, t<, rh@t_D := nh − Hsh@tD + eh@tD + ih@tDL; sv@t_D := nv − Hev@tD + iv@tDL; sol4 = With@8c = bi ∗ pvh<, NDSolve@8sh '@tD λ − Hn ∗ c ∗ iv@tD + µhL ∗ sh@tD, eh '@tD n ∗ c ∗ iv@tD ∗ sh@tD − Hτexh + µhL ∗ eh@tD, ih '@tD τexh ∗ eh@tD − Hτih + α + µhL ∗ ih@tD, ev '@tD c ∗ ih@tD ∗ H1 − ev@tD − iv@tDL − Hτexv + µvL ∗ ev@tD, iv '@tD τexv ∗ ev@tD − µv ∗ iv@tD, sh@0D sh0, eh@0D eh0, ih@0D ih0, ev@0D ev0, iv@0D iv0<, 8sh, eh, ih, ev, iv<, 8t, 0, T
H∗ Scipt Gambar 7 dan 8 ∗L Clear@sol6, plot3, plot4, x, y, z, µv, c, bi, R0, sh, eh, ih, rh, sv, ev, ivD rh@t_D := nh − Hsh@tD + eh@tD + ih@tDL; sv@t_D := nv − Hev@tD + iv@tDL; sol6@x_, µv_, bi_D := First@x ê. With@8c = bi ∗ pvh<, NDSolve@8sh '@tD λ − Hn ∗ c ∗ iv@tD + µhL ∗ sh@tD, eh '@tD n ∗ c ∗ iv@tD ∗ sh@tD − Hτexh + µhL ∗ eh@tD, ih '@tD τexh ∗ eh@tD − Hτih + α + µhL ∗ ih@tD, ev '@tD c ∗ ih@tD ∗ H1 − ev@tD − iv@tDL − Hτexv + µvL ∗ ev@tD, iv '@tD τexv ∗ ev@tD − µv ∗ iv@tD, sh@0D sh0, eh@0D eh0, ih@0D ih0, ev@0D ev0, iv@0D iv0<, 8sh, eh, ih, ev, iv<, 8t, 0, T
plot3@x_, y_, z_, xx_StringD := Plot@ Evaluate@Table@sol6@x, i, 0.40D, 8i, 0.01, 0.07, 0.02