MODEL SIRS PADA PROSES PENULARAN PENYAKIT INFLUENZA DENGAN POPULASI YANG TERINFEKSI VIRUS LIA MULYANAH

MODEL SIRS PADA PROSES PENULARAN PENYAKIT INFLUENZA DENGAN POPULASI YANG TERINFEKSI VIRUS

LIA MULYANAH

DEPARTEMEN MATEMATIKA FAKULTAS MATEMATIKA DAN ILMU PENGETAHUAN ALAM INSTITUT PERTANIAN BOGOR BOGOR 2008


LIA MULYANAH G54104018


ABSTRACT LIA MULYANAH. Application of SIRS Model in The Process of Outbreak of Influenza Disease in Infected Population. Supervisors : ALI KUSNANTO and ENDAR H. NUGRAHANI. Influenza is a viral infectious disease in which its infection was studied using the SIRS model. In this paper the analysis of the influenza infection in a virus infected population with different initial quantity of population already infected. Stability analysis of infection virus influence used SIRS Model produced 2 equilibrium points, in example one is called free of disease (T1) and the other is called an endemic case (T2). The parameter value is related to the basic reproduction rate R0 . If the basic reproduction rate was less than 1 ( R0 < 1 ), the equilibrium point T1 was a stable node and the equilibrium point T2 was an unstable spiral so virus was loss from population. If the basic reproduction rate was greater than 1 ( R 0 > 1 ), the equilibrium point T1 was a saddle and the equilibrium point T2 was a stable spiral so virus stayed in population. This work also compared the dynamic analysis of population with different parameter values, that is when initial population of infection I (0) is 1% and 10% from population. From the analysis, we have known that the 10% infected population was faster to stability than the 1% infected population.

ABSTRAK LIA MULYANAH. Model SIRS Pada Proses Penularan Penyakit Influenza dengan Populasi yang Terinfeksi Virus. Dibimbing oleh ALI KUSNANTO dan ENDAR H. NUGRAHANI. Penyakit influenza termasuk salah satu jenis penyakit menular yang disebabkan oleh virus influenza. Proses penularan penyakit ini dapat dimodelkan dengan menggunakan model SIRS. Dalam tulisan ini dianalisis penularan penyakit influenza pada populasi yang terinfeksi virus dengan jumlah individu yang terinfeksi awal berbeda. Analisis kestabilan pada proses penularan virus influenza dengan menggunakan model SIRS ini menghasilkan 2 titik tetap, yaitu titik tetap tanpa penyakit (T1) dan titik tetap endemik (T2) dengan beberapa nilai parameter. Nilai parameter ini ada kaitannya dengan konsep reproduksi dasar virus R0 . Ketika reproduksi dasar virus kurang dari 1 ( R0 < 1 ) titik tetap tanpa penyakit (T1) simpul stabil, dan titik tetap endemik (T2) spiral tak stabil sehingga virus akan hilang dari populasi. Ketika reproduksi dasar virus lebih besar dari 1 ( R0 > 1 ) titik tetap tanpa penyakit (T1) akan sadel, dan titik tetap endemik (T2) akan spiral stabil sehingga virus akan bertahan dalam populasi. Selain menganalisis kestabilan, dibandingkan juga analisis untuk dinamika populasi dengan nilai parameter yang berbeda, yaitu ketika jumlah populasi awal yang terinfeksi virus I (0) 1% dari populasi dan 10% dari populasi. Dari hasil analisis diperoleh bahwa ketika populasi awal yang terinfeksi 10% dari populasi, maka populasi akan lebih cepat menuju kestabilan daripada ketika populasi awal yang terinfeksi 1% dari populasi.


Skripsi Sebagai salah satu syarat untuk memperoleh gelar Sarjana Sains pada Fakultas Matematika dan Ilmu Pengetahuan Alam Institut Pertanian Bogor

Oleh: LIA MULYANAH G54104018


Judul: Model SIRS pada Proses Penularan Penyakit Influenza dengan Populasi yang Terinfeksi Virus Nama: Lia Mulyanah NIM : G54104018

Menyetujui

Pembimbing I,

Pembimbing II,

Drs. Ali Kusnanto, M.Si. NIP : 131913135

DR. Ir. Endar H. Nugrahani, MS NIP : 131842411

Mengetahui Dekan Fakultas Matematika dan Ilmu Pengetahuan Alam Institut Pertanian Bogor

Dr. Drh. Hasim, DEA NIP : 131578806

Tanggal Lulus :

KATA PENGANTAR Puji dan syukur penulis panjatkan kepada Allah SWT, yang telah melimpahkan rahmat dan karunia-Nya sehingga penulis dapat menyelesaikan karya ilmiah ini. Penyusunan karya ilmiah ini tidak terlepas dari bantuan berbagai pihak. Oleh karena itu, dalam kesempatan ini penulis mengucapkan terimakasih yang sebesar-besarnya kepada semua pihak yang telah membantu. Terima kasih kepada Bapak Drs. Ali Kusnanto, M.Si. selaku dosen pembimbing 1 yang telah dengan sabar membimbing penulis sehingga penulis bisa menyelesaikan skripsi ini. Begitu pun untuk dosen pembimbing 2 Ibu DR. Ir. Endar H. Nugrahani, MS. Terima kasih kepada Bapak Ir. N. K. Kutha Ardana, M. Sc. yang telah bersedia menjadi dosen penguji. Terima kasih penulis haturkan kepada semua dosen Departemen Matematika IPB yang telah memberikan ilmu dan motivasi kepada penulis. Terima kasih kepada semua staf Departemen Matematika IPB yang telah banyak membantu penulis, Bu Ade, Pak Yono, Pak Bono, Bu Marisi, Mas Deni, Bu Susi staf lainnya yang tidak bisa disebutkan semuanya. Terima kasih untuk keluarga tercinta, Ibu yang tiada henti slalu mendukung dan tetap berdoa dengan penuh kasih sayang. Bapak yang slalu memberi smangat, doa, kepercayaan dan kasih sayangnya selama ini. Kedua adik ku, Eva dan Ai yang slalu membangkitkan smangat untuk segera menyelesaikan skripsi ini. Kakek, Nenek, ua, bibi, paman dan semua sepupu yang slalu menjadi sumber smangat. Terima kasih untuk Supriatna dan keluarga yang slama ini slalu membantu. Teman-teman yang telah menjadi sahabat dan memberikan kenangan yang begitu indah dalam kebersamaan kita, Eli, Mukti, Rina, Janah, Rizul, Eva, Nidia, Ika, Sifa, Intan, Dian, Dee-dee, Liay, Mahar, Eni, Armi, Ayu, Ani, Tities, Tia, Fitri, Darwisah, Endit, Sita, Niken, Maryam, Peny, Uwi, Eci, Eneng, Reni, Roma, Ro’fah, Rangga M, Fariz, Adji, Fred, Mimin, Mahnur, Triyadi, Idris, Yaya, Mora, Iboy, Majid, Jali, Great dan semuanya yang tidak dapat disebutkan satu per satu. Kakak-kakak kelas Math 40: Teh Ulfa, Kak Abay, Kak Rusli, Teh Nci, Kak Aam yang telah membantu dan memberi smangat. Teman-teman sekostan yang slalu menyemangati, Mel, Wheny, Widia, Jo, Ana, Efi, Sinta, Ipo, Ima, Riska, Kak Marna, Popon, Eri, Uni, Ika dan teh Reni. Dan semua pihak yang telah membantu yang tidak dapat penulis sebutkan satu per satu. Semoga karya ilmiah ini bermanfaat bagi semua pihak, khususnya bagi pihak yang memiliki ketertarikan dalam bidang Matematika dan dapat menjadi inspirasi untuk penelitian-penelitian selanjutnya.

Bogor, Mei 2008

Lia Mulyanah

iv

RIWAYAT HIDUP Penulis dilahirkan pada tanggal 3 Oktober 1985 di Ciamis. Penulis merupakan anak pertama dari tiga bersaudara, anak dari pasangan Maman Kusmanto dan Ebah Nurjanah. Pada tahun 2001, penulis melanjutkan pendidikan ke jenjang pendidikan sekolah menengah umum di SMU Negeri 1 Ciamis. Pada tahun 2004, penulis lulus dari tingkat SMU dan melanjutkan pendidikan ke tingkat Perguruan Tinggi di Institut Pertanian Bogor (IPB). Penulis masuk IPB melalui jalur USMI dan mengambill program studi Matematika, Departemen Matematika, Fakultas matematika dan Ilmu Pengetahuan Alam, Institut Pertanian Bogor untuk tingkat Strata 1 (S1). Selama menjalani pendidikan di IPB, penulis memperoleh kesempatan untuk mendapatkan beasiswa PPA pada tahun 2005 dan 2006. Selain itu, penulis juga terlibat dalam beberapa kegiatan, antara lain mengikuti Program Guru Tambahan tingkat Sekolah Dasar (2004/2005), anggota Lingkar Muslim Matematika Himpro GUMATIKA (2005/2006, 2006/2007), Staf pengajar privat Kalkulus 1 Himpro GUMATIKA (2005/2006,2006/2007). Selain itu penulis juga aktif dalam berbagai kepanitiaan, diantaranya Penyambutan Mahasiswa Baru tahun 2005, Ramadhan In Action tahun 2005 dan 2006, Pesta Sains Nasional FMIPA IPB tahun 2006, Musyawarah Wilayah Ikahimatika tahun 2006, Welcome Ceremony Matematika tahun 2006 dan 2007.

v

DAFTAR ISI PENDAHULUAN ........................................................................................................ 1.1. Latar Belakang............................................................................................. 1.2. Tujuan..........................................................................................................

1 1 1

LANDASAN TEORI .................................................................................................... Sistem Persamaan Diferensial ............................................................................ Titik Tetap .......................................................................................................... Nilai Eigen dan Vektor Eigen............................................................................. Analisis Kestabilan Titik Tetap .......................................................................... Pelinearan ........................................................................................................... Model Primitif Penularan Virus .........................................................................

2 2 2 2 3 3 3

PEMBAHASAN ........................................................................................................... 3.1. Perumusan Model........................................................................................ 3.2. Titik Tetap ................................................................................................... 3.2.1. Titik Tetap Tanpa Penyakit ...................................................................... 3.2.2. Titik Tetap Endemik................................................................................. 3.3. Analisis Kestabilan Titik Tetap ketika Populasi Awal yang Terinfeksi 1 Orang ........................................................................................................ 3.3.1. Analisis Kestabilan Titik Tetap untuk R0 < 1........................................... 3.3.2. Analisis kestabilan titik tetap untuk R0 > 1 .............................................. 3.4. Dinamika Populasi Penularan Virus Influenza............................................ 3.4.1. Dinamika populasi untuk R0 < 1 .............................................................. 3.4.2. Dinamika populasi untuk R0 > 1 .............................................................. Perbandingan Dinamika Populasi ketika Populasi Awal yang Terinfeksi 1% dan 10% dari Populasi untuk R0 < 1................................................................... Perbandingan Dinamika Populasi ketika Populasi Awal yang Terinfeksi 1% dan 10% dari Populasi untuk R0 > 1...................................................................

5 5 6 6 7

KESIMPULAN .............................................................................................................

25

DAFTAR PUSTAKA ...................................................................................................

26

LAMPIRAN..................................................................................................................

27

vi

9 9 12 15 15 18 22 24

DAFTAR LAMPIRAN Penurunan persamaan (4), (5), (6) dari persamaan (1), (2), (3).................................... Pencarian titik tetap dari persamaan (4), (5) dan (6) ..................................................... Pencarian titik tetap dari persamaan (4), (5), (6) dengan menggunakan Mathematica . Penentuan nilai eigen titik tetap tanpa penyakit ............................................................ Penentuan nilai eigen titik tetap tanpa penyakit dengan Mathematica 6.0.................... Penentuan nilai eigen titik tetap endemic...................................................................... Bukti ............................................................................................................................ Gambar bidang fase ketika β = 0.0002 dan v = 0.9 dengan menggunakan Mathematica 5.2........................................................................................................... Gambar bidang fase ketika β = 0.0002 dan v = 0.8 dengan menggunakan Mathematica 5.2............................................................................................................ Gambar bidang fase ketika β = 0.0001 dan v = 0.9 dengan menggunakan Mathematica 5.2............................................................................................................ Gambar bidang fase ketika β = 0.0002 dan v = 0.5 dengan menggunakan Mathematica 5.2............................................................................................................ Gambar bidang fase ketika β = 0.0002 dan v = 0.4 dengan menggunakan Mathematica 5.2............................................................................................................ Gambar bidang fase ketika β = 0.0003 dan v = 0.4 dengan menggunakan Mathematica 5.2............................................................................................................ Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.9 ........... Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.8 ........... Program Mathematica untuk dinamika populasi ketika β = 0.0001 dan v = 0.9 ........... Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.5 ........... Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.4 ........... Program Mathematica untuk dinamika populasi ketika β = 0.0003 dan v = 0.4 ........... Program Mathematica ketika R0 < 1 dan populasi yang terinfeksi 1% dari total populasi ......................................................................................................................... Program Mathematica ketika R0 < 1 dan populasi yang terinfeksi 10% dari total populasi ......................................................................................................................... Program Mathematica ketika R0 > 1 dan populasi yang terinfeksi 1% dari total populasi ......................................................................................................................... Program Mathematica ketika R0 > 1 dan populasi yang terinfeksi 10% dari total populasi .........................................................................................................................

vii

28 29 31 31 32 32 37 38 39 40 42 43 44 46 46 47 47 48 49 49 50 50 51

DAFTAR GAMBAR Gambar 1 Model SIR .................................................................................................... Gambar 2 Model SIRS pada populasi tertutup.............................................................. Gambar 3 Model SIRS pada populasi terbuka .............................................................. Gambar 4 Diagram Model SIRS ................................................................................... Gambar 5 Orbit kestabilan pada bidang A, B, dan C ketika β = 0.0002, dan v = 0.9. ........... Gambar 6 Proyeksi bidang fase pada bidang A, B ketika β = 0.0002, dan v = 0.9. .................. Gambar 7 Proyeksi bidang fase pada bidang A, C ketika β = 0.0002, dan v = 0.9. ................ Gambar 8 Proyeksi bidang fase pada bidang B, C ketika β = 0.0002, dan v = 0.9. ................... Gambar 9 Orbit kestabilan pada bidang A, B, dan C ketika β = 0.0002, dan v = 0.8. ........... Gambar 10 Proyeksi bidang fase pada bidang A, B ketika β = 0.0002, dan v = 0.8. .............. Gambar 11 Proyeksi bidang fase pada bidang A, C ketika β = 0.0002, dan v = 0.8. ............. Gambar 12 Proyeksi bidang fase pada bidang B, C ketika β = 0.0002, dan v = 0.8. ............ Gambar 13 Orbit kestabilan pada bidang A, B, dan C ketika β = 0.0001, dan v = 0.9. ......... Gambar 14 Proyeksi bidang fase pada bidang A, B ketika β = 0.0001, dan v = 0.9. ............. Gambar 15 Proyeksi bidang fase pada bidang A, C ketika β = 0.0001, dan v = 0.9. ............ Gambar 16 Proyeksi bidang fase pada bidang B, C ketika β = 0.0001, dan v = 0.9. ............ Gambar 17 Bidang fase tiga dimensi (A, B, C) ketika β = 0.0002, dan v = 0.5 ................. Gambar 18 Proyeksi bidang fase pada bidang A, B ketika β = 0.0002, dan v = 0.5 ........... Gambar 19 Proyeksi bidang fase pada bidang A, C ketika β = 0.0002, dan v = 0.5 ................. Gambar 20 Proyeksi bidang fase pada bidang B, C ketika β = 0.0002, dan v = 0.5. ............... Gambar 21 Bidang fase tiga dimensi (A, B, C) ketika β = 0.0002, dan v = 0.4 ...................... Gambar 22 Proyeksi bidang fase pada bidang A, B ketika β = 0.0002, dan v = 0.4 ............ Gambar 23 Proyeksi bidang fase pada bidang A, C ketika β = 0.0002, dan v = 0.4 ............... Gambar 24 Proyeksi bidang fase pada bidang B, C ketika β = 0.0002, dan v = 0.4 ...................... Gambar 25 Bidang fase 3 dimensi (A, B, C) ketika β = 0.0003, dan v = 0.4 ...................... Gambar 26 Proyeksi bidang fase pada bidang A, B ketika β = 0.0003, dan v = 0.4 ................ Gambar 27 Proyeksi bidang fase pada bidang A, C ketika β = 0.0003, dan v = 0.4 ................ Gambar 28 Proyeksi bidang fase pada bidang B, C ketika β = 0.0003, dan v = 0.4 ................ Gambar 29 Dinamika populasi A terhadap waktu t............................................................. Gambar 30 Dinamika populasi B terhadap waktu t................................................................... Gambar 31 Dinamika populasi C terhadap waktu t ............................................................. Gambar 32 Dinamika populasi A, B, C terhadap waktu t.......................................................... Gambar 33 Dinamika populasi A terhadap waktu t ............................................................. Gambar 34 Dinamika populasi B terhadap waktu t ............................................................. Gambar 35 Dinamika populasi C terhadap waktu t ............................................................. Gambar 36 Dinamika populasi A, B, C terhadap waktu t ..................................................... Gambar 37 Dinamika populasi A terhadap waktu t ............................................................. Gambar 38 Dinamika populasi B terhadap waktu t ............................................................. Gambar 39 Dinamika populasi C terhadap waktu t ............................................................. Gambar 40 Dinamika populasi A, B, C terhadap waktu t ..................................................... Gambar 41 Dinamika populasi A terhadap waktu t ............................................................. Gambar 42 Dinamika populasi B terhadap waktu t ............................................................. Gambar 43 Dinamika populasi C terhadap waktu t ............................................................. Gambar 44 Dinamika populasi A, B, C terhadap waktu t ..................................................... Gambar 45 Dinamika populasi A terhadap waktu t ............................................................. Gambar 46 Dinamika populasi B terhadap waktu t ............................................................. Gambar 47 Dinamika populasi C terhadap waktu t ............................................................. Gambar 48 Dinamika populasi A, B, C terhadap waktu t ..................................................... Gambar 49 Dinamika populasi A terhadap waktu t ............................................................. Gambar 50 Dinamika populasi B terhadap waktu t .............................................................

viii

4 4 4 5 9 9 10 10 10 10 10 11 11 11 12 12 12 12 13

13 13 13 14 14 14 14 15 15 16 16 16 16

17 17 17 17 17 18 18 18 19 19 19 19 19 20 20 20 20 21

Gambar 51 Dinamika populasi C terhadap waktu t ............................................................. Gambar 52 Dinamika populasi A, B, dan C terhadap waktu t ............................................... Gambar 53 Dinamika populasi A terhadap waktu t ............................................................. Gambar 54 Dinamika populasi A terhadap waktu t ............................................................. Gambar 55 Dinamika populasi B terhadap waktu t ............................................................. Gambar 56 Dinamika populasi B terhadap waktu t ............................................................. Gambar 57 Dinamika populasi C terhadap waktu t ............................................................. Gambar 58 Dinamika populasi C terhadap waktu t ............................................................. Gambar 59 Dinamika populasi A, B dan C terhadap waktu t ................................................ Gambar 60 Dinamika populasi A, B dan C terhadap waktu t ................................................ Gambar 61 Dinamika populasi A terhadap waktu t ............................................................. Gambar 62 Dinamika populasi A terhadap waktu t ............................................................. Gambar 63 Dinamika populasi B terhadap waktu t ............................................................. Gambar 64 Dinamika populasi B terhadap waktu t ............................................................. Gambar 65 Dinamika populasi C terhadap waktu t ............................................................. Gambar 66 Dinamika populasi C terhadap waktu t ............................................................. Gambar 67 Dinamika populasi A, B dan C terhadap waktu t ................................................ Gambar 68 Dinamika populasi A, B dan C terhadap waktu t ................................................

ix

21 21 22 22 22 22 23 23 23 23 24 24 24 24 24 24 25 25

PENDAHULUAN 1.1. Latar Belakang Berbagai jenis penyakit saat ini semakin banyak. Salah satu penyebabnya, gaya hidup dan lingkungan yang semakin tidak sehat. Secara umum ada dua jenis penyakit, yaitu penyakit menular dan tidak menular. Dalam kelompok penyakit menular ada yang ringan dan ada yang berat. Yang ringan misalnya influenza dan diare. Sedangkan yang berat seperti HIV/AIDS, polio, demam berdarah, campak, TBC, malaria, flu burung, SARS, dan sederet penyakit lainnya. Menular atau tidaknya suatu penyakit tetap harus diwaspadai dan tidak boleh dianggap remeh. Sebab, ketika seseorang terkena suatu penyakit aktivitas kehidupannya akan terganggu. Apalagi jika penyakitnya sudah parah, bisa mengakibatkan kematian. Pada penyakit menular, medium penularannya bermacam-macam. Ada yang karena kontak langsung dengan penderita, lewat udara, kotoran, atau lewat perantara binatang. Virus influenza telah ada berabad-abad yang lalu dan telah menyebabkan kematian dan berbagai macam penyakit. Setiap tahun virus influenza menyerang 25-50 juta manusia, dan diperkirakan 20-40 ribu manusia di Amerika meninggal disebabkan oleh virus influenza. Ada tiga tipe serologi virus influenza, yaitu tipe A, B, dan C. Tipe A dapat menginfeksi manusia dan hewan, sedangkan tipe B dan C hanya menginfeksi manusia. Influenza adalah penyakit yang disebabkan oleh virus influenza pada saluran pernafasan dari hidung sampai trachea. Influenza (flu) yang merupakan suatu infeksi virus yang menyebabkan demam, hidung meler, sakit kepala, batuk, tidak enak badan (malaise) dan peradangan pada selaput lendir hidung dan saluran pernafasan. Influenza merupakan penyakit serius, tetapi sebagian besar penderita akan kembali sehat dalam waktu 710 hari. Selain menyerang manusia, ternyata virus influenza juga dapat ditemukan pada beberapa binatang, seperti unggas, babi, bebek, ikan paus, kuda, dan anjing laut. Jangka waktu seseorang terkena virus hingga munculnya gejala adalah satu sampai empat hari, dengan rata-rata dua hari. Sedangkan

periode seseorang dapat menularkan penyakitnya ke orang lain bervariasi untuk tiap usia. Penularan sudah mulai terjadi dari sebelum penderita merasa sakit, yang berlanjut hingga tiga sampai tujuh hari setelah timbul gejala pertama pada orang dewasa. Sedangkan pada anak-anak dapat lebih dari satu minggu. Untuk menganalisis dinamika penyebaran virus banyak metode yang dapat digunakan, diantaranya cellular automata yang pernah digunakan oleh Keeling dan Gillingan, forestfire lattices yang digunakan oleh Rhodes dkk., sedangkan dalam tulisan ini akan dibahas pemodelan matematika untuk suatu penyakit menular dengan memisahkan populasi menjadi beberapa kelas, yang lebih dikenal dengan pendekatan SIRS. Dalam tulisan ini akan dibahas model penyebaran penyakit menular yaitu influenza. Analisis kestabilan dan dinamika populasinya akan menjadi pembahasan. Pertama, ditentukan titik tetapnya yang terdiri dari titik tetap tanpa penyakit dan titik tetap endemik. Selanjutnya ditentukan matriks Jacobi dengan melakukan pelinearan setiap persamaan yang ada terhadap setiap variabel. Kemudian ditentukan nilai eigen dengan menyelesaikan persamaan karakteristik, nilai eigen tersebut akan digunakan untuk menganalisis kestabilan titik tetapnya. Selanjutnya laju reproduksi dasar virus ditentukan untuk mengetahui laju penularan penyakit tersebut dengan cara mengubah beberapa nilai parameter. Hal tersebut juga digunakan untuk menentukan bidang fasenya agar diketahui jenis kestabilannya. Dianalisis juga grafik dinamika populasi dengan beberapa nilai laju reproduksi dasar virus yang berbeda dengan mengubah beberapa parameter. Terakhir adalah membandingkan pengaruh dari jumlah populasi awal yang terinfeksi berbeda, yaitu 1% dan 10% dari total populasi. 1.2. Tujuan Tujuan dari penulisan karya ilmiah ini adalah menganalisa penyebaran penyakit influenza menggunakan model SIRS dan menentukan kestabilan serta dinamika populasi untuk beberapa kondisi yang berbeda.

2

LANDASAN TEORI Definisi 1 [Sistem Persamaan Diferensial Linear (SPDL)]

Definisi 5 [Titik Tetap Stabil Asimtotik Lokal]

Persamaan diferensial orde n dikatakan linear jika dapat ditulis dalam bentuk:

Titik x* dikatakan titik tetap stabil asimtotik lokal jika titik x* stabil dan terdapat

a0 ( t )

d nx dt n

+ a (t ) 0

d n −1 x dt n −1

+ ... + a

n −1

(t )

dx + dt

(1.1)

a ( t ) x = f (t )

ε > 0 sedemikian sehingga jika x* − x0 < ε maka lim x(t ) = x* . t →∞

(Szidarovszky & Bahill, 1998)

n

dengan a (t ) dan f (t ) adalah fungsi dari waktu t dan a (t ) ≠ 0 . (Farlow, 1994) Definisi 2 [Sistem Persamaan Diferensial Mandiri]

Misalkan diberikan suatu sistem persamaan diferensial orde 1 sebagai berikut: dx = x& = f ( x1 , x2 ,...) (1.2) dt f fungsi kontinu bernilai real, dengan laju perubahan dinyatakan dengan fungsi itu sendiri serta tidak berubah terhadap waktu, maka sistem (1.2) merupakan sistem persamaan differensial mandiri. (Verhulst, 1990) Definisi 3 [Titik tetap]

Misalkan diberikan sistem persamaan diferensial (SPD) sebagai berikut dx n = x& = f ( x1 , x2 ,...), ( x1 , x2 ,...) ∈ R (1.3) dt Suatu titik x* yang memenuhi f ( x* ) = 0 disebut titik keseimbangan atau titik tetap dari sistem (1.3). (Verhulst, 1990) Definisi 4 [Titik tetap stabil]

Titik x* adalah titik tetap sebuah SPD dan x (t ) adalah solusi yang memenuhi kondisi awal x(0) = x0 dan x0 ≠ x . Titik x dikatakan titik tetap stabil jika terdapat ε 0 > 0 , yang memenuhi sifat berikut: untuk setiap ε1 , 0 < ε1 < ε 0 , terdapat ε > 0 *

*

sedemikian sehingga jika x* − x0 < ε maka

x* − x(t ) < ε1 , untuk setiap t > t0 . (Szidarovszky & Bahill, 1998)

Definisi 6 [Titik Tetap Stabil Asimtotik Global]

Titik x* dikatakan titik tetap stabil asimtotik global jika titik x* stabil dan x0 ∈ Ω ⊆ R n , lim x(t ) = x* . t →∞

(Szidarovszky & Bahill, 1998) Definisi 4 menyatakan bahwa titik x* stabil jika seluruh orbit (lintasan kurva dan yang menggambarkan solusi x (t ) ) berada pada radius ε1 , jika nilai awal ( x0 ) yang dipilih cukup dekat dengan x* . Pada Definisi 5 (titik tetap stabil asimtotik lokal), dipilih nilai awal ( x0 ) yang cukup dekat dengan x* sehingga solusi x (t ) adalah x* untuk t → ∞ . Sedangkan pada definisi 6 (titik tetap stabil asimtotik global), dipilih nilai awal ( x0 ) di luar radius ε 0 sehingga solusi x (t ) adalah x* untuk t → ∞ . Definisi 7 [Titik Tetap Tidak stabil]

Misalkan x* adalah titik tetap sebuah SPD dan x (t ) adalah solusi dengan nilai awal x(0) = x0 dengan x0 ≠ x* . Titik x* dikatakan titik tetap tidak stabil jika terdapat ε 0 > 0 , yang memenuhi sifat berikut: untuk setiap ε > 0 , 0 < ε < ε 0 , sedemikian sehingga jika

x* − x0 < ε maka

x* − x(t ) < ε 0 , untuk

setiap t > t0 . (Szidarovszky & Bahill, 1998) Definisi 8 [Nilai Eigen dan Vektor Eigen]

Misalkan A adalah matriks n × n , maka suatu vektor taknol x di dalam R n disebut vektor eigen dari A , jika untuk suatu skalar λ , yang disebut nilai eigen dari A , berlaku:

3

Ax = λ x , (1.4) vektor x disebut vektor eigen yang bersesuaian dengan nilai eigen λ . Untuk mencari nilai eigen dari matriks A yang berukuran n × n , maka persamaan (1.4) dapat dituliskan sebagai berikut: (1.5) (A − λ I)x = 0 , dengan I matriks identitas. Persamaan (1.4) mempunyai solusi tak nol jika dan hanya jika det (A − λ I) = 0 (1.6) Persamaan (1.6) disebut persamaan karakteristik. (Leon, 1998)

a. Setiap nilai eigen real adalah negatif ( λi < 0 untuk setiap i). b. Setiap komponen nilai eigen kompleks lebih kecil atau sama dengan nol, ( Re(λi ≤ 0) untuk setiap i). 2. Tak stabil jika a. Setiap nilai eigen real adalah positif ( λi > 0 untuk setiap i). b. Setiap komponen nilai eigen kompleks lebih besar dari nol, ( Re(λi > 0) untuk setiap i). 3. Sadel jika perkalian dua buah nilai eigen real sembarang adalah negatif ( λi λ j < 0 untuk i dan j sembarang).

Analisis Kestabilan Titik Tetap

Misalkan diberikan matriks A berukuran n × n sebagai berikut: ⎛ a11 L a1n ⎞ ⎟ ⎜ A = ⎜ M O M ⎟ . Dengan persamaan ⎜a L a ⎟ nn ⎠ ⎝ n1 karakteristik det (A − λ I) = 0 , dan I adalah matriks identitas, maka persamaan karakteristiknya menjadi: a1n ⎞ ⎛ a11 − λ L ⎟ ⎜ sedemikian O M ⎟=0, det ⎜ M ⎜ a ⎟ L ann − λ ⎠ ⎝ n1 sehingga akan diperoleh nilai eigen dari matriks tersebut. Analisis kestabilan titik tetap dilakukan untuk setiap nilai eigen yang diperoleh, maka akan ada 4 kasus sebagai berikut: 1. Jika nilai eigennya real dan berbeda tanda, maka titik tetap bersifat sadel. 2. Jika semua nilai eigennya real dan bertanda sama maka titik tetapnya merupakan simpul tak sejati (nodes). Jika bertanda positif maka nodes tak stabil. Jika bertanda negatif, maka nodes stabil. 3. Jika nilai eigennya merupakan complex conjugate dengan bagian real yang positif maka titik tetap bersifat spiral tak stabil. Jika bagian realnya negatif maka titik tetap bersifat spiral stabil. 4. Jika nilai eigen merupakan imajiner murni, maka titik tetap bersifat center yang selalu stabil. (Huntley & Johnson, 1983) Berdasarkan uraian di atas dapat disimpulkan bahwa kestabilan titik tetap mempunyai 3 perilaku sebagai berikut: 1. Stabil jika

Pelinearan

Untuk suatu SPD taklinear, analisis kestabilannya dilakukan melalui pelinearan. Misalkan diberikan SPD taklinear sebagai berikut: x& =f(x):U ⊂ R n → R n . (1.7) Dengan menggunakan ekspansi Taylor untuk *

suatu titik tetap x , maka persamaan (1.7) dapat ditulis sebagai berikut: x& = Ax + ϕ ( x) . (1.8) Persamaan tersebut merupakan SPD taklinear dengan A adalah matriks Jacobi, A = Df ( x * ) = Df ( x ) x = x * ⎡ ∂f1 ⎢ ∂x ⎢ 1 =⎢ M ⎢ ⎢ ∂f n ⎢⎣ ∂x1 ⎡ a11 ⎢ =⎢ M ⎢⎣ an1

∂f1 ⎤ ∂xn ⎥ ⎥ O M ⎥ ⎥ ∂f n ⎥ L ∂xn ⎥⎦ x = x* L a1n ⎤ ⎥ O M⎥ L ann ⎥⎦ L

dan ϕ ( x ) suku berorde tinggi yang bersifat lim ϕ ( x) = 0. Selanjutnya Ax pada x →0

persamaan (1.8) disebut pelinearan dari sistem taklinear persamaan (1.7) yang didapatkan dalam bentuk x& = Ax . (Verhulst, 1990) Model primitif penularan virus

Misalkan x = Populasi koloni manusia, y = Populasi virus, dengan asumsi sebagai berikut:

4

1. Laju kelahiran manusia α adalah konstan, 2. Infeksi virus menyebabkan peningkatan kematian terhadap penyakit, akibatnya g ( y ) > 0 , dengan g adalah fungsi sebarang, partikel virus 3. Perkembangbiakkan tergantung pada kehadiran manusia, 4. Dalam ketiadaan koloni manusia, partikel virus “mati“ atau menjadi tidak mampu bertahan hidup dengan laju γ. Persamaannya adalah sebagai berikut dx = [α − g ( y )]x, dt dy = β xy − γ y. dt Dalam Model ini, populasi dibagi kedalam kelas-kelas yang berbeda sesuai dengan tingkat kesehatan anggotanya. Pembagian kelas-kelas tersebut terdiri dari kelas individu yang rentan terserang penyakit S, yang terinfeksi I, dan kelas individu yang telah pulih R dari individu yang terjangkit penyakit tidak dalam jangka waktu lama, karena mereka telah pulih dan kebal, sudah diisolasi atau sudah mati. Jika penyakit memberikan sistem kekebalan sementara kepada penderitanya, atau penderita yang telah pulih dapat terjangkit kembali, maka individu dapat juga bergerak dari kelas ketiga ke kelas pertama. Jangka waktu epidemik dapat bervariasi dari minggu sampai tahun. Model klasik tentang teori epidemik telah dibuat oleh Kermack dan McKendrick (1927). Salah satu kasus khusus yang dipelajari digambarkan di bawah ini.

β S

v

I

R

Gambar 1 model SIR

γ β S

I

v

R

Gambar 2 model SIRS pada populasi tertutup

dt dI

S

δ

β I

v

δ

Gambar 3 model SIRS pada populasi terbuka

= δ NS − β SI − δ S + γ R

= β SI − δ I − vI dt dR = vI − γ R − δ R dt

γ δN

Gambar di atas menjelaskan laju perpindahan antara ketiga kelas dengan parameter β, yaitu laju penularan penyakit dan laju pemulihan v. Diasumsikan bahwa setiap pembagian kelas terdiri dari individu yang sama sehatnya atau sama sakitnya, dan tidak ada kelahiran atau kematian dalam populasi. (Dalam beberapa terminologi, model yang digambarkan pada gambar 1 disebut model SIR tanpa dinamika yang lebih rumit, karena perpindahan hanya dari kelas S ke I dan kemudian ke R). Total populasi N dibagi menjadi kelas yang rentan terserang infeksi (S), yang terinfeksi (I), dan yang pulih (sembuh) dari penyakit (R). Laju penularan penyakit, pemulihan, dan laju hilangnya sistem kekebalan masing-masing β, v, dan γ konstan. Persamaan yang dibuat Kermack dan McKendrick untuk penyakit yang diilustrasikan dalam gambar 1 adalah: dS = − β SI dt dI = β SI − vI dt dR = vI dt Dapat ditunjukkan bahwa total populasi N = S + I + R tidak berubah. Persamaan tersebut tak linear. S, I, dan R adalah fungsi dari waktu t. Misalkan sistem kekebalan tubuh dapat hilang sehingga individu yang telah pulih menjadi terserang penyakit kembali (gambar 2). Persamaan dapat dituliskan menjadi: dS = − β SI + γ R dt dI = β SI − vI dt dR = vI − γ R dt Model di atas disebut model SIRS pada populasi tertutup. Sedangkan persamaan pada gambar 3 yang merupakan model SIRS pada populasi terbuka adalah dS

R

δ

R0 merupakan rata-rata jumlah infeksi sekunder yang disebabkan oleh datangnya individu terinfeksi tunggal kedalam populasi

5

Analisis lebih lanjut mengenai model ini dapat dibuat kedalam suatu perhitungan, faktanya adalah bahwa total populasi: N =S+I+R tidak pernah berubah. Hal tersebut menjelaskan bahwa satu variabel, misalkan R dapat dihilangkan, sehingga model menjadi hanya dua persamaan.

N yang rentan terserang penyakit, atau bisa juga dikatakan R0 merupakan reproduksi dasar virus. Berikut adalah analisis untuk nilai R0: ¾ R0 < 1 : Virus tidak dapat bertahan hidup didalam populasi. ¾ R0 > 1 : Virus dapat bertahan hidup di dalam populasi.

PEMBAHASAN rentan sampai sembuh dan imunitas permanen diperoleh dinotasikan dengan R , atau dapat juga ditulis:

3.1. Perumusan Model

Model ini dideskripsikan oleh Waltman (1974). Pada model ini total populasi N (t ) yang merupakan fungsi dari waktu t dibagi menjadi tiga kelas yaitu; kelas individu yang rentan terserang penyakit dinotasikan dengan S , kelas individu yang terinfeksi yaitu individu yang mampu menularkan penyakit ke individu lain dinotasikan dengan I , dan kelas yang sembuh yaitu individu yang terserang penyakit dan mati, atau terserang penyakit dan sembuh, atau individu yang diisolasi dari

N (t ) = S (t ) + I (t ) + R(t )

Pada populasi, penyebaran penyakit menular mengikuti model dinamik SIRS (Susceptible > Infected -> Removed -> Susceptible) dan proses analisis ini telah banyak dilakukan pada manusia dan mamalia (Anderson and May, 1978; Giesecko, 1994; Cairns, 1995). Proses penggambarannya diperlihatkan pada gambar 5.

Birth

Susceptible (S)

β

Infected (I) α+ b

b

Dead

v

Remove (R) b

Dead

Dead

Gambar 4 Diagram Model SIRS

Model ini menggunakan beberapa asumsi, diantaranya: Total populasi N (t ) konstan. Perpindahan penyakit terjadi secara langsung, yang akan menghasilkan kepekatan S individu per satuan waktu dan kepekatan I individu per satuan waktu.

Adanya kepastian proporsi populasi, yaitu S dan sisanya I . Persatuan waktu, proporsi individu I membuat β I secara potensial merupakan interaksi penularan. Tingkat transmisi penularan penyakit ( β ) proporsional terhadap tingkat

6

pertemuan individu yang rentan terkena penyakit dengan penular yang dimodelkan oleh perkalian β SI , dimana

β

adalah koefisien transmisi penularan. Persamaan dideskripsikan mengikuti model deterministik yang memasukkan nilai konstan untuk setiap parameternya. Setiap individu pada populasi memiliki kesempatan yang sama untuk berinteraksi dengan individu yang terinfeksi atau penular. Individu yang sembuh dari infeksi dengan tingkat v , mereka telah memiliki imunitas untuk infeksi selanjutnya pada periode waktu tertentu. Semua parameter dan variable yang digunakan tak negatif. Sebagaimana dideskripsikan Anderson (1982), koefisien transmisi penularan β ditentukan oleh 2 faktor, yaitu frekuensi interaksi yang terjadi secara langsung untuk kepekatan hewan yang rentan terserang penyakit dan kemungkinan adanya kontak pada transmisi patogen. Banyaknya individu pada waktu t untuk setiap kelas dinyatakan sebagai S , I , dan R yang mana ketiganya merupakan fungsi dari t . Pada saat t diobservasi akan mengikuti persamaan:

dS = S0 + (−bS − β SI + vR) dt dI = β SI − (α + b + v) I dt dR = vI − bR dt

(1)

S0 adalah laju awal populasi yang rentan terserang penyakit. Selanjutnya digunakan skala berikut: S = N 0 A, I = N 0 B, dan R = N 0C. Sekarang, setiap fungsi merupakan fraksi dari populasi awal N 0 . Misal, A(t ) merupakan fraksi dari populasi awal N 0 yang rentan terserang penyakit. Dengan mensubstitusikan skala di atas ke dalam persamaan (1), (2) dan (3), maka akan diperoleh:

dA = A0 + (−bA − β N 0 AB + vC ) dt dB = β N 0 AB − (α + b + v) B dt dC = vB − bC dt

(4) (5) (6)

Selanjutnya akan dicari titik tetap untuk model (4), (5) dan (6) yang kemudian akan dianalisis kestabilan di sekitar titik tetap tersebut. 3.2. Titik tetap

Analisis titik tetap pada SPD sering digunakan untuk menentukan solusi yang tidak berubah terhadap waktu, yaitu pada saat

dx dy dz = 0, = 0 , dan = 0 . Titik tetap dt dt dt dari persamaan (4), (5) dan (6) akan diperoleh

(2) (3)

dimana S adalah banyaknya individu dalam populasi yang rentan terserang penyakit. I adalah banyaknya individu dalam populasi yang terinfeksi dan dapat menularkan penyakit. R adalah banyaknya individu dalam populasi yang sembuh dari penyakit. β adalah koefisien transmisi (laju) penularan penyakit. α adalah tingkat mortalitas individu yang terinfeksi (penular). v adalah laju penyembuhan. b adalah tingkat mortalitas alami. N adalah total populasi.

dengan menentukan

dC = 0, sehingga dt

dA dB = 0, = 0, dan dt dt diperoleh

persamaan

dibawah ini:

A0 + (−bA − β N 0 AB + vC ) = 0 β N 0 AB − (α + b + v) B = 0 vB − bC = 0

(7) (8) (9)

Dari hasil analisis (lihat lampiran), akan diperoleh dua jenis titik tetap, yaitu titik tetap tanpa penyakit dan titik tetap endemik. 3.2.1. Titik tetap tanpa penyakit

Titik tetap tanpa penyakit diperoleh ketika B = 0 . Dengan mensubstitusikan nilai tersebut kedalam persamaan (7), (8), dan (9) maka akan diperoleh nilai A dan C .

7

Sehingga diperoleh titik tetap dari persamaan

⎛ A0 ⎞ , 0, 0 ⎟ ⎝ b ⎠

tetap ini akan stabil dan sebaliknya jika

β N 0 A0

> (α + b + v)

(4), (5) dan (6) yaitu T1 = ⎜

λ3 > 0

Misalkan persamaan (4), (5) dan (6) dituliskan sebagai berikut:

maka titik tetapnya akan sadel. Kondisi stabil yang dipenuhi ketika

x( A, B, C ) = A0 + (−bA − β N 0 AB + vC ) y ( A, B, C ) = β N 0 AB − (α + b + v) B

z ( A, B, C ) = vB − bC Dengan melakukan pelinearan pada system persamaan di atas, maka akan diperoleh matriks Jacobi:

⎡ dx ⎢ dA ⎢ dy J =⎢ ⎢ dA ⎢ ⎢ dz ⎣⎢ dA

dx ⎤ dC ⎥ ⎥ dy ⎥ dC ⎥ ⎥ dz ⎥ dC ⎥⎦

dx dB dy dB dz dB

⎡ −b − β N 0 B = ⎢⎢ β N 0 B ⎢⎣ 0 Pelinearan

tetap

T1

v⎤ 0 ⎥⎥ −b ⎥⎦ akan

b

karena semua parameter tak negatif, maka λ1 = λ2 = −b < 0 sehingga kestabilan dititik ini tergantung pada nilai eigen

λ3 .

Jika

yang mana kondisi ini akan dipenuhi

b

merupakan

reproduksi dasar virus dalam populasi R0 . ketika

β N 0 A0 <1 b(α + b + v)

atau

β N 0 A0 >1 b(α + b + v)

atau

R0 > 1 maka populasi tidak stabil karena

λ1 = λ2 = −b β N 0 A0 λ3 = − (α + b + v)

β N 0 A0

β N0 A0 β N 0 A0 < 1. b(α + b + v) b(α + b + v)

Sebaliknya, ketika

sehingga akan diperoleh nilai eigen untuk matriks J0, yaitu:

ketika

< (α + b + v) dapat ditulis dalam

virus tidak dapat bertahan di dalam populasi.

menghasilkan matriks jacobi sebagai berikut: − β N 0 A0 ⎤ ⎡ v⎥ ⎢ −b b ⎥ ⎢ β N 0 A0 J0 = ⎢ 0 − (α + b + v) 0 ⎥ ⎢ ⎥ b ⎢ ⎥ v −b ⎥ ⎢0 ⎢⎣ ⎥⎦ Untuk memperoleh nilai eigen digunakan persamaan karakteristik det( J 0 − λ I ) = 0 ,

λ3 < 0

b

b

R0 < 1 yang merupakan kondisi stabil maka

v titik

β N 0 A0

Sehingga

−β N0 A β N 0 A − (α + b + v )

pada

yang berarti

< (α + b + v) maka titik

virus akan bertahan dalam populasi. Bisa disimpulkan bahwa titik tetap tanpa penyakit bersifat stabil asimtotik lokal jika dan hanya jika R0 < 1 dan tidak stabil jika dan hanya jika R0 > 1 . 3.2.2. Titik tetap endemik

Titik tetap endemic yaitu ketika ada individu yang terinfeksi dalam populasi. Dari persamaan (7), (8), dan (9) diperoleh titik tetap endemik T2 = ( x, y , z ) , dengan

x=

α +b+v β N0

b(α + b + v) ) β N0 y= b(α + b + v) − v 2 v(b(α + b + v) − β N 0 A0 ) z= β N 0 (v 2 − b(α + b + v)) b( A0 −

Pelinearan

pada

titik

tetap

T2

akan

menghasilkan matriks jacobi sebagai berikut:

8

⎡ ⎤ β bA0 N 0 − b 2 (α + b + v) −(α + b + v) v ⎥ ⎢ −b − 2 2 + + − b bv b v α ⎢ ⎥ ⎢ β bA0 N 0 − b 2 (α + b + v) ⎥ J1 = ⎢ 0 0⎥ 2 2 + + − b bv b v α ⎢ ⎥ v 0 −b ⎥ ⎢ ⎢ ⎥ ⎣ ⎦ Untuk memperoleh nilai eigen digunakan persamaan karakteristik −b −

β bA0 N 0 − b (α + b + v) − λ −(α + b + v) b 2 + bv + bα − v 2 β bA0 N 0 − b 2 (α + b + v) 0−λ b 2 + bv + bα − v 2

det( J1 − λ I ) = 0 ,

2

v

0

v =0

0 −b − λ

Sehingga diperoleh nilai eigen J1, yaitu:

λ1 → −

3 [d + d 2 + 4(− a 2 + 3bc)3 ]1/ 3 a 2(− a 2 + 3bc) − + 3b 3b[d + d 2 + 4(− a 2 + 3bc )3 ]1/ 3 3 3 2b

λ2 → −

(1 − i 3)[d + d 2 + 4(− a 2 + 3bc)3 ]1/ 3 a (1 + i 3)(− a 2 + 3bc ) + − 3b 3 3 4b[d + d 2 + 4(− a 2 + 3bc)3 ]1/ 3 6 3 2b

λ3 → −

(1 + i 3)[d + d 2 + 4(− a 2 + 3bc)3 ]1/ 3 a (1 − i 3)(− a 2 + 3bc) + − 3b 3 3 4b[d + d 2 + 4(− a 2 + 3bc)3 ]1/ 3 6 3 2b

dengan nilai a, b, c , dan d dalam lampiran. Berikut adalah tabel kondisi kestabilan dari dua titik tetap yang diperoleh.

β N 0 A0 b

< (α + b + v)

atau R0 < 1

β N 0 A0 b

Kondisi

> (α + b + v)

T1

T2

Simpul stabil

Spiral tak stabil

Sadel

Spiral stabil

atau R0 > 1 Tabel 1

Dari tabel di atas terlihat bahwa kondisi kestabilan dari dua titik tetap yang diperoleh saling bertentangan, meskipun jenis kestabilannya berbeda. Ketika titik tetap yang pertama stabil, titik tetap yang kedua tidak stabil dan ketika titik tetap yang pertama sadel, titik tetap yang kedua stabil. Selain itu, terlihat ada korelasi antara R0 dengan titik tetap endemik, sehingga titik tetap endemik dapat dinyatakan dalam R0 yaitu sebagai berikut. x=

A0 bR0

1 − 1) R0 y= 2 v − b(α + b + v) b(

1 − 1) R0 z= 2 v − b(α + b + v ) v(

Jika R0 < 1 maka x > A0 , karena x ,

y,

dan z positif sehingga akan mengakibatkan b(α + b + v) − v 2 > 0 . Artinya virus akan hilang dari populasi jika populasi yang rentan terserang penyakit melebihi populasi awal yang rentan terserang penyakit. Dan

9

sebaliknya jika

R0 > 1 maka

x < A0 ,

sehingga b(α + b + v) − v 2 < 0 , artinya virus akan bertahan dalam populasi ketika populasi yang rentan terserang penyakit kurang dari populasi awal yang rentan terserang penyakit. 3.3. Analisis kestabilan titik tetap ketika populasi awal yang terinfeksi 1 orang

a. Ketika R0 = 0.76

Kondisi R0 = 0.76 dipenuhi ketika β = 0.0002, dan v = 0.9, sehingga diperoleh dan titik tetap T1=(2.499,0,0) Orbit kestabilannya T2=(3.3,0.45,1.02). diperlihatkan di bawah ini. 8 α , β , b, v, p< = 8

se4

Di sini akan digambarkan bidang fase dari persamaan (4), (5), (6). Hal ini membutuh nilai awal untuk masing-masing parameter dan variabel . Misal populasi awal N (0) kita sebut N 0 . Misalkan hanya ada satu individu yang terinfeksi ( I (0) = 1 ), sehingga populasi lainnya rentan terserang penyakit Kondisi awal lainnya S (0) = N 0 − 1 . R(0) = 0 , sehingga dipenuhi:

N (0) = S (0) + I (0) + R(0) . Jika dituliskan secara menyeluruh, kondisi awal dari model ini adalah

S (0) = N (0) − 1 I (0) = 1 R(0) = 0 Dengan

S = N 0 A,

menggunakan skala berikut I = N 0 B, dan R = N 0C.

Persamaan kondisi awal di atas akan menjadi

A(0) = 1 −

1 1 , B(0) = , C (0) = 0. N0 N0

Pada proses penggambaran ini diambil nilai awal dari populasi N 0 = 2000 . Selain kondisi awal di atas, dibutuhkan juga nilai dari parameter yang belum diketahui. 3.3.1. Analisis kestabilan titik tetap untuk R0 < 1

Berikut adalah gambar bidang fase tiga dimensi dari persamaan (4), (5), dan (6) yang menunjukkan orbit kestabilannya. Dengan menggunakan Dynpac dari DynSys10. 71 pada Mathematica 5.2 ketika tingkat mortalitas alami b = 0.4, dan α = 0.02, 1 A(0) = 1 − , dengan nilai awal 2000 1 B(0) = , C (0) = 0. 2000

0.02 ,

0.00 ,

0.40 ,

0.90 ,

2000.00

0.0002 c 0.0001 0.0004 0 1

0.0002B 1.5 2 A

2.5

0

Gambar 5 Orbit kestabilan pada bidang A, B, dan C ketika β = 0.0002, dan v = 0.9. se4

8Bα , β , b, v, p< = 8

0.02 ,

0.00 ,

0.40 ,

0.90 , 2000.00

0005

0004

0003

0002

0001

A 0.5

1

1.5

2

2.5

3

3.5

Gambar 6 Proyeksi bidang fase pada bidang A, B ketika β = 0.0002, dan v = 0.9.

Dari gambar 6 terlihat bahwa orbit menuju titik T1 = (2.499, 0) sehingga T1 = (2.499, 0) berada pada bidang fase tersebut, sedangkan titik T2 = (3.3, 0.45) sangat jauh dari bidang fase. Hal ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node, dan T2 tak stabil.

10

se4 8cα, β, b, v, p<=8

0.02 ,

0.00 ,

0.40 ,

0.90 , 2000.00 8 α , β , b, v, p < = 8

se2

0.02 ,

0.00 ,

0.40 ,

0.80 ,

2000.00

0025

0002

0015

0.0002

0001

0.0001

c

0.0004

0005

0 1

0.0002B 1.5

A 0.5

1

1.5

2

2.5

3

2 A

3.5

Gambar 7 Proyeksi bidang fase pada bidang A, C ketika β = 0.0002, dan v = 0.9.

Gambar 7 yang merupakan proyeksi atau pencerminan orbit kestabilan dari gambar 5 terhadap sumbu A dan C menunjukkan bahwa T1 = (2.499, 0) bersifat stabil karena dituju oleh orbit kestabilan sehingga T1 = (2.499, 0) berada pada bidang fase. Jenis kestabilannya adalah node dan T2 = (3.3,1.02) tak stabil karena jauh dari bidang fase. se1 8c α , β , b, v, p< = 8 0035

0.02 ,

0.00 ,

0.40 ,

2.5

0

Gambar 9 Orbit kestabilan pada bidang A, B, dan C ketika β = 0.0002, dan v = 0.8. se2

8Bα , β , b, v, p < = 8

0.02 ,

0.00 ,

0.40 ,

0.80 , 2000.00

0005

0004

0003

0002

0.90 , 2000.00

0001 0003

A 0.5

0025

1

1.5

2

2.5

3

3.5


0002

0015

0001

0005

B 0.0001

0.0002

0.0003

0.0004

0.0005

Gambar 8 Proyeksi bidang fase pada bidang B, C ketika β = 0.0002, dan v = 0.9.

Gambar 8 merupakan proyeksi atau pencerminan dari orbit kestabilan tiga dimensi gambar 5 terhadap bidang B, C menunjukkan bahwa T1 = (0, 0) dituju oleh bidang fase. Hal ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node dan T2 tak stabil karena jauh dari bidang fase.

Dari gambar terlihat bahwa titik tetap T1 = (2.499, 0) dituju oleh bidang fase sehingga berada pada bidang fase dan titik tetap T2 = (3.05, 0.58) jauh dari bidang fase. Ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node dan T2 tak stabil. Selain itu terlihat pula bahwa dengan menurunkan nilai parameter v kurva menjadi lebih curam jika dibandingkan dengan kurva sebelumnya. se2

8c α , β , b, v, p< = 8

0.02 ,

0.00 ,

0.40 ,

0.80 , 2000.00

2.5

3

0025

0002

0015

b. Ketika R0 = 0.82 0001

Berikut akan dianalisis untuk nilai R0 yang berbeda dengan mengubah parameter v menjadi 0.8 sehingga diperoleh nilai R0 yang lebih besar yaitu R0 = 0.82 dengan titik tetap T1=(2.499,0,0) dan T2=(3.05,0.58,1.16).

0005

A 0.5

1

1.5

2

3.5


11

Dari gambar terlihat bahwa titik tetap T1 = (2.499, 0) dituju oleh bidang fase sehingga berada pada bidang fase dan titik tetap T2 = (3.05,1.16) jauh dari bidang fase. Ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node dan T2 tak stabil. Selain itu terlihat pula bahwa jika dibandingkan dengan kurva sebelumnya, dengan menurunkan nilai parameter v yang mengakibatkan nilai R0 menjadi lebih besar gambar bidang fase menjadi semakin mampat. se2 0035

8c α , β , b, v, p< = 8

se1

8 α , β , b, v, p < = 8

0.02 ,

0.00 ,

0.40 ,

0.90 ,

2000.00

0.0002 0.00015 c 0.0001 0.00005 0.0004 0 1

B 0.0002 1.5

0.02 ,

0.00 ,

0.40 ,

2

0.80 , 2000.00

A

2.5

0

0003

Gambar 13 Orbit kestabilan pada bidang A, B, dan C ketika β = 0.0001, dan v = 0.9.

0025

0002 se1

0015

8Bα , β , b, v, p < = 8

0.02 ,

0.00 ,

0.40 ,

0.90 ,

2000.00

0005

0001

0004

0005

B 0.0001

0.0002

0.0003

0.0004

0.0005


Dari gambar terlihat bahwa titik tetap T1 = (0, 0) dituju oleh bidang fase sehingga berada pada bidang fase dan titik tetap T2 = (0.58,1.16) jauh dari bidang fase. Ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node dan T2 tak stabil. Selain itu terlihat pula bahwa dengan menurunkan nilai parameter v yang mengakibatkan nilai R0 menjadi lebih besar, gambar bidang fase menjadi semakin sempit terlihat dari ujung bidang fase yang menurun. c.

Ketika R0 = 0.38

Berikut akan dianalisis untuk nilai R0 yang berbeda dengan mengubah nilai parameter β menjadi 0.0001, sehingga nilai R0 menjadi lebih kecil. Dengan perubahan nilai β tersebut mengakibatkan nilai titik tetapnya juga berubah menjadi T1=(2.499,0,0) dan T2=(6.6,2.33,5.24).

0003

0002

0001

A 0.5

1

1.5

2

2.5

3

3.5


Dari gambar terlihat bahwa titik tetap T1 = (2.499, 0) dituju oleh bidang fase sehingga berada pada bidang fase dan titik tetap T2 = (6.6, 2.33) jauh dari bidang fase. Ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node dan T2 tak stabil. Selain itu terlihat pula bahwa dengan menurunkan nilai parameter β kurva menjadi lebih landai dibandingkan dengan kurva sebelumnya.

12

se1 8c α, β, b, v, p<=8

0.02 ,

0.00 ,

0.40 ,

0.90 , 2000.00

3.3.2. Analisis kestabilan titik tetap untuk R0 > 1

Berikut akan dianalisis orbit kestabilan untuk nilai R0 > 1 dengan tingkat mortalitas alami b = 0.4, α = 0.02, dengan nilai awal 1 A(0) = 1 − , yang sama yaitu 2000 1 B(0) = , dan C (0) = 0. 2000

00025

0002

00015

0001

00005

a. Ketika R0 = 1.09 A 0.5

1

1.5

2

2.5

3

3.5


Dari gambar terlihat bahwa titik tetap T1 = (2.499, 0) dituju oleh bidang fase sehingga berada pada bidang fase dan titik tetap T2 = (6.6,5.24) jauh dari bidang fase. Ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node dan T2 tak stabil. Selain itu terlihat pula bahwa dengan menurunkan nilai parameter β puncak kurva menjadi lebih rendah dibandingkan dengan kurva dengan nilai parameter β = 0.0002, dan v = 0.9. se1 8cα, β, b, v, p<=8 00035

0.02 ,

0.00 ,

0.40 ,

Kondisi R0 = 1.09 dipenuhi ketika β = 0.0002, dan v = 0.5, dengan titik tetap dan T2=(2.3,0.27,0.34). T1=(2.499,0,0) Diperoleh orbit kestabilan yang digambarkan dengan bidang fase pada gambar di bawah. 8α , β, b, v, p< =8

se1

0.02 ,

0.00 ,

0.40 ,

0.50 ,

2000.00

0.004

c 0.002

0.004

0.90 , 2000.00

0

0.003 0.002 B

1 1.5

0.001

2

0003

A

2.5

00025

Gambar 17 Bidang fase tiga dimensi (A, B, C) ketika β = 0.0002, dan v = 0.5

0002

se1 B8 α , β , b, v, p < = 8 006

00015

0001

0.02 ,

0.00 ,

0.40 ,

0.50 , 2000.00

005

00005 004

B 0.0001

0.0002

0.0003

0.0004

0.0005


003

002

001

Dari gambar terlihat bahwa titik tetap T1 = (0, 0) dituju oleh bidang fase, sedangkan titik tetap T2 = (2.33,5.24) jauh dari bidang fase. Ini menunjukkan bahwa T1 stabil dengan jenis kestabilan node dan T2 tak stabil. Selain itu terlihat pula bahwa dengan menurunkan nilai parameter β ujung dari bidang fase menjadi lebih tinggi dibandingkan dengan bidang fase dengan nilai parameter β = 0.0002, dan v = 0.9.

A 0.5

1

1.5

2

2.5

Gambar 18 Proyeksi bidang fase pada bidang A, B ketika β = 0.0002, dan v = 0.5

Gambar di atas merupakan proyeksi atau pencerminan dari orbit kestabilan tiga dimensi pada gambar 17 terhadap sumbu A, B. Bidang fase bergerak dari titik awal (0.9995, 0.0005) mendekati titik tetap T1 = (2.499, 0) tapi tidak berhenti di T1, karena bidang fase tersebut

13

berbelok menuju ke suatu nilai dari B, mungkin menuju ke titik tetap T2 = (2.3, 0.27) . se1 c8 α , β , b, v, p < = 8 006

0.02 ,

0.00 ,

0.40 ,

0.50 ,

2000.00

005

004

b. Ketika R0 = 1.22

Berikut dilakukan analisis untuk nilai v = 0.4, sehingga menghasilkan nilai R0 yang berbeda yang lebih besar. Kondisi tersebut menghasilkan titik tetap yang berbeda pula yaitu T1=(2.499,0,0) dan T2=(2.05,0.43,0.43). Di bawah ini adalah gambar bidang fase dalam bentuk 3 dimensi.

003 8 α , β , b, v, p < = 8

se1

0.02 ,

0.00 ,

0.40 ,

0.40 ,

2000.00

002

001

A 0.5

1

1.5

2

2.5

Gambar 19 Proyeksi bidang fase pada bidang A, C ketika β = 0.0002, dan v = 0.5

0.2 c 0.1

Gambar di atas merupakan proyeksi atau pencerminan dari orbit kestabilan tiga dimensi pada gambar 17 terhadap sumbu A, C. Bidang fase bergerak dari titik awal (0.9995, 0) mendekati titik tetap T1 = (2.499, 0) tapi tidak berhenti di T1, karena bidang fase tersebut berbelok menuju ke suatu nilai dari B, mungkin menuju ke titik tetap T2 = (2.3, 0.34) . se1 c8 α , β , b, v, p < = 8 006

0.2

0 1

B 0.1

1.5 2 A

2.5

0

Gambar 21 Bidang fase tiga dimensi (A, B, C) ketika β = 0.0002, dan v = 0.4

se1 B8 α , β , b, v, p < = 8

0.02 ,

0.00 ,

0.40 ,

0.40 ,

2000.00

0.3

.25

0.02 ,

0.00 ,

0.40 ,

0.50 ,

2000.00 0.2

005 .15

004

0.1

.05

003

A 0.5

002

1

1.5

2

2.5


001

0.001

0.002

0.003

0.004

0.005

B 0.006


Gambar di atas merupakan proyeksi atau pencerminan dari orbit kestabilan tiga dimensi pada gambar 17 terhadap sumbu B, C. Bidang fase bergerak dari titik awal (0.0005, 0) menjauhi titik tetap T1 = (0, 0) , dan mungkin menuju ke titik tetap T2 = (0.27, 0.34) . Dari ketiga gambar di atas terlihat bahwa bidang fase menjauhi titik tetap T1, dan mungkin menuju ke titik tetap T2. Hal ini menunjukkan bahwa pada kondisi ini T1 sadel.

Gambar di atas merupakan proyeksi atau pencerminan dari orbit kestabilan tiga dimensi pada gambar 21. Bidang fase bergerak dari titik awal (0.9995, 0.0005) mendekati titik tetap T1 = (2.499, 0) tapi tidak berhenti di T1. Bidang fase terus bergerak naik menuju ke titik tetap T2 = (2.05, 0.43) . Ini menunjukkan bahwa T1 sadel. Selain itu, dengan berkurangnya nilai v mengakibatkan bidang fase semakin rapat dengan sumbu absis.

14

se1 c8 α , β , b, v, p < = 8

0.02 ,

0.00 ,

0.40 ,

0.40 ,

2000.00

c. Ketika R0 = 1.83

Berikut akan dilakukan analisis untuk perubahan nilai β yaitu dengan menaikkan nilai β menjadi 0.0003 dengan nilai v = 0.4, dan parameter yang lain tetap, menghasilkan nilai R0 yang lebih besar yaitu R0 = 1.83 , dan titik tetapnya menjadi T1=(2.499,0,0) dan T2=(1.37,1.08,1.08). Di bawah ini adalah gambar bidang fase dalam 3 dimensi.

0.3

.25

0.2

.15

0.1

.05

A 0.5

1

1.5

2

2.5


Gambar di atas merupakan proyeksi atau pencerminan dari orbit kestabilan tiga dimensi pada gambar 21 terhadap sumbu A, C. Bidang fase bergerak dari titik awal (0.9995, 0) mendekati titik tetap T1 = (2.499, 0) tapi tidak berhenti di T1. Bidang fase terus bergerak naik menuju ke titik kestabilan T2 = (2.05, 0.43) . Ini menunjukkan bahwa T1 sadel. Selain itu, dengan berkurangnya nilai v mengakibatkan bidang fase semakin rapat dengan sumbu absis. se1 c8 α , β , b, v, p < = 8

0.02 ,

0.00 ,

0.40 ,

0.40 ,

2000.00

8 α , β , b, v, p < = 8

se1

0.00 ,

0.40 ,

0.40 ,

2000.00

1 0.75 c 0.5 0.25

1 0.75

0 1

0.5 B 1.5

0.25 2

A

0

Gambar 25 Bidang fase 3 dimensi (A, B, C) ketika β = 0.0003, dan v = 0.4 se1 B 8 α , β , b, v, p< = 8

0.3

0.02 ,

0.02 ,

0.00 ,

0.40 ,

0.40 , 2000.00

.25 1

0.2 0.8

.15 0.6

0.1 0.4

.05

B 0.05

0.1

0.15

0.2

0.25

0.2

0.3

Gambar 24 Proyeksi bidang fase pada bidang B, C ketika β = 0.0002, dan v = 0.4

Gambar di atas merupakan proyeksi atau pencerminan dari orbit kestabilan tiga dimensi pada gambar 21 terhadap sumbu B, C. Bidang fase bergerak dari titik awal (0.0005, 0) menjauhi titik kestabilan T1 = (0, 0) menuju titik kestabilan di T2 = (0.43, 0.43) . Selain itu, terlihat bahwa gambar bidang fase menjadi lebih landai. Dari ketiga gambar di atas terlihat bahwa bidang fase menjauhi titik tetap T1, dan mungkin menuju ke titik tetap T2. Hal ini menunjukkan bahwa pada kondisi ini T1 sadel. Selain itu, dengan penurunan nilai v mengakibatkan peningkatan untuk nilai B dan C.

A 0.5

1

1.5

2

2.5


Gambar di atas merupakan proyeksi dari gambar 25 terhadap sumbu A, B. Pada gambar di atas titik awal pergerakkan bidang fase tidak terlihat dengan jelas, hal ini karena skala nilai yang dihasilkan lebih besar daripada ketika β = 0.0002. Di sini jelas terlihat bahwa bidang fase menuju titik tetap T2 di (1.37,1.08).

15

se1 c 8 α , β , b, v, p < = 8

0.02 ,

0.00 ,

0.40 ,

0.40 ,

2000.00

1

0.8

0.6

0.4

0.2

A 0.5

1

1.5

2

2.5


Gambar di atas merupakan proyeksi dari gambar 25 terhadap sumbu A, C. Pada gambar di atas titik awal pergerakkan bidang fase tidak terlihat dengan jelas, hal ini karena skala nilai yang dihasilkan lebih besar daripada ketika β = 0.0002. Di sini jelas terlihat bahwa bidang fase menuju titik tetap T2 di (1.37,1.08). se1 c 8 α , β , b, v, p< = 8

0.02 ,

0.00 ,

0.40 ,

0.40 , 2000.00

1

0.8

0.6

0.4

0.2

B 0.2

0.4

0.6

0.8

1

Gambar 28 Proyeksi bidang fase pada bidang B, C ketika β = 0.0003, dan v = 0.4

Gambar di atas merupakan proyeksi dari gambar 25 terhadap sumbu B, C. Pada gambar di atas titik awal pergerakkan bidang fase terlihat seolah-olah dari titik (0,0) padahal sebenarnya bidang fase bergerak dari titik (0.0005,0), hal ini karena skala nilai yang dihasilkan lebih besar daripada ketika β = 0.0002.

Dari ketiga proyeksi di atas dapat dilihat bahwa dengan bertambahnya nilai β titik tetapnya semakin besar, begitu pula dengan nilai R0. Dengan semakin besarnya nilai titik tetap mengakibatkan gambar bidang fase lebih jelas pada ujungnya tapi kurang jelas kondisi awalnya.

3.4. Dinamika populasi penularan virus influenza

Untuk mengamati pengaruh masuknya patogen atau virus ke dalam populasi pada waktu tertentu maka diperlukan kurva yang menunjukkan pengaruh masuknya patogen kedalam populasi hubungannya dengan periode waktu. Hal ini membutuh nilai awal untuk masing-masing parameter dan variabel . Pada proses penggambaran ini diambil nilai awal dari populasi N 0 = 2000 . Selain kondisi awal di atas, dibutuhkan juga nilai dari parameter yang belum diketahui. Dalam karya tulis ini dianalisis dinamika populasi untuk dua kondisi, yaitu ketika R0 < 1 dimana populasi akan stabil karena virus akan hilang dari populasi dan ketika R0 > 1 dimana virus akan bertahan dalam populasi. 3.4.1. Dinamika populasi untuk R0 < 1

Kondisi R0 < 1 adalah kondisi dimana populasi akan stabil menuju musnahnya virus dari populasi. Dalam proses penggambarannya saya menggunakan Mathematica 6.0 yang dievaluasi ketika tingkat mortalitas alami b = 0.4, dan tingkat mortalitas populasi yang terinfeksi α = 0.02, 1 A(0) = 1 − , dengan nilai awal 2000 1 B(0) = , C (0) = 0. 2000 Untuk nilai R0 < 1 ini akan dianalisis untuk tiga kondisi dengan mengubah nilai v dan β sehingga diperoleh nilai R0 yang berbeda yaitu R0=0.76, R0=0.81 dan R0 = 0.38 . a. Ketika R0 = 0.76 Kondisi R0=0.76 dipenuhi ketika β = 0.0002, dan v = 0.9 , sehingga diperoleh dan titik tetap T1=(2.499,0,0) T2=(3.3,0.45,1.02).

16

A 2.5

2.0

1.5

1.0

0.5

0.0

10

5

0

15

20

t

25

Gambar 29 Dinamika populasi A terhadap waktu t

Gambar di atas menyatakan bahwa populasi yang rentan terserang penyakit mengalami peningkatan yang signifikan pada waktu yang cukup singkat, tapi kemudian stabil pada angka 2.499 setelah 8 hari. Peningkatan tersebut terjadi karena pada kondisi ini x > A0 .

Gambar di atas menyatakan bahwa pada hari ke-1 laju kesembuhan meningkat tajam, kemudian pada hari ke-2 peningkatan yang terjadi tidak signifikan seperti hari pertama. Pada hari ketiga, laju kesembuhan mengalami penurunan yang tajam menuju kestabilan pada angka nol setelah hari ke-18. Jumlah populasi yang sembuh meningkat karena ketika banyak populasi yang terinfeksi maka populasi yang sembuh pun akan banyak hal ini karena laju kesembuhan yang cukup besar yaitu mendekati 1. Namun setelah mencapai tingkat tertentu populasi yang sembuh akan menurun, hal ini seiring dengan menurunnya jumlah populasi yang juga menurun. A,B,C 2.5

2.0

1.5

1.0

B 0.0006

0.5

0.0005 0.0

0.0004

0.0003

0.0002

0.0001

0.0000

0

5

10

15

20

25

t

Gambar 30 Dinamika populasi B terhadap waktu t

Gambar di atas menyatakan bahwa jumlah populasi yang terinfeksi mengalami penurunan yang sangat signifikan pada waktu yang sangat singkat. Kemudian mencapai kestabilan pada angka 0 dihari ke-11. Penurunan ini terjadi karena pada kondisi ini virus tidak dapat bertahan dalam populasi. c 0.0004

0.0003

0.0002

0.0001

0.0000

0

5

10

15

20

25

t

Gambar 32 Dinamika populasi A, B, C terhadap waktu t

0

5

10

15

20

25

t

Gambar 31 Dinamika populasi C terhadap waktu t

Gambar di atas merupakan gabungan kurva A, B, dan C. Pada gambar di atas terlihat bahwa kurva B dan C berupa garis horizontal yang berada disumbu absis dan nilainya mendekati nol. Hal ini berarti jumlah populasi yang terinfeksi sedikit sehingga jumlah populasi yang sembuh pun sedikit. Inilah yang menyebabkan virus tidak dapat bertahan dalam populasi atau virus akan musnah dari populasi. b. Ketika R0 = 0.81

Berikut adalah dinamika populasi ketika koefisien penularan penyakit β tetap dan laju kesembuhan v diturunkan menjadi 0.8. Hal ini menaikkan nilai R0 menjadi 0.81 dengan dan titik tetap T1=(2.499,0,0) T2=(3.05,0.58,1.16).

17

A

Populasi yang terinfeksi meningkat pada hari pertama secara tajam, kemudian pada hari kedua peningkatan yang terjadi hanya sedikit sebelum kemudian menurun secara tajam menuju kestabilan di nol pada hari ke-24. Kurva ini lebih landai dibandingkan kurva yang sama untuk nilai R0 yang lebih kecil. Hal ini menunjukkan bahwa kecepatannya menurun sehingga lebih lambat menuju kestabilan.

2.5

2.0

1.5

1.0

0.5

A,B,c 0.0

0

5

10

15

20

25

t

2.5


Kurva A di atas lebih landai dari kurva A sebelumnya. Hal ini menunjukkan bahwa kecepatan menuju kestabilan lebih kecil, meskipun populasi mencapai kestabilan dengan tingkat yang sama yaitu pada 2.499.

2.0

1.5

1.0

B 0.0006

0.5 0.0005

0.0 0.0004

0.0003

0.0002

0.0001

0.0000

0

5

10

15

20

t

25


0

5

10

15

20

25

t


Populasi penular terus mengalami penurunan menuju kestabilan pada angka nol. Tetapi untuk kondisi ini populasi lebih lambat menuju kestabilan. Dimana untuk nilai R0=0.76 populasi yang terinfeksi mencapai kestabilan pada hari ke-11, tetapi untuk R0=0.81 populasi yang terinfeksi mencapai kestabilan pada hari ke-14. Kurva di atas lebih landai daripada kurva B sebelumnya, ini menunjukkan bahwa kecepatan menuju kestabilan menurun.

Gambar di atas merupakan dinamika populasi untuk A, B, dan C. Kurva B dan C berupa garis horizontal karena nilainya yang jauh lebih kecil dari A dan keduanya memiliki nilai yang tidak jauh berbeda. Karena jumlah populasi yang terinfeksi sedikit sehingga populasi yang sembuh juga kecil, hal ini akan mengakibatkan virus musnah dari populasi. c. Ketika R0 = 0.38

Berikut adalah dinamika populasi ketika koefisien penularan penyakit β diturunkan menjadi 0.0001 dan laju kesembuhan v tetap yaitu 0.9. Hal ini menurunkan nilai R0 menjadi 0.38 sehingga diperoleh titik tetap T1=(2.499,0,0) dan T2=(6.6,2.33,5.24).

c 0.0004

A 2.5

0.0003

2.0

1.5

0.0002

1.0

0.0001 0.5

0.0000

0

5

10

15

20

25

t


0.0

0

5

10

15

20

25

t


18

Kurva A di atas lebih curam dari kurva A ketika β = 0.0002, dan v = 0.9 . Hal ini menunjukkan bahwa kecepatan menuju kestabilan lebih besar sehingga populasi akan lebih cepat menuju kestabilan, meskipun pada tingkat kestabilan yang sama yaitu 2.499. B 0.0006

A,B,c 2.5

2.0

1.5

1.0

0.0005

0.5

0.0004

0.0

0

5

10

15

20

25

t


0.0003

0.0002

0.0001

0.0000

0

5

10

15

20

25

t


Kurva B di atas lebih curam dari kurva B ketika β = 0.0002, dan v = 0.9 . Hal ini menunjukkan bahwa kecepatan menuju kestabilan meningkat sehingga populasi akan lebih cepat menuju kestabilan. Ketika β = 0.0002, dan v = 0.9 populasi yang terinfeksi mencapai kestabilan pada hari ke11, sedangkan pada kondisi β = 0.0003, populasi yang terinfeksi mencapai kestabilan pada hari ke-6. c 0.0004

0.0003

Gambar di atas merupakan dinamika populasi A, B, dan C. Kurva B dan C berupa garis horizontal karena nilainya yang jauh lebih kecil dari A dan keduanya memiliki nilai yang tidak jauh berbeda. Karena jumlah populasi yang terinfeksi sedikit sehingga populasi yang sembuh juga kecil, hal ini akan mengakibatkan virus musnah dari populasi. 3.4.2. Dinamika populasi untuk R0 > 1

Kondisi R0 > 1 adalah kondisi ketika virus akan bertahan dalam populasi. Dalam proses penggambarannya saya menggunakan Mathematica 6.0 yang dievaluasi ketika tingkat mortalitas alami b = 0.4, dan 1 , α = 0.02, dengan nilai awal A(0) = 1 − 2000 1 B(0) = , C (0) = 0. 2000 Untuk nilai R0 > 1 ini akan dianalisis

0.0002

untuk tiga kondisi yang berbeda, yaitu ketika R0=1.09, R0=1.22 dan R0=1.83.

0.0001

a. Ketika R0=1.09

0.0000

0

5

10

15

20

25

t


Kurva C di atas lebih rendah dari kurva C ketika β = 0.0002, dan v = 0.9 . Hal ini menunjukkan bahwa kecepatan menuju kestabilan meningkat sehingga populasi akan lebih cepat menuju kestabilan. Ketika β = 0.0002, dan v = 0.9 populasi yang sembuh mencapai kestabilan pada hari ke-20, sedangkan pada kondisi β = 0.0003, populasi yang sembuh mencapai kestabilan pada hari ke-15.

Berikut adalah dinamika populasi ketika R0=1.09, yang mana kondisi ini dipenuhi ketika nilai parameter β = 0.0002, dan v = 0.5, sehingga diperoleh titik tetap T1=(2.499,0,0) dan T2=(2.3,0.27,0.34).

19

A

Pada gambar di atas, jumlah populasi yang sembuh pada awal pengamatan adalah nol dan terus mengalami peningkatan. Peningkatan yang signifikan terjadi setelah hari ke-70. Kemudian populasi yang sembuh stabil pada 0.34 setelah hari ke-190.

2.5

2.0

1.5

A,B,c

1.0

2.5 0.5

2.0 0.0

0

50

100

150

200

t


1.5

Pada grafik di atas, populasi yang rentan pada awal pengamatan mengalami kenaikan yang sangat tajam. Kemudian stabil pada hari ke-17 di 2.499, tapi itu tidak bertahan lama. Grafik kemudian menurun dan mencapai kestabilan kembali di 2.3 yang merupakan titik kestabilan yang ke-2 setelah hari ke-160.

1.0

B 0.35

0.0

0

50

100

150

200

t


Gambar di atas merupakan dinamika populasi A, B, dan C terhadap waktu t. Kurva A mengalami peningkatan yang signifikan pada awal pengamatan sebelum kemudian stabil. Sementara kurva B dan C terus mengalami kenaikkan tapi dengan rentan yang sangat kecil. Dari ketiga kurva, dapat diketahui bahwa populasi menuju kestabilan titik tetap T2.

0.30 0.25 0.20 0.15 0.10 0.05 0.00

0.5

b. Ketika R0 = 1.22 0

50

100

150

200

t


Pada gambar di atas, populasi yang terinfeksi pada awal pengamatan sangat kecil bahkan terlihat seperti nol. Tapi terus mengalami peningkatan dan peningkatan yang signifikan terjadi setelah hari ke-60. Kemudian populasi yang terinfeksi stabil pada 0.27 setelah hari ke-170. c

Dengan mengubah nilai parameternya, yaitu v menjadi 0.4, dengan nilai β yang tetap maka diperoleh nilai R0 = 1.22 dan titik tetapnya adalah T1=(2.499,0,0) dan T2=(2.05,0.43,0.43). A 2.5

2.0

0.35 0.30

1.5 0.25

1.0

0.20 0.15

0.5

0.10

0.0 0.05 0.00

0

50

100

150

200

t

Gambar 45 Dinamika populasi A terhadap waktu t 0

50

100

150

200

t


20

Pada grafik di atas, populasi yang rentan pada awal pengamatan mengalami kenaikan yang sangat tajam. Kemudian stabil pada hari ke-11 di 2.499, tapi itu tidak bertahan lama. Grafik kemudian menurun dan mencapai kestabilan kembali di 2.05 yang merupakan titik kestabilan yang ke-2 setelah hari ke-60. Grafik ini lebih cepat menuju kestabilan dibandingkan dengan ketika nilai v=0.5, tapi dengan kestabilan yang lebih kecil. 0.5

B

curam dari grafik yang sejenis. Ini menunjukkan bahwa kecepatan menuju kestabilannya meningkat. A,B,c 2.5

2.0

1.5

1.0

0.5

0.4

0.0

0.3

0

50

100

150

200

t


0.2

0.1

0.0

0

50

100

150

200

t


Gambar di atas memperlihatkan bahwa populasi penular terus mengalami peningkatan, dan peningkatan yang signifikan terjadi setelah hari ke-30. Grafik kemudian stabil pada 0.43 setelah hari ke-100. Dibandingkan dengan grafik sebelumnya, grafik ini lebih curam sehingga populasi lebih cepat menuju kestabilannya. c

Gambar di atas menunjukkan hubungan antara A, B, dan C terhadap waktu t. Kurva A mengalami peningkatan yang signifikan pada awal pengamatan sebelum kemudian stabil. Sementara kurva B dan C terus mengalami kenaikkan tapi dengan rentan yang sangat kecil. Dari ketiga kurva, dapat diketahui bahwa populasi menuju kestabilan titik tetap T2. Selain itu terlihat bahwa ketiga kurva lebih curam dari kurva sebelumnya. Jarak antara kurva populasi yang rentan dengan kedua kurva lainnya lebih dekat daripada kurva sebelumnya. c.Ketika R0 = 1.83

Dengan mengubah nilai parameternya, yaitu β menjadi 0.0003, dan v tetap 0.4, maka diperoleh nilai R0 = 1.83 dan titik tetapnya dan adalah T1=(2.499,0,0) T2=(1.37,1.08,1.08).

0.5

0.4

0.3 A 2.5

0.2 2.0

0.1 1.5

0.0

0

50

100

150

200

t

1.0

Gambar 47 Dinamika populasi C terhadap waktu t 0.5

Gambar di atas memperlihatkan bahwa jumlah populasi yang sembuh pada awal pengamatan adalah nol dan terus mengalami peningkatan. Peningkatan yang signifikan terjadi setelah hari ke-30 sampai hari ke-100. Kemudian stabil pada 0.43 setelah hari ke110. Selain itu terlihat bahwa grafik ini lebih

0.0

0

50

100

150

200

t


Pada grafik di atas, populasi yang rentan pada awal pengamatan mengalami kenaikan yang sangat tajam mencapai titik

21

kesetimbangan yang pertama yaitu 2.499. Tapi kemudian menurun dan stabil setelah hari ke-40 di 1.37 yang merupakan titik kestabilan yang ke-2. Grafik ini lebih cepat menuju kestabilan dibandingkan dengan ketika nilai β = 0.0002, tapi dengan kestabilan T2 yang lebih kecil.

A,B,c 2.5

2.0

1.5

1.0 B

1.4

0.5

1.2

0.0

1.0

0

50

100

150

200

t

Gambar 52 Dinamika populasi A, B, dan C terhadap waktu t

0.8

0.6 0.4 0.2 0.0

0

50

100

150

200

t


Gambar di atas memperlihatkan bahwa populasi penular terus mengalami peningkatan, dan peningkatan yang signifikan terjadi setelah hari ke-10. Grafik kemudian stabil pada 1.08 setelah hari ke-50. Dibandingkan dengan grafik sebelumnya, grafik ini lebih curam sehingga populasi lebih cepat menuju kestabilan. c 1.4 1.2 1.0

0.8 0.6 0.4 0.2

0.0

0

50

100

150

200

t


Gambar di atas memperlihatkan bahwa jumlah populasi yang sembuh pada awal pengamatan adalah nol dan terus mengalami peningkatan. Peningkatan yang signifikan terjadi setelah hari ke-14. Kemudian stabil pada 1.08 setelah hari ke-50. Selain itu terlihat bahwa grafik ini lebih curam dari grafik yang sejenis. Ini menunjukkan bahwa kecepatan menuju kestabilannya meningkat.

Gambar di atas menunjukkan hubungan antara A, B, dan C terhadap waktu t. Kurva A mengalami peningkatan yang signifikan pada awal pengamatan sebelum kemudian stabil. Sementara kurva B dan C terus mengalami kenaikkan sebelum kemudian mencapai kondisi yang konstan. Dari ketiga kurva, dapat diketahui bahwa populasi menuju kestabilan titik tetap T2. Selain itu terlihat bahwa ketiga kurva lebih curam dari kurva sebelumnya yang berarti bahwa kecepatan menuju kestabilannya bertambah. Jarak antara kurva populasi yang rentan dengan kedua kurva lainnya lebih dekat daripada kurva sebelumnya, ini menunjukkan bahwa banyak populasi yang terinfeksi sehingga populasi yang rentan terserang penyakit berkurang meskipun populasi yang sembuh pada kondisi kesetimbangan sama dengan populasi yang terinfeksi. Selain itu, dapat pula dikatakan semakin dekat jarak ketiga kurva populasi semakin cepat menuju kesetimbangan T2. Berikut akan dibandingkan ketika populasi awal yang terinfeksi 1% dari populasi yaitu sebanyak 20 ( I (0) = 20 ), dengan ketika populasi awal yang terinfeksi 10% dari populasi yaitu sebanyak 200 ( I (0) = 200 ). Ketika R0 < 1 parameter α = 0.02, lainnya yaitu, b = 0.4, β = 0.0002, dan v = 0.9. Ketika R0 > 1

α = 0.02, parameter lainnya b = 0.4, β = 0.0003, dan v = 0.5. Ketika populasi yang terinfeksi 1% dari populasi maka populasi lainnya rentan terserang penyakit S (0) = N 0 − 20 . Kondisi awal lainnya R(0) = 0 , sehingga dipenuhi: N (0) = S (0) + I (0) + R(0) . Jika dituliskan secara menyeluruh, kondisi awal dari model ini adalah

22

S (0) = N (0) − 20 I (0) = 20 R (0) = 0

Dengan

menggunakan skala berikut dan R = N 0C.

S = N 0 A,

I = N 0 B,

Persamaan kondisi awal di atas akan menjadi 20 B(0) = 20 , C (0) = 0. Dengan A(0) = 1 − , N0 N0 mengambil N 0 = 2000 maka kondisi awalnya 20 adalah A(0) = 1 − = 0.99, 2000 20 B(0) = = 0.01, dan C (0) = 0. 2000 Ketika populasi yang terinfeksi 10% dari populasi maka populasi lainnya rentan terserang penyakit S (0) = N 0 − 200 . Kondisi awal lainnya R(0) = 0 , sehingga dipenuhi:

N (0) = S (0) + I (0) + R(0) . Jika dituliskan secara menyeluruh, kondisi awal dari model ini adalah S (0) = N (0) − 200 I (0) = 200

R(0) = 0 Dengan

S = N 0 A,

menggunakan skala berikut I = N 0 B, dan R = N 0C.

Persamaan kondisi awal di atas akan menjadi 200 200 , C (0) = 0. A(0) = 1 − , B(0) = N0 N0 Dengan mengambil N 0 = 2000 maka kondisi 200 awalnya adalah A(0) = 1 − = 0.9, 2000 B(0) =

200 = 0.1, C (0) = 0. 2000

Perbandingan dinamika populasi ketika populasi awal yang terinfeksi 1% dan 10% dari populasi untuk R0 < 1 Ketika yang terinfeksi 1% dari populasi

Ketika yang terinfeksi 10% dari populasi A

A 2.5

2.5

2.0

2.0

1.5

1.5

1.0

1.0

0.5

0.5

0.0

0

5

10

15

20

25

0.0

t

Gambar 53 Dinamika populasi A terhadap waktu t B

0

5

15

20

25

t


B

0.010

0.10

0.008

0.08

0.06

0.006

0.04

0.004

0.02

0.002

0.000

10

0

5

10

15

20

25

t


0.00

0

5

10

15

20

25

t


23

c

c

0.006

0.06

0.005

0.05

0.004

0.04

0.003

0.03

0.02

0.002

0.01

0.001

0.000

0

5

10

15

20

25

0.00

t


A,B,c

10

15

20

25

t

A,B,c 2.5

2.0

2.0

1.5

1.5

1.0

1.0

0.5

0.5

0

5


2.5

0.0

0

5

10

15

20

25

Gambar 59 Dinamika populasi A, B dan C terhadap waktu t

Pada gambar 53 dan 54 terlihat bahwa dengan bertambahnya populasi yang terinfeksi, kurva lebih curam yang berarti bahwa kurva semakin cepat menuju kestabilan tapi dengan tingkat kestabilan yang lebih kecil. Hal ini dikarenakan banyaknya populasi yang terinfeksi yang akan mengurangi jumlah populasi yang rentan karena total populasi konstan. Begitu pun pada gambar 55-58, pada populasi yang terinfeksi 10% kurva lebih curam yang berarti lebih cepat menuju kestabilan.

t

0.0

0

5

10

15

20

25

t


Pada gambar 59, populasi yang terinfeksi dan populasi yang sembuh hanya berupa garis horizontal, ini karena nilainya yang kecil. Sedangkan pada gambar 60, bisa dilihat bahwa kurva populasi yang terinfeksi pada awalnya berada di atas kurva populasi yang sembuh. Kemudian berpotongan, dan setelah itu kurva populasi yang sembuh berada di atas kurva populasi yang terinfeksi. Namun, pada kondisi kesetimbangan keduanya berada pada kondisi yang sama yaitu nol.

24

Perbandingan dinamika populasi ketika jumlah yang terinfeksi 1% dan 10% dari populasi ketika R0 > 1

Ketika yang terinfeksi 1% dari populasi

Ketika yang terinfeksi 10% dari populasi

A

A

2.5

2.5

2.0

2.0

1.5

1.5

1.0

1.0

0.5

0.5

0.0

0

20

40

60

80

100

t


0.0

40

60

80

100

t

B

1.4

1.4

1.2

1.2

1.0

1.0

0.8

0.8

0.6

0.6

0.4

0.4 0.2

0.2

0

20

40

60

80

100

0.0

t


20

40

60

80

100

t

c

1.5

1.5

1.0

1.0

0.5

0.5

0

0


c

0.0

20


B

0.0

0

20

40

60

80

100

t


0.0

0

20

40

60

80

100

t


25

A,B,c

A,B,c

2.5

2.5

2.0

2.0

1.5

1.5

1.0

1.0

0.5

0.5

0.0

0

20

40

60

80

100

t


Pada gambar 61, kurva meningkat pada awal pengamatan lalu menurun dan konstan pada 1.53. Begitupun pada gambar 62, kurva mengalami peningkatan, tapi peningkatan yang terjadi tidak sebesar pada populasi awal yang terinfeksi 1% dari populasi. Yang sama dari kedua kurva ini adalah keduanya konstan pada kondisi yang sama. Pada gambar 63 dan 64, selain kedua kurva berawal dari nilai yang berbeda, keduanya juga konstan pada nilai yang berbeda. Pada populasi awal yang terinfeksi 1% dari populasi, kurva konstan pada 1.28. Sedangkan pada populasi awal yang terinfeksi 10% dari populasi, kurva konstan pada 0.97. Pada gambar 65 dan 66, kedua kurva berawal dari nilai yang sama yaitu nol kemudian meningkat dan mencapai keadaan yang konstan. Tapi peningkatan yang terjadi tidak sama, dimana pada populasi awal yang terinfeksi 1% dari populasi peningkatan yang

0.0

0

20

40

60

80

100

t


terjadi lebih tinggi daripada pada populasi awal yang terinfeksi 10% dari populasi. Hal ini mengakibatkan kondisi konstan yang terjadi berbeda. Pada populasi awal yang terinfeksi 1% dari populasi, kurva konstan pada 1.6. Sedangkan pada populasi awal yang terinfeksi 10% dari populasi, kurva konstan pada 1.2. Pada gambar 67 terlihat bahwa kurva populasi yang sembuh lebih tinggi daripada kurva populasi yang rentan terserang penyakit. Sedangkan pada gambar 68, kurva populasi yang sembuh berada dibawah kurva populasi yang rentan terserang penyakit. Selain itu, pada gambar 67 ketiga kurva lebih berdekatan daripada pada gambar 68. Ini berarti bahwa dengan bertambahnya populasi awal yang terinfeksi akan menurunkan tingkat kestabilan populasi yang terinfeksi yang berakibat menurunnya populasi yang sembuh.

KESIMPULAN

Dalam karya tulis ini dianalisis model kestabilan SIRS pada proses penularan penyakit influenza. Analisis kestabilan penyebaran penyakit influenza pada suatu populasi menghasilkan dua titik tetap, yaitu titik tetap tanpa penyakit dan titik tetap endemik. Dari gambar bidang fase yang menunjukkan orbit kestabilan, dapat disimpulkan bahwa kedua titik tidak pernah stabil secara bersamaan. Ketika titik tetap pertama stabil titik tetap kedua tak stabil dan begitu sebaliknya. Titik tetap pertama berada dalam kestabilan ketika reproduksi dasar virus kurang dari 1 dan titik tetap kedua berada

dalam kestabilan ketika reproduksi dasar virus lebih besar dari 1. Dinamika populasi ketika reproduksi dasar virus kurang dari satu ( R0 < 1 ), saat laju kesembuhan (v) diperkecil populasi yang sembuh akan semakin lambat menuju kestabilan yang berakibat populasi pada kelas yang rentan dan terinfeksi pun lambat menuju kestabilan. Sedangkan ketika koefisien transmisi penularan penyakit (β) diperkecil populasi yang terinfeksi akan semakin cepat menuju kestabilan yang berakibat populasi yang rentan dan sembuh pun semakin cepat menuju kestabilan. Selain itu, titik balik

26

kesembuhannya lebih kecil daripada ketika β lebih besar. Dinamika populasi ketika reproduksi dasar virus lebih besar dari satu ( R0 > 1 ), saat laju kesembuhan (v) diperkecil populasi yang sembuh dan populasi yang terinfeksi semakin cepat menuju kestabilan meskipun dengan tingkat yang lebih besar. Begitupun ketika koefisien penularan penyakit (β) diperbesar, populasi yang terinfeksi semakin cepat menuju kestabilan tapi dengan tingkat yang lebih besar. Hal ini terjadi juga pada populasi yang sembuh. Sedangkan pada populasi yang rentan terserang penyakit, meskipun semakin cepat menuju kondisi tidak stabil tapi dengan jumlah populasi yang lebih kecil. Dari dinamika populasi juga dapat dismpulkan bahwa semakin kecil nilai R0 untuk R0 < 1 maka populasi semakin cepat

menuju kestabilan. Dan begitu sebaliknya, semakin besar nilai R0 untuk R0 > 1 maka populasi semakin cepat menuju ketidakstabilan yang berarti virus dapat bertahan dalam populasi. Semakin banyak populasi awal yang terinfeksi pada R0 < 1 akan menyebabkan populasi menjadi lebih cepat menuju kestabilan. Sedangkan pada R0 > 1 semakin banyak populasi awal yang terinfeksi menyebabkan populasi yang terinfeksi dan sembuh pada kondisi konstan menurun. Selain itu terlihat pula bahwa banyaknya populasi yang sembuh proporsional terhadap banyaknya populasi yang terinfeksi.

DAFTAR FUSTAKA Abel, M. L. dan J. P. Braselton. 1994. Mathematica by Example. Revised Ed. AP Professional, USA. Anderson R. M. 1982. The Population Dynamics of Infectious Diseases: Theory and Applications. Chapman and Hall, London. Anderson R. M. & May R. M. 1978. Population Biology of Infectious Diseases Part I, Nature, 280: 361-367. Cairns, A. 1995. Primary Component of Epidemiological Model. D. Mollison (Ed.), Epidemic models: Their Structure and Relation to Data, Newton Inst. Public, 350-371. Farlow, S. J. 1994. An Introduction to Differential Equation and Their Application. Mc. Graw-Hill, New York. Giesecko, J. 1994. Mathematical Models for Epidemics. Modern Infectious Disease Epidemiology. Edward Arnold, London, 109-123. Huntley, I. and R. M. Johnson. 1983. Linear and Nonlinear Differential Equations. Ellis Horword Limited, Chichester.

Keshet, L. E. 1987. Mathematical Models in Biology. Random House, New York. Kutha, A. N. K. 2004. Panduan Penggunaan Mathematica. Bogor. Mulchy and Pascho, 1984. Leon, S. J. 1998. Aljabar Linear dan Aplikasinya, Edisi ke-5. Terjemahan Alit Bondan. Erlangga. Ogut, H. 2001. Modeling of Fish Disease Dynamics: A New Approach to an Old Problem. Turkish Journal of Fisheries and Aquatic Sciences, Turkey, 67-74. Szidarovszky, F. & A. T. Bahill. 1998. Linear System Theory. CRC Press. Florida. Verhulst, F. 1990. Nonlinear Differential Equation and Dynamical System. Springer. Verlag. Heidelberg. Germany. Waltman, P. 1974. Lecture Notes in Biomathematics, Deterministic Threshold Models inTthe Theory of Epidemics. Springer-Verlag, New York.

27

LAMPIRAN

28

Lampiran 1 Penurunan persamaan (4), (5), (6) dari persamaan (1), (2), (3)

Diketahui persamaan (1), (2), (3) adalah

dS = S0 + (−bS − β SI + vR) dt dI = β SI − (α + b + v) I dt dR = vI − bR dt

(1) (2) (3)

Dan didefinisikan skala sebagai berikut:

S = N0 A I = N0 B R = N 0C Kemudian substitusikan skala diatas kedalam persamaan (1), (2), (3) sehingga diperoleh persamaan (4), (5), (6). •

dS = S0 + (−bS − β SI + vR) dt d ( N 0 A) = N 0 A0 + (−bN 0 A − β N 0 AN 0 B + vN 0C ) dt N dA ⇔ 0 = N 0 ( A0 + (−bA − β N 0 AB + vC )) dt dA ⇔ = A0 + (−bA − β N 0 AB + vC ) dt ∴ diperoleh persamaan (4) di bawah ini: dA = A0 + (−bA − β N 0 AB + vC ) dt ⇔

•

dI = β SI − (α + b + v) I dt d ( N 0 B) = β N 0 AN 0 B − (α + b + v) N 0 B dt N dB ⇔ 0 = N 0 ( β N 0 AB − (α + b + v) B) dt dB ⇔ = β N 0 AB − (α + b + v) B dt ∴ diperoleh persamaan (5) di bawah ini: dB = β N 0 AB − (α + b + v) B dt ⇔

•

dR = vI − bR dt d ( N0C ) = vN 0 B − bN 0C dt N dC ⇔ 0 = N 0 (vB − bC ) dt ⇔

29

dC = (vB − bC ) dt ∴ diperoleh persamaan (6) di bawah ini: dC = (vB − bC ) dt ⇔

Lampiran 2 Pencarian titik tetap dari persamaan (4), (5) dan (6)

Untuk menentukan titik tetap dari persamaan (4), (5) dan (6) maka persamaan tersebut dibuat sama dengan nol seperti dalam persamaan (7), (8) dan (9).

A0 + (−bA − β N 0 AB + vC) = 0 β N 0 AB − (α + b + v) B = 0

(8)

vB − bC = 0

(9)

(7)

Persamaan (8) dapat disederhanakan agar diperoleh nilai B.

β N0 AB − (α + b + v) B = 0 ⇔ B( β N 0 A − (α + b + v)) = 0 ⇔ B = 0 V β N 0 A − (α + b + v) = 0 Ketika B = 0 maka akan diperoleh titik tetap tanpa penyakit. Sekarang kita substitusikan B = 0 kedalam persamaan (9)

⇔ v.0 − bC = 0 ⇔ −bC = 0 ⇔b=0 V C =

0

Karena b ≠ 0 , maka dipilih C = 0 . Sekarang kita substitusikan nilai B dan C yang sudah didapat ke persamaan (7).

⇔ A0 + (−bA − β N0 A.0 + v.0) = 0

⇔ A0 − bA = 0 ⇔ A0 = bA

⇔ A=

A0 b

⎛A ⎞ ∴ diperoleh titik tetap tanpa penyakit yaitu T0 = ⎜ 0 , 0, 0 ⎟ . ⎝ b ⎠ Dari persamaan (8) tadi diperoleh β N 0 A − (α + b + v) = 0 atau A =

α +b+v . β N0

substitusikan nilai A kedalam persamaan (7) sehingga diperoleh titik tetap endemik. α +b+v α +b+v ⇔ A0 + (−b − β N0 B + vC ) = 0 β N0 β N0 ⇔ A0 + (−b ⇔ A0 − b

α +b+v − (α + b + v) B + vC ) = 0 β N0

α +b+v − (α + b + v) B + vC = 0 β N0

Kemudian

30

α +b+v = (α + b + v) B − vC β N0 Persamaan diatas eliminasikan dengan persamaan (9) α +b+v (α + b + v) B − vC = A0 − b β N0 ⇔ A0 − b

vB − bC = 0 ⇔ (α + b + v) B − vC = A0 − b

α +b+v β N0

⇔ vB − bC = 0

xv

⇔ b(α + b + v) B − bvC = bA0 − b 2 v 2 B − bvC = 0

⇔

xb

α +b+v β N0 -

α +b+v β N0 α +b+v ⇔ B (b(α + b + v) − v 2 ) = bA0 − b 2 β N0 α +b+v ) b( A0 − b β N0 ⇔B= b(α + b + v) − v 2 Selanjutnya substitusikan nilai B ke persamaan (9) untuk memperoleh nilai C. ⇒ b(α + b + v) B − v 2 B = bA0 − b 2

vB − bC = 0 b( A0 − b ⇔ v(

α +b+v ) β N0

b(α + b + v ) − v 2 b( A0 − b

⇔ v(

) − bC = 0

α +b+v ) β N0

) = bC b(α + b + v) − v 2 α +b+v b( A0 − b ) β N0 v( ) b(α + b + v) − v 2 ⇔C = b v( A0 − b(α + b + v )) ⇔C= β N 0 (b(α + b + v) − v 2 )

∴ diperoleh titik tetap endemik yaitu T1 = ( A, B, C ) , dengan:

A=

α +b+v β N0 b( A0 − b

B=

C=

α +b+v ) β N0

b(α + b + v) − v 2

v ( A0 − b(α + b + v)) β N 0 (b(α + b + v ) − v 2 )

31

Lampiran 3 Pencarian titik tetap dari persamaan (4), (5), (6) dengan menggunakan Mathematica

Program Mathematica untuk menentukan titik tetap dari persamaan (4), (5), (6) adalah Solve[{A0-b*A+v*C-β*B*A*N0 0,β*B*A*N0-(α+b+v)*B 0,v*B-b*C 0},{A,B,C}]

Diperoleh hasil sebagai berikut untuk titik tetap tanpa penyakit

A0 ⎛ ⎞ ⎜ A → , B → 0, C → 0 ⎟ b ⎝ ⎠ Dan titik tetap endemik adalah

α +b+v ⎛ ⎞ b( A0 − b ) ⎜ v( A0 − b(α + b + v)) ⎟ β N0 α +b+v ⎜A→ ⎟ ,B → , C → b(α + b + v) − v 2 β N0 β N 0 (b(α + b + v) − v 2 ) ⎟ ⎜ ⎜ ⎟ ⎝ ⎠ Lampiran 4 Penentuan nilai eigen titik tetap tanpa penyakit

Diketahui matriks Jacobi adalah

⎡ dx ⎢ dA ⎢ dy J =⎢ ⎢ dA ⎢ ⎢ dz ⎢⎣ dA

dx ⎤ dC ⎥ ⎡ −b − β N B v⎤ −β N0 A ⎥ 0 dy ⎥ ⎢ = β N0 B β N 0 A − (α + b + v) 0 ⎥⎥ dC ⎥ ⎢ v b ⎥⎦ 0 ⎥ ⎢ dz ⎥ ⎣ dC ⎥⎦ A ⎛ ⎞ Pelinearan pada titik tetap tanpa penyakit ⎜ A → 0 , B → 0, C → 0 ⎟ , diperoleh matriks J0. b ⎝ ⎠ − β N 0 A0 ⎡ ⎤ v⎥ ⎢ −b b ⎢ ⎥ N A β ⎢ 0 0 − (α + b + v) 0 ⎥ J0 = 0 ⎢ ⎥ b ⎢ ⎥ −b ⎥ v ⎢0 ⎣⎢ ⎦⎥ dx dB dy dB dz dB

Kemudian dicari nilai eigennya dengan menggunakan persamaan karakteristik det(A-λI)=0.

−b − λ 0 0

β N 0 A0 b

− β N 0 A0 b

v

− (α + b + v) − λ

0

v

−b − λ

=0

32

− β N 0 A0 − (α + b + v)) − λ ) = 0 b − β N 0 A0 ⇒ (−b − λ )2 = 0 atau ( − (α + b + v)) − λ = 0 b β N 0 A0 − (α + b + v) Jadi nilai eigennya adalah λ1 = λ2 = −b dan λ3 = b ⇒ (−b − λ )2 ((

Lampiran 5 Penentuan nilai eigen titik tetap tanpa penyakit dengan Mathematica 6.0

Nilai eigen dari titik tetap tanpa penyakit dengan menggunakan Mathematica 6.0 : Eigenvalues[{{-b,(-β*N0*A(0))/b,v},{0,(β*N0*A(0))/b-(α+b+v),0},{0,v,-b}}]

Sehingga diperoleh nilai eigen yaitu:

β N 0 A0 ⎧ ⎫ − (α + b + v) ⎬ ⎨−b, −b, b ⎩ ⎭ Lampiran 6 Penentuan nilai eigen titik tetap endemic

Diketahui matriks Jacobi untuk titik tetap endemik

α +b+v ⎛ ⎞ b( A0 − b ) ⎜ β N0 v( A0 − b(α + b + v)) ⎟ α +b+v ⎜A= ⎟ adalah J1. = ,B = , C β N0 β N 0 (b(α + b + v) − v 2 ) ⎟ b(α + b + v) − v 2 ⎜ ⎜ ⎟ ⎝ ⎠ ⎡ ⎤ β bA0 N 0 − b 2 (α + b + v) −(α + b + v) v ⎥ ⎢ −b − 2 2 b + bv + bα − v ⎢ ⎥ ⎢ β bA0 N 0 − b 2 (α + b + v) ⎥ J1 = ⎢ 0 0⎥ 2 2 b + bv + bα − v ⎢ ⎥ v −b ⎥ 0 ⎢ ⎢ ⎥ ⎣ ⎦ Kemudian dicari nilai eigennya dengan menggunakan persamaan karakteristik det(A-λI)=0. β bA0 N 0 − b 2 (α + b + v) v −b − − λ −(α + b + v) b 2 + bv + bα − v 2 β bA0 N 0 − b 2 (α + b + v) 0−λ 0 =0 b 2 + bv + bα − v 2 v 0 −b − λ 2 β bN 0 A0 − b 2 (α + b + v) 2 β bN 0 A0 − b (α + b + v ) ) − λ )( − λ )( − b − λ ) + v b 2 + bv + bα − v 2 b 2 + bv + bα − v 2 β bN 0 A0 − b 2 (α + b + v) +( )(α + b + v)(−b − λ ) = 0 b 2 + bv + bα − v 2

⇒ ((−b −

Perhitungan ini sulit dilakukan secara manual sehingga dilakukan perhitungan dengan menggunakan Mathematica dengan program sebagai berikut:

33

Solve[((-b-(β*b*N0*A 0-b2 (α+b+v))/(b2+b*v+b*α-v2))-λ) (λ) (b+λ)+v2 ((β*b*N0*A 0-b2 (α+b+v))/(b2+b*v+b*α-v2))+((β*b*N0*A 0-b2 (α+b+v))/(b2+b*v+b*αv2)) (α+b+v) (-b-λ) 0,λ]

Sehingga diperoleh nilai eigen dibawah ini.

λ1 → −

−b3 − b 2 v + 2bv 2 − b 2α − β bA0 N 0 − (21/3 (−(−b3 − b 2 v + 2bv 2 − b 2α − β bA0 N 0 ) 2 + 3(−b 2 − bv + v 2 − bα ) 3(−b 2 − bv + v 2 − bα )

(b 4 + 2b3v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 ))) (3(−b 2 − bv + v 2 − bα )(−16b9 − 66b8 v − 51b 7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v5 − 34b3v 6 + 27b 2 v 7 − 66b8α − 204 b 7 v α − 54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3v5α + 27b 2 v 6α −102b7α 2 − 210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v3α 2 − 63b3v 4α 2 − 70b 6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b 6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 −3b3 v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3 v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A20 N 02 − 21b 4 vβ 2 A20 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v 3 β 2 A20 N 02 − 21b 4αβ 2 A20 N 02 − 18b3vαβ 2 A2 0 N 02 +9b 2 v 2αβ 2 A2 0 N 02 − 9b3α 2 β 2 A20 N 02 + 2b3 β 3 A30 N 03 + ((−16b9 − 66b8 v − 51b7 v 2 + 86b 6 v 3 + 75b5 v 4 − 75b 4 v5 − 34b3v 6 + 27b 2 v 7 − 66b8α − 204 b7 v α − 54b6 v 2α + 246b5 v 3α + 33b 4 v 4α − 126b3v5α + 27b 2 v 6α − 102b 7α 2 − 210b 6 vα 2 + 51b5 v 2α 2 + 162b 4 v3α 2 − 63b3v 4α 2 − 70b 6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b6 v β A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 − 3b3 v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A20 N 02 − 21b 4 vβ 2 A20 N 02 − 3b3v 2 β 2 A2 N 02 +9b 2 v3 β 2 A20 N 02 − 21b 4αβ 2 A20 N 02 − 18b3 vαβ 2 A20 N 02 + 9b 2 v 2αβ 2 A2 0 N 02 − 9b3α 2 β 2 A20 N 02 + 2b3 β 3 A30 N 03 ) 2 + 4(−(−b3 − b 2 v + 2bv 2 − b 2α − bβ A0 N 0 ) 2 + 3(−b 2 − bv + v 2 − bα )(b 4 + 2b3v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 ))3 ))1/3 ) + (−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v5 −34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b 7 v α − 54b 6 v 2α + 246b5 v 3α + 33b 4 v 4α − 126b3 v 5α + 27b 2 v 6α − 102b 7α 2 − 210b 6 vα 2 + 51b5 v 2α 2 + 162b 4 v 3α 2 − 63b3 v 4α 2 − 70b 6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b 7 β A0 N 0 + 75b 6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v3 β A0 N 0 − 3b3v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b 6αβ A0 N 0

+156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 v β 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3 vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 + ((−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b7 v α − 54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3 v5α + 27b 2 v 6α − 102b7α 2 − 210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v3α 2 − 63b3v 4α 2 − 70b 6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 − 3b3v 4 β A0 N 0 + 63b 2 v 5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0

34

+81b 4 vα 2 β A0 N 0 − 63b3 v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A30 N 03 ) 2 + 4(−(−b3 − b 2 v + 2bv 2 − b 2α − bβ A0 N 0 ) 2 + 3(−b 2 − bv + v 2 − bα )(b 4 + 2b3 v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 ))3 ))1/3 (3(21/ 3 )(−b 2 − bv + v 2 − bα ))

λ2 → −

−b3 − b 2 v + 2bv 2 − b 2α − β bA0 N 0 + ((1 + i 3)(−(−b3 − b 2 v + 2bv 2 − b 2α − β bA0 N 0 ) 2 + 3(−b 2 − bv + v 2 − bα ) 3(−b 2 − bv + v 2 − bα )

(b 4 + 2b3v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 ))) (3(22 / 3 )(−b 2 − bv + v 2 − bα )(−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v5 − 34b3v 6 + 27b 2 v 7 − 66b8α − 204 b 7 v α − 54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3v 5α + 27b 2 v 6α − 102b7α 2 − 210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v 3α 2 − 63b3 v 4α 2 − 70b6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b 7 β A0 N 0 + 75b 6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 − 3b3 v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b 6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3 v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 + ((−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v 5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b7 vα −54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3v5α + 27b 2 v 6α − 102b 7α 2 − 210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v3α 2 − 63b3v 4α 2 − 70b6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 +75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v3 β A0 N 0 − 3b3v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3v 2 β 2 A2 N 02 +9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 )2 + 4(−(−b3 − b 2 v + 2bv 2 − b 2α − bβ A0 N 0 )2 + 3( −b 2 − bv + v 2 − bα )(b 4 + 2b3 v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bv β A0 N 0 − bαβ A0 N 0 ))3 ))1/3 ) − ((1 − i 3)(−16b9 − 66b8 v − 51b 7 v 2 + 86b 6 v 3 + 75b5 v 4 − 75b 4 v 5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b 7 v α − 54b 6 v 2α + 246b5 v 3α + 33b 4 v 4α − 126b3 v 5α + 27b 2 v 6α − 102b 7α 2 − 210b 6 vα 2 + 51b5 v 2α 2 + 162b 4 v 3α 2 − 63b3 v 4α 2 − 70b 6α 3 −72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 − 3b3 v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0

+156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 v β 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3 vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 + ((−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b7 v α − 54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3 v5α + 27b 2 v 6α − 102b7α 2 − 210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v3α 2 − 63b3v 4α 2 − 70b 6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 −

35

3b3v 4 β A0 N 0 + 63b 2 v 5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 +81b 4 vα 2 β A0 N 0 − 63b3 v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 ) 2 + 4(−(−b3 − b 2 v + 2bv 2 − b 2α − bβ A0 N 0 ) 2 + 3(−b 2 − bv + v 2 − bα )(b 4 + 2b3v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 ))3 ))1/ 3 ) (6(21/ 3 )(−b 2 − bv + v 2 − bα ))

λ3 → −

−b3 − b 2 v + 2bv 2 − b 2α − β bA0 N 0 + ((1 − i 3)(−(−b3 − b 2 v + 2bv 2 − b 2α − β bA0 N 0 ) 2 + 3(−b 2 − bv + v 2 − bα ) 2 2 3(−b − bv + v − bα )

(b 4 + 2b3v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 ))) (3(22 / 3 )(−b 2 − bv + v 2 − bα )(−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v5 − 34b3v 6 + 27b 2 v 7 − 66b8α − 204 b 7 v α − 54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3v 5α + 27b 2 v 6α − 102b7α 2 − 210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v 3α 2 − 63b3 v 4α 2 − 70b6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b 7 β A0 N 0 + 75b 6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 − 3b3 v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b 6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3 v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 + ((−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v 5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b7 vα −54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3v5α + 27b 2 v 6α − 102b 7α 2 − 210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v3α 2 − 63b3v 4α 2 − 70b6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 +75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v3 β A0 N 0 − 3b3v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3v 2 β 2 A2 N 02 +9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 )2 + 4(−(−b3 − b 2 v + 2bv 2 − b 2α − bβ A0 N 0 )2 + 3( −b 2 − bv + v 2 − bα )(b 4 + 2b3 v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bv β A0 N 0 − bαβ A0 N 0 ))3 ))1/3 ) − ((1 + i 3)(−16b9 − 66b8 v − 51b 7 v 2 + 86b 6 v 3 + 75b5 v 4 − 75b 4 v 5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b 7 v α − 54b 6 v 2α + 246b5 v 3α + 33b 4 v 4α − 126b3 v 5α + 27b 2 v 6α − 102b 7α 2 − 210b 6 vα 2 + 51b5 v 2α 2 + 162b 4 v 3α 2 − 63b3 v 4α 2 − 70b 6α 3 −72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 − 3b3 v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0

+156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 v β 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3 vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 + ((−16b9 − 66b8 v − 51b7 v 2 + 86b6 v3 + 75b5 v 4 − 75b 4 v5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b7 v α − 54b6 v 2α + 246b5 v3α + 33b 4 v 4α − 126b3 v5α + 27b 2 v 6α − 102b7α 2 −

36

210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v3α 2 − 63b3v 4α 2 − 70b 6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v 3 β A0 N 0 − 3b3v 4 β A0 N 0 + 63b 2 v 5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 +81b 4 vα 2 β A0 N 0 − 63b3 v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3vαβ 2 A02 N 02 + 9b 2 v 2αβ 2 A02 N 02 − 9b3α 2 β 2 A02 N 02 + 2b3 β 3 A03 N 03 ) 2 + 4(−(−b3 − b 2 v + 2bv 2 − b 2α − bβ A0 N 0 ) 2 + 3(−b 2 − bv + v 2 − bα )(b 4 + 2b3v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 ))3 ))1/ 3 ) (6(21/ 3 )(−b 2 − bv + v 2 − bα ))

Nilai diatas dapat disederhanakan menjadi λ1 → −

3 [d + d 2 + 4(− a 2 + 3bc)3 ]1/ 3 a 2(− a 2 + 3bc) − + 3b 3b[d + d 2 + 4(−a 2 + 3bc)3 ]1/ 3 3 3 2b

λ2 → −

(1 − i 3)[d + d 2 + 4(−a 2 + 3bc)3 ]1/ 3 a (1 + i 3)(− a 2 + 3bc) + − 3b 3 3 4b[d + d 2 + 4(−a 2 + 3bc)3 ]1/ 3 6 3 2b

λ3 → −

(1 + i 3)[d + d 2 + 4(− a 2 + 3bc)3 ]1/ 3 a (1 − i 3)(− a 2 + 3bc) + − 3b 3 3 4b[d + d 2 + 4(−a 2 + 3bc)3 ]1/ 3 6 3 2b

Dengan nilai a, b, c, dan d sebagai berikut: a = −b3 − b 2 v + 2bv 2 − b 2α − β bA0 N 0 b = −b 2 − bv + v 2 − bα c = b 4 + 2b3 v + 2b 2 v 2 + 2b3α + 2b 2 vα + b 2α 2 − 2b 2 β A0 N 0 − bvβ A0 N 0 − bαβ A0 N 0 d = −16b9 − 66b8 v − 51b 7 v 2 + 86b 6 v 3 + 75b5 v 4 − 75b 4 v 5 − 34b3 v 6 + 27b 2 v 7 − 66b8α − 204 b 7 v α − 54b 6 v 2α + 246b5 v 3α + 33b 4 v 4α − 126b3 v 5α + 27b 2 v 6α − 102b 7α 2 −210b6 vα 2 + 51b5 v 2α 2 + 162b 4 v 3α 2 − 63b3 v 4α 2 − 70b6α 3 − 72b5 vα 3 + 54b 4 v 2α 3 − 18b5α 4 + 24b7 β A0 N 0 + 75b6 vβ A0 N 0 + 27b5 v 2 β A0 N 0 − 87b 4 v3 β A0 N 0

−3b3 v 4 β A0 N 0 + 63b 2 v5 β A0 N 0 − 27bv 6 β A0 N 0 + 75b6αβ A0 N 0 + 156b5 vαβ A0 N 0 − 33b 4 v 2αβ A0 N 0 − 126b3v3αβ A0 N 0 + 63b 2 v 4αβ A0 N 0 + 78b5α 2 β A0 N 0 + 81b 4 vα 2 β A0 N 0 − 63b3v 2α 2 β A0 N 0 + 27b 4α 3 β A0 N 0 − 12b5 β 2 A02 N 02 − 21b 4 vβ 2 A02 N 02 − 3b3 v 2 β 2 A2 N 02 + 9b 2 v 3 β 2 A02 N 02 − 21b 4αβ 2 A02 N 02 − 18b3 vαβ 2 A03 N 02 + 9b 2 v 2αβ 2 A03 N 02 − 9b3α 2 β 2 A03 N 02 + 2b3 β 3 A03 N 03

37

Lampiran 7 Bukti

Titik tetap tanpa penyakit dari model ini bersifat stabil asimtotik lokal jika dan hanya jika R0 < 1 dan tidak stabil jika dan hanya jika R0 > 1. Karena λ1 = λ2 < 0, maka titik tetap tanpa penyakit akan bersifat asimtotik lokal jika λ3 < 0 dan tak stabil jika λ3 > 0. Akan dibuktikan bahwa nilai λ3 < 0 jika dan hanya jika R0 < 1 dan λ3 > 0 jika dan hanya jika R0 > 1. • Nilai λ3 < 0 jika dan hanya jika R0 < 1 Bukti arah ke kanan.

⇒ λ3 < 0 β N 0 A0 ⇔ − (α + b + v) < 0 b β N 0 A0 ⇔ −1 < 0 b(α + b + v) β N 0 A0 ⇔ <1 b(α + b + v) ⇔ R0 < 1

Bukti arah ke kiri.

⇒ R0 < 1 β N 0 A0 ⇔ <1 b(α + b + v) β N 0 A0 ⇔ −1 < 0 b(α + b + v) β N 0 A0 − (α + b + v) < 0 b ⇔ λ3 < 0 • Nilai λ3 > 0 jika dan hanya jika R0 > 1 Bukti arah ke kanan.

⇒ λ3 > 0 β N 0 A0 ⇔ − (α + b + v) > 0 b β N 0 A0 ⇔ −1 > 0 b(α + b + v) β N 0 A0 ⇔ >1 b(α + b + v) ⇔ R0 > 1

Bukti arah ke kiri.

⇒ R0 > 1 β N 0 A0 ⇔ >1 b(α + b + v) β N 0 A0 ⇔ −1 > 0 b(α + b + v)

38

⇔

β N 0 A0

b ⇔ λ3 > 0

− (α + b + v) < 0

Terbukti bahwa λ3 < 0 jika dan hanya jika R0 < 1 dan λ3 > 0 jika dan hanya jika R0 > 1. Lampiran 8

Gambar bidang fase ketika

β = 0.0002 dan v = 0.9 dengan menggunakan Mathematica 5.2

Untuk memperoleh bidang fase digunakan Dynpac dari DynSys10. 71 pada Mathematica 5.2 dengan program dibawah ini. sysid Mathematica 5.2 . 0 , DynPac intreset True plotreset 340 setstate[{A,B,c}]

10.71 ,

3 / 17 / 2008

setparm[{α,β,b,v,p}] {α,β,b,v,p}

slopevec= 9J :

1999 N − b∗A−β∗ p∗A∗ B+ v∗c, β∗ p∗ A∗ B−Hα+ b+ vL∗ B,v∗ B− b∗c= 2000

1999 −Ab+ cv− ABpβ, −BHb+ v+αL +ABpβ, −bc+ Bv> 2000

parmval={0.02,0.0002,0.4,0.9,2000} {0.02,0.0002,0.4,0.9,2000} sysname="fase1" fase1

initvec= 9 :

1999 1 ,0= , 2000 2000

1999 1 , ,0> 2000 2000

t0=0.0 0. window={{0,3},{0,0.5},{0,0.5}} {{0,3},{0,0.5},{0,0.5}} s=sugtimestep[window] 0.0225494 nsteps=2000 2000 sol=integrate[initvec,t0,s,nsteps]; Gambar hasil simulasi 3D graf1=phaseplot3D[sol,1,2,3]

se1 B8α, β, b, v, p<=8 0.02 , 0.0004 20 1

0.00 ,

0.40 ,

0.90 , 2000.00

1.5 A

2 1 2.50.00c02

Graphics3D boxrat={1,1,1} {1,1,1} graf2=phaseplot3D[sol,1,2,3] Mentransformasikan gambar hasil simulasi 3D ke dalam 2D

6D2dalamgambarhasilkeMentransformasikansimulasi

39

asprat=1.0 1. plrange={{0,3.5},{0,0.00055}} {{0,3.5},{0,0.00055}} arrowflag=True True

1 arrowvec= 9 = 3 : >

1 3

grafxyproj1=phaseplot[sol,1,2] plrange={{0,3.5},{0,0.0003}} {{0,3.5},{0,0.0003}} grafxzproj2=phaseplot[sol,1,3] plrange={{0,0.00055},{0,0.00035}} {{0,0.00055},{0,0.00035}} grafyzproj3=phaseplot[sol,2,3] Lampiran 9



Untuk memperoleh bidang fase digunakan Dynpac dari DynSys10. 71 pada Mathematica 5.2 dengan program seperti dibawah ini. sysid Mathematica 5.2 . 0 , DynPac intreset True plotreset 340 setstate[{A,B,c}]

10.71 ,

3 / 17 / 2008


slopevec= 9J :




initvec= 9 :

1999 1 ,0= , 2000 2000

1999 1 , ,0> 2000 2000


40 se2 B8α, β, b, v, p<=8 0.02 , 0.0004 20 1

0.00 ,

0.40 ,

0.80 , 2000.00

1.5 A

2 1 2.50.00c02


6D2dalamgambarhasilkeMentransformasikansimulasi asprat=1.0 1. plrange={{0,3.5},{0,0.00055}} {{0,3.5},{0,0.00055}} arrowflag=True True

1 arrowvec= 9 = 3 : >

1 3


Gambar bidang fase ketika β = 0.0001 dan v = 0.9 dengan menggunakan Mathematica 5.2


10.71 ,

3 / 17 / 2008


slopevec= 9J :




initvec= 9

1999 1 ,0= , 2000 2000

41

:

1999 1 , ,0> 2000 2000


se1 B8α, β, b, v, p<=8 0.02 , 0.0004 20 1

0.00 ,

0.40 ,

0.90 , 2000.00

1.5 A

2 05 1 2.50.00c02



1 arrowvec= 9 = 3 : >

1 3

grafxyproj1=phaseplot[sol,1,2] plrange={{0,3.5},{0,0.0003}} {{0,3.5},{0,0.0003}} grafxzproj2=phaseplot[sol,1,3] plrange={{0,0.00055},{0,0.00035}} {{0,0.00055},{0,0.00035}} grafyzproj3=phaseplot[sol,2,3]

42

Lampiran 11 Gambar bidang fase ketika β = 0.0002 dan v = 0.5 dengan menggunakan Mathematica 5.2


10.71 ,

3 / 19 / 2008


slopevec= 9J :




initvec= 9 :

1999 1 ,0= , 2000 2000

1999 1 , ,0> 2000 2000

t0=0.0 0. window={{0,1},{0,1},{0,1}} {{0,1},{0,1},{0,1}} s=sugtimestep[window] 0.0238366 nsteps=2000 2000 sol=integrate[initvec,t0,s,nsteps]; Gambar hasil simulasi 3D graf1=phaseplot3D[sol,1,2,3]

se1 B8α, β, b, v, p<=8 0.02 , 0.004 3 1 2 1

0.00 ,

0.40 ,

0.50 , 2000.00

1.5 A

2 0.004 2.5 c2



43

1 arrowvec= 9 = 3 : >

1 3





10.71 ,

3 / 17 / 2008


slopevec= 9J :




initvec= 9 :

1999 1 ,0= , 2000 2000

1999 1 , ,0> 2000 2000

t0=0.0 0. window={{0,1},{0,1},{0,1}} {{0,1},{0,1},{0,1}} s=sugtimestep[window] 0.0261755 nsteps=2000 2000 sol=integrate[initvec,t0,s,nsteps]; Gambar hasil simulasi 3D graf1=phaseplot3D[sol,1,2,3]

44 se1 8α, β, b, v, p<=8 0.02 ,

0.00 ,

0.40 ,

0.40 , 2000.00

B 0.2 0.1 0

1 0.2 0.1 c 0

1.5 A

2 2.5



1 arrowvec= 9 = 3 : >

1 3





10.71 ,

3 / 18 / 2008


slopevec= 9J :



parmval={0.02,0.0003,0.4,0.4,2000} {0.02,0.0003,0.4,0.4,2000}

45

sysname="fase1" fase1

initvec= 9 :

1999 1 ,0= , 2000 2000

1999 1 , ,0> 2000 2000

t0=0.0 0. window={{0,1},{0,1},{0,1}} {{0,1},{0,1},{0,1}} s=sugtimestep[window] 0.0248331 nsteps=2000 2000 sol=integrate[initvec,t0,s,nsteps]; Gambar hasil simulasi 3D graf1=phaseplot3D[sol,1,2,3] 8 α , β, b, v, p< = 8

se1

0.02 ,

0.00 ,

0.40 ,

0.40 , 2000.00

1 0.75 c 0.5 0.25 1

0

0.75

1

0.5 B

1.5 A

0.25

2 0



1 arrowvec= 9 = 3 : >

1 3

grafxyproj1=phaseplot[sol,1,2] plrange={{0,2.6},{0,1.2}} {{0,2.6},{0,1.2}} grafxzproj2=phaseplot[sol,1,3] plrange={{0,1.2},{0,1.2}} {{0,1.2},{0,1.2}} grafyzproj3=phaseplot[sol,2,3] Lampiran 14 Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.9

46

Untuk menggambar dinamika populasi digunakan Mathematica 6.0, dengan program seperti dibawah ini. α=0.02; b=0.4; v=0.9; β=0.0002; a=1999/2000; n=2000; sol=NDSolve[{A'[t] 1999/2000-b*A[t]-β*n*A[t]*B[t]+v*c[t], B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1999/2000, B[0] 1/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} ran1=Table[p,{t,0,1000}]; ran2=Table[q,{t,0,1000}]; ran3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,2.6}},ImageSize→300] gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.0006}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.0004}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,PlotRange→{{0,25},{0,2.6}},ImageSize→300] Lampiran 15 Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.8

Untuk menggambar dinamika populasi digunakan Mathematica 6.0, dengan program seperti dibawah ini. α=0.02; b=0.4; v=0.8; β=0.0002; a=1999/2000; n=2000; sol=NDSolve[{A'[t] 1999/2000-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]b*c[t],A[0] 1999/2000,B[0] 1/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} ran1=Table[p,{t,0,1000}]; ran2=Table[q,{t,0,1000}]; ran3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,2.6}},ImageSize→300]

47

gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.0006}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.0004}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,C"},DisplayFunction→$DisplayFuncti on,PlotRange→{{0,25},{0,2.6}},ImageSize→300] Lampiran 16 Program Mathematica untuk dinamika populasi ketika β = 0.0001 dan v = 0.9

Untuk menggambar dinamika populasi digunakan Mathematica 6.0, dengan program seperti dibawah ini. α=0.02; b=0.4; v=0.9; β=0.0001; a=1999/2000; n=2000; sol=NDSolve[{A'[t] 1999/2000-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1999/2000,B[0] 1/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} ran1=Table[p,{t,0,1000}]; ran2=Table[q,{t,0,1000}]; ran3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,2.6}},ImageSize→300] gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.0006}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.0004}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,PlotRange→{{0,25},{0,2.6}},ImageSize→300] Lampiran 17 Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.5

Untuk menggambar dinamika populasi digunakan Mathematica 6.0, dengan program seperti dibawah ini. α=0.02; b=0.4; v=0.5; β=0.0002; a=1999/2000; n=2000; sol=NDSolve[{A'[t] 1999/2000-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1999/2000,B[0] 1/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t],

48

B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} ran1=Table[p,{t,0,1000}]; ran2=Table[q,{t,0,1000}]; ran3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,200},{0,2.6}},ImageSize→300] gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange->{{0,200},{0,0.35}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,200},{0,0.35}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,ImageSize→300] Lampiran 18 Program Mathematica untuk dinamika populasi ketika β = 0.0002 dan v = 0.4

Untuk menggambar dinamika populasi digunakan Mathematica 6.0, dengan program seperti dibawah ini. α=0.02; b=0.4; v=0.4; β=0.0002; a=1999/2000; n=2000; sol=NDSolve[{A'[t] 1999/2000-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1999/2000,B[0] 1/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} ran1=Table[p,{t,0,1000}]; ran2=Table[q,{t,0,1000}]; ran3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,200},{0,2.6}},ImageSize→300] gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange->{{0,200},{0,0.5}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,200},{0,0.5}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,ImageSize→300] Lampiran 19 Program Mathematica untuk dinamika populasi ketika β = 0.0003 dan v = 0.4

Untuk menggambar dinamika populasi digunakan Mathematica 6.0, dengan program seperti dibawah ini. α=0.02;

49

b=0.4; v=0.4; β=0.0003; a=1999/2000; n=2000; sol=NDSolve[{A'[t] 1999/2000-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1999/2000,B[0] 1/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} ran1=Table[p,{t,0,1000}]; ran2=Table[q,{t,0,1000}]; ran3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,200},{0,2.6}},ImageSize→300] gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange->{{0,200},{0,1.5}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,200},{0,1.5}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,ImageSize→300] Lampiran 20 Program Mathematica ketika R0 < 1 dan populasi yang terinfeksi 1% dari total populasi α=0.02; β=0.0002; b=0.4; v=0.9; a=1980/2000; n=2000; sol=NDSolve[{A'[t] (1980/2000)-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1980/2000,B[0] 20/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} run1=Table[p,{t,0,1000}]; run2=Table[q,{t,0,1000}]; run3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,2.6}},ImageSize→300] gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.011}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.006}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,PlotRange→{{0,25},{0,2.6}},ImageSize→300]

50

Lampiran 21 Program Mathematica ketika R0 < 1 dan populasi yang terinfeksi 10% dari total populasi α=0.02; β=0.0002; b=0.4; v=0.9; a=1800/2000; n=2000; sol=NDSolve[{A'[t] (1800/2000)-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1800/2000,B[0] 200/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} ran1=Table[p,{t,0,1000}]; ran2=Table[q,{t,0,1000}]; ran3=Table[s,{t,0,1000}]; gambr1=ListPlot[ran1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,2.6}},ImageSize→300] gambr2=ListPlot[ran2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.11}},ImageSize→300] gambr3=ListPlot[ran3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,25},{0,0.06}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,PlotRange→{{0,25},{0,2.6}},ImageSize→300] Lampiran 22 Program Mathematica ketika R0 > 1 dan populasi yang terinfeksi 1% dari total populasi α=0.02; β=0.0003; b=0.4; v=0.5; a=1980/2000; n=2000; sol=NDSolve[{A'[t] (1-20/2000)-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1-20/2000,B[0] 20/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} run1=Table[p,{t,0,1000}]; run2=Table[q,{t,0,1000}]; run3=Table[s,{t,0,1000}]; gambr1=ListPlot[run1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,100},{0,2.5}},ImageSize→300]

51

gambr2=ListPlot[run2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange→{{0,100},{0,1.5}},ImageSize→300] gambr3=ListPlot[run3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,100},{0,1.7}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,PlotRange→{{0,100},{0,2.6}},ImageSize→300] Lampiran 23 Program Mathematica ketika R0 > 1 dan populasi yang terinfeksi 10% dari total populasi α=0.02; β=0.0003; b=0.4; v=0.5; a=1800/2000; n=2000; sol=NDSolve[{A'[t] (1800/2000)-b*A[t]-β*n*A[t]*B[t]+v*c[t],B'[t] β*n*A[t]*B[t]-(α+b+v) B[t],c'[t] v*B[t]-b*c[t],A[0] 1800/2000,B[0] 200/2000,c[0] 0},{A[t],B[t],c[t]},{t,0,1000}, DependentVariables→{A,B,c},MaxSteps→∞] {{A[t]→InterpolatingFunction[{{0.,1000.}},<>][t], B[t]→InterpolatingFunction[{{0.,1000.}},<>][t], c[t]→InterpolatingFunction[{{0.,1000.}},<>][t]}} {p,q,s}={A[t],B[t],c[t]}/.Flatten[sol] {InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t], InterpolatingFunction[{{0.,1000.}},<>][t]} run1=Table[p,{t,0,1000}]; run2=Table[q,{t,0,1000}]; run3=Table[s,{t,0,1000}]; gambr1=ListPlot[run1,PlotStyle→{RGBColor[1,0,0],Thickness[0.01]},AxesLabel→{"t","A"}, Joined→True,AspectRatio→1,PlotRange→{{0,100},{0,2.5}},ImageSize→300] gambr2=ListPlot[run2,PlotStyle→{RGBColor[0,1,0],Thickness[0.01]},AxesLabel→{"t","B"}, Joined→True,AspectRatio→1,PlotRange→{{0,100},{0,1.5}},ImageSize→300] gambr3=ListPlot[run3,PlotStyle→{RGBColor[0,0,1],Thickness[0.01]},AxesLabel→{"t","c"}, Joined→True,AspectRatio→1,PlotRange→{{0,100},{0,1.7}},ImageSize→300] Show[gambr1,gambr2,gambr3,AxesLabel→{"t","A,B,c"},DisplayFunction→$DisplayFunctio n,PlotRange→{{0,100},{0,2.6}},ImageSize→300]

MODEL SIRS PADA PROSES PENULARAN PENYAKIT INFLUENZA DENGAN POPULASI YANG TERINFEKSI VIRUS LIA MULYANAH

Recommend Documents