Jurnal Matematika, Sains, dan Pembelajarannya
Vol.8, No.2, Juli 2022, pp. 197~203
ISSN 2407-4098 (print)
DOI:10.31605/saintifik.v8i2.381
ISSN 2622-8904 (online)

## Analisis Model Dua Mangsa Satu Pemangsa dengan Pertahanan Kelompok dan Pemanenan Linear pada Mangsa
Fardinah* 1 , Darma Ekawati 2 , Hikmah 3 , Ahmad Ansar 4 1,2,3,4 Universitas Sulawesi Barat e-mail: *1 fardinah@unsulbar.ac.id , 2 darmaekawati@unsulbar.ac.id , 3 hikmah@unsulber.ac.id , 4 ahmad.ansar@unsulbar.ac.id

## Abstrak
Dalam  lingkungan  mangsa  pemangsa,  terutama  yang  terdiri  atas  lebih  dari  satu  spesies  mangsa, terdapat beberapa perilaku pemangsa dalam berburu mangsa. Pada umumnya pemangsa lebih suka berburu di  suatu  habitat  untuk  beberapa  waktu  dan  kemudian  merubah  kemauannya  ke  habitat  lain.  Situasi  ini sangat  berbeda  ketika  spesies  mangsa  terdiri  dari  individu-individu  dalam  jumlah  yang  lebih  besar  dan dengan  ukuran  tubuh  yang  lebih  besar  dari  pemangsa  serta  memiliki  kemampuan  komunikasi  untuk membentuk pertahanan kelompok sehingga mangsa tersebut dapat bertahan dan melawan pemangsa. Selain itu  dapat  juga  dijumpai  dalam  suatu  lingkungan  bahwa  terdapat  spesies  yang  dapat  dipanen  untuk memenuhi  kebutuhan  manusia.  Penelitian  ini  bertujuan  untuk  menganalisis  kestabilan  model  dua  mangsa satu  pemangsa  dengan  pertahanan  kelompok  dan  pemanenan  linear  pada  mangsa  yang  terdiri  dari  tiga subpopulasi  yaitu  dua  jenis  spesies  mangsa  dan  satu  pemangsa.  Jenis  kestabilan  ditentukan  berdasarkan karakteristik  nilai  eigen  yang  diperoleh  dengan  menggunakan  kriteria  Routh-Hurwitz.  Dari  penelitian  ini diperoleh  bahwa  kepunahan  subpopulasi  pemangsa  dan  eksistensi  semua  subpopulasi  dapat  terjadi  jika memenuhi kondisi yang disyaratkan

## Kata kunci: Model Mangsa Pemangsa, Pertahanan Kelompok, Pemanenan Linear,  Kriteria Routh-Hurwitz

## 1. PENDAHULUAN
Model mangsa pemangsa pertama kali diperkenalkan oleh Alfred Lotka dan Vito Volterra pada tahun 1926 yang selanjutnya dikenal dengan model Lotka-Voltera. Dalam lingkungan mangsa pemangsa, terutama yang terdiri atas lebih dari satu spesies mangsa, terdapat beberapa perilaku pemangsa dalam berburu mangsa. Pada umumnya  pemangsa  lebih  suka  berburu  di  suatu  habitat  untuk  beberapa  waktu  dan  kemudian  merubah kemauannya  ke habitat lain. Salah satu alasan beralih dari pemangsa adalah semakin berkurangnya spesies mangsa di habitat pemangsa sedang mencari mangsa.
Situasi  ini  sangat  berbeda  ketika  spesies  mangsa  terdiri  dari  individu-individu  dalam  jumlah  yang lebih besar dan dengan ukuran tubuh yang lebih besar dari pemangsa. Selain itu, individu-individu mangsa tersebut  juga  memiliki  kemampuan  untuk  berkomunikasi  dan  bekerjasama  baik  antar  sesama  mangsa maupun  dengan mangsa  dari  spesies yang lain.  Kemampuan  berkomunikasi  dan bekerjasama  antar  sesama mangsa sangat berguna dalam hal mendeteksi keberadaan dan pergerakan pemangsa serta untuk membentuk pertahanan  kelompok  sehingga  mangsa  tersebut  dapat  bertahan  dan  melawan  pemangsa.  Sedangkan kemampuan  berkomunikasi  dan  bekerjasama  dengan  spesies  mangsa  lain  sangat  berguna  ketika  terjadi serangan  oleh  pemangsa  maka  mangsa  dari  spesies  lain  tersebut  dapat  membantu  pertahanan  kelompok dalam menghadapi pemangsa.
Dalam  perkembangannya,  model  Lotka-Voltera  ini  kemudian  mengalami  banyak  modifikasi. Pertumbuhan  mangsa  dan  pemangsa  menggunakan  fungsi  yang  lebih  kompleks  seperti  fungsi  logistik dimana pertumbuhan mangsa dibatasi oleh kapasitas batas lingkungan. Selain bentuk interaksi antara mangsa dan  pemangsa  yang  dikenal  dengan  fungsi  respon  juga  mengalami  perkembangan.  Fungsi  respon  dalam
https://jurnal.unsulbar.ac.id/index.php/saintifik
ekologi  diartikan  sebagai  jumlah  makanan  yang  dikonsumsi  oleh  pemangsa  sebagai  fungsi  kepadatan makanan  (Panja,  et  al.,  2017).  Tahun  1953  Holling  memperkenalkan  fungsi  respon  yang  dikenal  sebagai fungsi respon Holling tipe I, tipe II dan tipe III ( Mortoja, et al., 2018). Salah satu fungsi respon yang sering digunakan  adalah  fungsi  respon  Holling  Tipe  II  khususnya  fungsi  respon  yang  hanya  bergantung  pada populasi mangsa.
Disamping  itu,  model  mangsa  pemangsa  juga  dikembangkan  dengan  penelitian  model  interaksi antara  lebih  dari  satu  mangsa  dengan  pemangsa.  Penelitian  tersebut  salah  satunya  dilakukan  oleh    Tansky (Q.J.A  Khan,  2001)  yang  meneliti  model  matematika  dari  sistem  dua  mangsa  dan  satu  pemangsa  dengan sifat  pemangsa  yang beralih.  Penelitian lain  dilakukan  oleh  Freedman  dan  Wolkowicz  (Q.J.A  Khan, 2001) yang  mengembangkan  dengan  kasus  mangsa  memiliki  pertahanan  kelompok  yaitu  dengan  memberikan fungsi  mekanisme  peralihan  pemangsa.  Selain  itu,  model  mangsa  pemangsa  juga  dikembangkan  dengan melakukan pemanenan pada populasi mangsa atau pemangsa. Penelitian ini diantaranya telah dilakukan oleh Benjamin  Leard,  dkk  (2008)  yang  telah  mempelajari  dinamika  populasi  dengan  pemanenan  dari  satu  atau beberapa spesies sehingga dihasilkan berbagai model matematika dengan melibatkan pemanenan pada satu atau beberapa spesies populasi.
Artikel  ini  mengkaji  kembali  model  yang  telah  dipaparkan  oleh  Q.J.A  Khan  (2001)  dan memodifikasi model tersebut dengan menambahkan pemanenan pada populasi mangsa. Meskipun modifikasi yang  dilakukan  tidak  terlalu  signifikan  namun  hasil  analisis  memperlihatkan  bahwa  terjadi  perubahan signifikan pada dinamika dari model dibandingkan model rujukan.

## 2. METODE PENELITIAN
Metode  yang  digunakan  dalam  penelitian  ini  adalah  studi  literatur  dengan  prosedur  penelitian diberikan sebagai berikut:
1. Membangun model berdasarkan asumsi yang diperoleh dari literatur yang telah dikaji
2. Menentukan titik kesetimbangan model
3. Menganalisis titik kesetimbangan model
4. Melakukan simulasi model

## 3. HASIL DAN PEMBAHASAN

## 3.1 Pembentukan Model Dua Mangsa Satu Pemangsa dengan Pemanenan Linear pada Mangsa
Dalam  pembentukan  model  dua  mangsa  satu  pemangsa  dengan  pemanenan  linear  pada  mangsa diperlukan  asumsi-asumsi  yang  berlaku  dalam  ekosistemnya.  Asumsi-asumsi  tersebut  diberikan  sebagai berikut:

## 3.1.1. Asumsi Model
Dalam  pembentukan  model  dua  mangsa  satu  pemangsa  dengan  pemanenan  linear  pada  mangsa, ekosistemnya dianggap memenuhi asumsi-asumsi sebagai berikut:
1. Terdapat tiga spesies yaitu mangsa 𝑥 1 , mangsa 𝑥 2 dan pemangsa y dengan mangsa 𝑥 1 dan mangsa 𝑥 2 tersebut hidup di lingkungan berbeda dan memiliki ukuran tubuh yang lebih besar dari pemangsa.
2. Dalam keadaan tanpa pemangsa, populasi mangsa tumbuh secara logistik.
3. Mangsa memiliki makanan yang cukup.
4. Mangsa dapat berpindah diantara dua lingkungan mangsa.
5. Laju kematian pemangsa konstan dan mangsa merupakan makanan utama bagi pemangsa.
6. Mangsa 𝑥 1 dan  mangsa 𝑥 2 memiliki  kemampuan  berkomunikasi  dan  bekerjasama  antar  sesama mangsa  dalam  membentuk  pertahanan  kelompok.  Pertahanan  kelompok  ini  sangat  berguna  dalam mendeteksi keberadaan pemangsa dan melakukan perlawanan terhadap pemangsa.
7. Fungsi respon yang digunakan adalah fungsi respon Holling Tipe II
https://jurnal.unsulbar.ac.id/index.php/saintifik
8. Pemanenan hanya dilakukan pada mangsa 𝑥 1 dengan laju pemanenan linear.

## 3.1.2. Model Dasar Mangsa Pemangsa Lotka Volterra
Model  dalam  penelitian  ini  dibentuk  dari  model  dasar  mangsa  pemangsa  yaitu  model  Lotka Volterra yaitu,
𝑑𝑥 𝑑𝑡 = 𝑔𝑥 − 𝛽𝑥𝑦
𝑑𝑦
𝑑𝑡 = 𝛿𝛽𝑥𝑦 − 𝜇𝑦,
dengan x adalah kepadatan populasi mangsa, y adalah kepadatan populasi pemangsa, 𝑔 adalah koefisien kelahiran  mangsa, 𝛽 adalah  koefisien  pemangsaan, 𝜇 koefisien  kematian  pemangsa,  dan 𝛿 koefisien konversi mangsa ke pemangsa.
Populasi mangsa tidak selamanya meningkat tetapi dapat terjadi jika subpopulasi naik maka angka pertumbuhan  cenderung  turun.  Bahkan  untuk  subpopulasi  yang  cukup  besar,  bukan  mustahil  angka pertumbuhan  negatif.  Fenomena  ini  disebabkan  area  dan  fasilitas  hidup  terbatas.  Oleh  karena  itu, persamaan (1) kurang relevan karena populasi mangsa tidak selamanya naik atau turun. Populasi mangsa memiliki  kapasitas  batas  lingkungan  (disimbolkan k )  Jika  populasi  mangsa  meningkat  maka  daya tampung  atau  kapasitas  batas  akan  berkurang  (disimbolkan 𝑘 − 𝑥 ).  Akibatnya,  masih  tersisa 𝑘−𝑥 𝑘 atau 1 − 𝑥 𝑘 bagian  kapasitas  lingkungan.  Oleh  karena  itu,  pertumbuhan  populasi  mangsa  akan berbanding lurus dengan sisa kapasitas lingkungan, sehingga diperoleh persamaan:
𝑑𝑥 𝑑𝑡 = 𝑔𝑥 (1 − 𝑥 𝑘 ) − 𝛽𝑥𝑦
𝑑𝑦 𝑑𝑡 = 𝛿𝛽𝑥𝑦 − 𝜇𝑦
3.1.3.  Modifikasi  Model  Dasar  Mangsa  Pemangsa  Lotka  Volterra  menjadi  Model  Dua  Mangsa  Satu Pemangsa
Berdasarkan  asumsi  (4)  pada  pembahasan  (3.1.1)  diperoleh  bahwa  angka  pertumbuhan  populasi mangsa 𝑥 1 bertambah dengan adanya perpindahan mangsa 𝑥 2 dan angka pertumbuhan populasi mangsa 𝑥 2 bertambah dengan adanya perpindahan mangsa 𝑥 1 . Perpindahan  mangsa 𝑥 2 ke  lingkungan 𝑥 1 dipengaruhi  oleh  koefisien  perpindahan 𝑥 2 (disimbolkan 𝑒 2 ) dan peluang kesuksesan perpindahan 𝑥 2 (disimbolkan 𝑝 21 ). Dengan demikian pertambahan mangsa 𝑥 1 adalah  jumlah  mangsa 𝑥 2 yang  pindah  yaitu 𝑒 2 𝑝 21 𝑥 2 .  Dengan  prinsip  yang  sama  diperoleh pertambahan mangsa 𝑥 2 adalah jumlah mangsa 𝑥 1 yang pindah yaitu 𝑒 1 𝑝 12 𝑥 1 . Subtitusi 𝑒 2 𝑝 21 𝑥 2 dan 𝑒 1 𝑝 12 𝑥 1 ke  persamaan  (2)  diperoleh  model  dua  mangsa  satu  pemangsa sebagai berikut:
𝑑𝑥 1 𝑑𝑡 = 𝑥 1 𝑔 1 (1 − 𝑥 1 𝑘 1 ) +   𝑒 2 𝑝 21 𝑥 2 −   𝛽 1 𝑥 1 𝑥 2 𝑑𝑥 2 𝑑𝑡 = 𝑥 2 𝑔 2 (1 − 𝑥 2 𝑘 2 ) +   𝑒 1 𝑝 12 𝑥 1 −   𝛽 2 𝑥 1 𝑥 2 𝑑𝑦 𝑑𝑡 = 𝛿 1 𝛽 1 𝑥 1 𝑦 + 𝛿 2 𝛽 2 𝑥 2 𝑦 − 𝜇𝑦.
3.1.4.   Modifikasi  Model Dua  Mangsa  –  Satu  Pemangsa  menjadi  Model  Dua Mangsa  Satu  Pemangsa dengan Pertahanan Kelompok dan Pemanenan Linear pada Mangsa
Berdasarkan  asumsi  (6)  diperoleh  fungsi  pemangsaan 𝐻(𝑥) = 𝛽𝑥 1+𝛽𝑡 ℎ 𝑥 .  Dalam  hal  ini  diasumsikan koefisien  pemangsaan  yaitu 𝛽 dan  waktu  yang  dibutuhkan  untuk  mengkonsumsi  mangsa  yaitu 𝑡 ℎ berbanding terbalik yang dapat dituliskan sebagai 𝛽 = 1 𝑡 ℎ sehingga diperoleh 𝛽𝑡 ℎ = 1, sehingga peroleh
(3) (2) (1)
https://jurnal.unsulbar.ac.id/index.php/saintifik
fungsi  respon  untuk  mangsa 𝑥 1 yaitu 𝐻(𝑥 1 ) = 𝛽𝑥 1 1+𝑥 1 dan  fungsi  respon  untuk  mangsa 𝑥 2 yaitu 𝐻(𝑥 2 ) = 𝛽𝑥 2 1+𝑥 2 .
Fakta  menunjukkan  bahwa  fungsi  respon  untuk  mangsa 𝑥 1 juga  bergantung  pada  mangsa 𝑥 2 dan fungsi respon untuk mangsa 𝑥 2 juga bergantung pada mangsa 𝑥 1 ketika pemangsa mencari makan. Dengan demikian diperoleh fungsi pemangsaan untuk mangsa 𝑥 1 yaitu:
𝐻 ( 𝑥 1 𝑥 2 ) = 𝛽 1 ( 𝑥1 𝑥2 )𝑥 2 𝑦 1+( 𝑥1 𝑥2 ) atau 𝐻 ( 𝑥 1 𝑥 2 ) = 𝛽 1 𝑥 1 𝑥 2 𝑦 𝑥 1 +𝑥 2
dan fungsi pemangsaan untuk mangsa 𝑥 2 yaitu:
𝐻 ( 𝑥 2 𝑥 1 ) = 𝛽 2 ( 𝑥2 𝑥1 )𝑥 1 𝑦 1+( 𝑥2 𝑥1 ) atau 𝐻 ( 𝑥 2 𝑥 1 ) = 𝛽 1 𝑥 1 𝑥 2 𝑦 𝑥 1 +𝑥 2 .
Berdasarkan asumsi (7) diperoleh bahwa angka pertumbuhan populasi mangsa 𝑥 1 berkurang dengan adanya  pemanenan  linear  pada  mangsa 𝑥 1 .  Berkurangnya  populasi  mangsa 𝑥 1 ini  dipengaruhi  oleh koefisien  pemanenan  (dilambangkan ℎ̂). Dengan  demikian  mangsa 𝑥 1 yang  berkurang  akibat pemanenan linear yaitu ℎ̂𝑥 1 .
Berdasarkan penjelasan di atas, diperoleh persamaan diferensial model dua mangsa satu  pemangsa dengan pertahanan kelompok dan pemanenan  linear pada  mangsa 𝑥 1 yaitu: 𝑑𝑥 1 𝑑𝑡 = 𝑥 1 𝑔 1 (1 − 𝑥 1 𝑘 1 ) +   𝑒 2 𝑝 21 𝑥 2 − 𝛽 1 𝑥 1 𝑥 2 𝑦 𝑥 1 +𝑥 2 − ℎ̂𝑥 1 𝑑𝑥 2 𝑑𝑡 = 𝑥 2 𝑔 2 (1 − 𝑥 2 𝑘 2 ) +   𝑒 1 𝑝 12 𝑥 1 − 𝛽 2 𝑥 1 𝑥 2 𝑦 𝑥 1 +𝑥 2 𝑑𝑦 𝑑𝑡 = (−𝜇 + 𝛿 1 𝛽 1 𝑥 1 𝑥 2 𝑥 1 +𝑥 2 + 𝛿 2 𝛽 2 𝑥 1 𝑥 2 𝑥 1 +𝑥 2 ) 𝑦,
dengan: 𝑔 1 : koefisien kelahiran mangsa 𝑥 1 𝑔 2 : koefisien kelahiran mangsa 𝑥 2 𝑘 1 : kapasitas batas lingkungan untuk mangsa 𝑥 1 𝑘 2 : kapasitas batas lingkungan untuk mangsa 𝑥 2 𝑒 1,2 : koefisien perpindahan mangsa 𝑥 1 ke habitat 𝑥 2 𝑒 2,1 : koefisien perpindahan mangsa 𝑥 2 ke habitat 𝑥 1 𝑝 1,2 : peluang kesuksesan peralihan mangsa 𝑥 1 ke habitat 𝑥 2 𝑝 2,1 : peluang kesuksesan peralihan mangsa 𝑥 2 ke habitat 𝑥 1 𝛽 1 : koefisien pemangsaan pada mangsa 𝑥 1 𝛽 2 : koefisien pemangsaan pada mangsa 𝑥 2 𝛿 1 : koefisien konversi mangsa 𝑥 1 ke pemangsa 𝛿 2 : koefisien konversi mangsa 𝑥 2 ke pemangsa 𝜇 : koefisien kematian pemangsa
ℎ̂ : koefisien pemanenan.
Koefisien 𝑔 1 , 𝑔 2 , 𝑘 1 , 𝑘 2 , 𝑒 1,2 , 𝑒 2,1 , 𝑝 1,2 , 𝑝 2,1 , 𝛽 1 , 𝛽 2 , 𝛿 1 , 𝛿 2 , 𝜇 dan ℎ̂ adalah konstanta-konstanta positif.
Dengan  memisalkan 𝑔 1 =  𝑔, 𝛽 1 = 𝛽, 𝑒 1 𝑝 1,2 = 𝑒 , 𝑒𝑡 = 𝜏 , 𝑝 1 = 𝑔 𝑒 ,  𝛽 ′ = 𝛽 𝑒 , 𝑔 2 = 𝑔, 𝛽 2 = 𝛽 , 𝑒 2 𝑝 2,1 = 𝑒 , 𝛽 1 = 𝛽 2 = 𝛽,  𝛿 1 = 𝛿 2 = 𝛿 , 𝜇 1 = 𝜇 𝑒 dan 𝛿 ′ = 2𝛿 maka dari persamaan (4) diperoleh:
𝑑𝑥 1 𝑑𝜏 = 𝑝 1 𝑥 1 (1 − 𝑥 1 𝑘 ) +   𝑥 2 − 𝛽 ′ 𝑥 1 𝑥 2 𝑦 𝑥 1 + 𝑥 2 𝑑𝑥 2 𝑑𝜏 = 𝑝 1 𝑥 2 (1 − 𝑥 2 𝑘 ) + 𝑥 1 − 𝛽 ′ 𝑥 1 𝑥 2 𝑦 𝑥 1 +𝑥 2 (5) 𝑑𝑦 𝑑𝜏 = (−𝜇 1 + 𝛿 ′ 𝛽 ′ 𝑥 1 𝑥 2 𝑥 1 + 𝑥 2 ) 𝑦.

## 3.2. Kestabilan Titik Kesetimbangan
(4)
https://jurnal.unsulbar.ac.id/index.php/saintifik
Titik  kesetimbangan  persamaan  (5)  diperoleh  berdasarkan  (Verhulst,  1989)  dengan  menentukan solusi  yang  memenuhi 𝑑𝑥 1 𝑑𝜏 = 0 , 𝑑𝑥 2 𝑑𝜏 = 0, dan 𝑑𝑦 𝑑𝜏 = 0 .  Dari  persamaan  tersebut  diperoleh  titik kesetimbangan sebagai berikut:
(i) Titik  kesetimbangan 𝐸 1 (𝑥̅ 1 , 𝑥̅ 2 , 0) dengan 𝑥̅ 1 = −𝑥̅ 2 𝑝 1 (1 − 𝑥̅ 2 𝑘 ) dan 𝑥̅ 2 adalah  solusi  positif  dari 𝑝 1 3 𝑘 3 𝑥 2 3 − 2𝑝 1 3 𝑘 2 𝑥 2 2 − ( 𝑝 1 2 𝑘 − 𝑝 1 3 𝑘 − ℎ𝑝 1 𝑘 )𝑥 2 − (−𝑝 1 2 + ℎ𝑝 1 + 1) = 0 (ii) Titik kesetimbangan 𝐸 2 (𝑥̅ 1 , 𝑥̅ 2 , 𝑦̅) dengan 𝑥̅ 1 adalah solusi positif dari 𝛽 ′2 𝛿 ′2 𝑝 1 𝑥 1 3 + (𝛽 ′2 𝛿 ′2 𝑘 − 𝛽 ′2 𝛿 ′2 𝑘𝑝 1 − 2𝛽 ′ 𝛿 ′ 𝜇 1 𝑝 1 +𝛽 ′2 𝛿 ′2 𝑘ℎ)𝑥 1 2 + (3𝛽 ′ 𝛿 ′ 𝑘𝜇 1 𝑝 1 − 3𝛽 ′ 𝛿 ′ 𝑘𝜇 1 − 2𝛽 ′ 𝛿 ′ 𝑘𝜇 1 ℎ)𝑥 1 + 2𝑘𝜇 1 2 + 𝑘ℎ𝜇 1 2 − 2𝑘𝜇 1 2 𝑝 1 = 0, 𝑥̅ 2 = 𝜇 1 𝑥̅ 1 𝛿 ′ 𝛽 ′ 𝑥̅ 1 −𝜇 1 dan 𝑦̅ = 𝛿 ′ 𝑥̅ 1 𝜇 1 ( 𝜇 1 𝑝 1 𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 − 𝜇 1 2 𝑝 1 𝑥̅ 1 𝑘(𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 ) 2 + 1) , dengan 𝑥̅ 1 > 𝜇 1 𝛿 ′ 𝛽 ′ .
Selanjutnya  akan  dipelajari  dinamika  di  sekitar  titik-titik  kesetimbangan  tersebut.  Untuk mengidentifikasi  kestabilan  lokal  dari  titik  kesetimbangan,  dilakukan  pelinearan  terhadap  model  (5) berdasarkan (Perko,  1991) .  Hasil  pelinearan  memberikan  matriks  Jacobian  dari  model  (5)  yaitu  sebagai berikut:
𝑱 =
[ 𝑝 1 − 2𝑝 1 𝑥 1 𝑘 − 𝛽 ′ 𝑥 2 2 𝑦 (𝑥 1 + 𝑥 2 ) 2 − ℎ 1 − 𝛽 ′ 𝑥 1 2 𝑦 (𝑥 1 + 𝑥 2 ) 2 − 𝛽 ′ 𝑥 1 𝑥 2 𝑥 1 + 𝑥 2 1 − 𝛽 ′ 𝑥 2 2 𝑦 (𝑥 1 + 𝑥 2 ) 2 𝑝 1 − 2𝑝 1 𝑥 2 𝑘 − 𝛽 ′ 𝑥 1 2 𝑦 (𝑥 1 + 𝑥 2 ) 2 − 𝛽 ′ 𝑥 1 𝑥 2 𝑥 1 + 𝑥 2 𝛿 ′ 𝛽 ′ 𝑥 2 2 𝑦 (𝑥 1 + 𝑥 2 ) 2 𝛿 ′ 𝛽 ′ 𝑥 1 2 𝑦 (𝑥 1 + 𝑥 2 ) 2 −𝜇 1 + 𝛿 ′ 𝛽 ′ 𝑥 1 𝑥 2 𝑥 1 + 𝑥 2 ] Dengan mensubstitusi 𝐸 1 (𝑥̅ 1 , 𝑥̅ 2 , 0) ke matriks Jacobian J diperoleh:
𝐽 1 = [ 𝑝 1 − 2𝑝 1 𝑥̅ 1 𝑘 − ℎ 1 − 𝛽 ′ 𝑥̅ 1 𝑥̅ 2 𝑥̅ 1 + 𝑥̅ 2 1 𝑝 1 − 2𝑝 1 𝑥̅ 2 𝑘 − 𝛽 ′ 𝑥̅ 1 𝑥̅ 2 𝑥̅ 1 + 𝑥̅ 2
0 0 −𝜇 1 + 𝛿 ′ 𝛽 ′ 𝑥̅ 1 𝑥̅ 2 𝑥̅ 1 + 𝑥̅ 2 ] diperoleh persamaan karakteristik 𝜆 2 + (2𝑝 1 − 2𝑝 1 𝑘 (𝑥̅ 1 + 𝑥̅ 2 ) − ℎ)𝜆 + (𝑝 1 − 2𝑝 1 𝑥̅ 1 𝑘 − ℎ)(𝑝 1 − 2𝑝 1 𝑥̅ 2 𝑘 ) = 0.
Dengan kriteria Routh-Hurwitz (Gantmacher,1959) diperoleh bahwa titik kesetimbangan 𝐸 1 akan stabil dengan syarat 𝜇 1 > 𝛿 ′ 𝛽 ′ 𝑥̅ 1 𝑥̅ 2 𝑥̅ 1 +𝑥̅ 2 , ℎ > 2𝑝 1 − 2𝑝 1 𝑘 (𝑥̅ 1 + 𝑥̅ 2 ) dan (𝑝 1 − 2𝑝 1 𝑥̅ 1 𝑘 − ℎ)(𝑝 1 − 2𝑝 1 𝑥̅ 2 𝑘 ) > 0 .
Dengan mensubstitusi 𝐸 2 (𝑥̅ 1 , 𝑥̅ 2 , 𝑦̅) ke matriks Jacobian J diperoleh:
𝐽 2 =
[ 𝑀 1 − 𝐵 − 𝜇 1 𝛿 ′ 1 − 𝐶 𝐿 − 𝜇 1 𝛿 ′ 𝛿 ′ 𝐶 𝛿 ′ 𝐵 0 ]
dengan:
𝐿 = 𝑝 1 − 2𝜇 1 𝑝 1 𝑥̅ 1 𝑘(𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 ) − (𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 ) 2 𝛽 ′ 𝛿 ′ 𝜇 1 𝑥̅ 1 ( 𝜇 1 𝑝 1 𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 − 𝜇 1 2 𝑝 1 𝑥̅ 1 𝑘(𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 ) 2 + 1)
𝑀 = 𝑝 1 − 2𝑝 1 𝑥̅ 1 𝑘 − 𝜇 1 𝛽 ′ 𝛿 ′ 𝑥̅ 1 ( 𝜇 1 𝑝 1 𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 − 𝜇 1 2 𝑝 1 𝑥̅ 1 𝑘(𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 ) 2 + 1) − ℎ
𝐶 = 𝜇 1 𝛽 ′ 𝛿 ′ 𝑥̅ 1 ( 𝜇 1 𝑝 1 𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 − 𝜇 1 2 𝑝 1 𝑥̅ 1 𝑘(𝛽 ′ 𝛿 ′ 𝑥 1 −𝜇 1 ) 2 + 1) 𝐵 = (𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 ) 2 𝛽 ′ 𝛿 ′ 𝜇 1 𝑥̅ 1 ( 𝜇 1 𝑝 1 𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 − 𝜇 1 2 𝑝 1 𝑥̅ 1 𝑘(𝛽 ′ 𝛿 ′ 𝑥̅ 1 −𝜇 1 ) 2 + 1).
https://jurnal.unsulbar.ac.id/index.php/saintifik
Diperoleh persamaan karakteristik 𝜆 3 + 𝑎 1 𝜆 2 + 𝑎 2 𝜆 + 𝑎 3 = 0 dengan 𝑎 1 = −𝐿 − 𝑀 , 𝑎 2 = 𝐿𝑀 + 𝜇 1 𝐶 + 𝜇 1 𝐵 − 1 + 𝐵 + 𝐶 − 𝐵𝐶 dan 𝑎 3 = 𝜇 1 𝐶 + 𝜇 1 𝐵 − 𝜇 1 𝐵𝑀 − 𝜇 1 𝐶𝐿 − 2𝜇 1 𝐵𝐶 . Dengan kriteria Routh-Hurwitz diperoleh bahwa titik kesetimbangan 𝐸 2 akan stabil dengan syarat dengan syarat −𝐿 > 𝑀 dan 𝐿𝑀 + 𝜇 1 𝐶 + 𝜇 1 𝐵 − 1 + 𝐵 + 𝐶 − 𝐵𝐶 > 0 .

## 3.3. Simulasi Numerik
Untuk memperkuat hasil analisis yang diberikan sebelumnya, dilakukan simulasi numerik terhadap model  (5).  Karena  keterbatasan  data  yang  ada,  pemilihan  parameter  disesuaikan dengan  kondisi  kestabilan yang diberikan pada hasil sebelumnya.
Simulasi titik kesetimbangan E 1 dipilih parameter sebagai berikut. p = 0.1, k = 0.7, β = 0.2, δ = 0.02, μ = 0.4 dan h = 0.04 Diperoleh gambar sebagai berikut:
Gambar 1. Simulasi numerik disekitar titik kesetimbangan E 1
Berdasarkan  analisis  kestabilan  yang  dilakukan  sebelumnya  diperoleh  bahwa  titik  kesetimbangan E 1 stabil asimtotik dengan syarat seperti yang telah disebutkan pada bagian 3.2. Dari gambar tersebut terlihat bahwa subpopulasi pemangsa menuju kepunahan sedangkan subpopulasi mangsa akan eksis. Simulasi titik kesetimbangan E 2 dipilih parameter sebagai berikut. p = 0.05, k = 0.5, β = 0.2, δ = 0.1, μ = 0.04 dan h = 0.1 Diperoleh gambar sebagai berikut:
Gambar 2. Simulasi numerik disekitar titik kesetimbangan E 2
https://jurnal.unsulbar.ac.id/index.php/saintifik
Berdasarkan  analisis  kestabilan  yang  dilakukan  sebelumnya  diperoleh  bahwa  titik  kesetimbangan E 2 stabil asimtotik dengan syarat seperti yang telah disebutkan pada bagian 3.2. Dari gambar tersebut terlihat bahwa kehidupan populasi akan eksis.

## 4. KESIMPULAN
Model model dua mangsa satu pemangsa dengan pemanenan  linear pada  mangsa seperti ditunjukkan dalam model  (4).  Hasil  analisis  menunjukkan  bahwa  terdapat  dua  titik  kesetimbangan  yaitu  titik  kesetimbangan kepunahan  pemangsa  dan  titik  kesetimbangan  eksistensi  semua  subpopulasi.  Dari  analisis  titik kesetimbangan  yang  dilakukan  diperoleh  bahwa    kepunahan  subpopulasi  pemangsa  dan  eksistensi  semua subpopulasi  dapat  terjadi  jika  memenuhi  kondisi  yang  disyaratkan.  Di  akhir  artikel  telah  ditunjukkan simulasi  numerik  berupa  gambar  deret  waktu  di  sekitar  setiap  titik  kesetimbangan  serta  interpretasi biologisnya.

## DAFTAR PUSTAKA
Benjamin,  L.,Catherine,  L.,&  Jorge,  R.,  2008,  Dynamics  Of  Ratio-Dependent  Predator-Prey  Models  With Nonconstant Harvesting, Discrete And Continuous Dynamical Systems Series, Vol 1, pp. 303–315.
Gantmacher, F. R., 1959, Applications of The Theory of Matrices , Intersciences Publishers, New-York, USA.
Khan,  Q.J.A.  dan  Balakrishnan,  B.,  2001,  Analysis  of  a  Predator  –  Prey  System  with  Predator  Switching , Bulletin of Mathematical Biology , Vol. 38 , pp. 449-458 .
Mortoja, S. G., Panja, P., & Mondal, S. K., 2018, Dynamics of a Predator-Prey Model with Stage Structure on Both Species and Anti-Predator Behavior, Informatics in Medicine Unlocked, Vol 10, pp. 50–57.
Panja, P., Mondal, S. K., & Chattopadhyay, J., 2017, Dynamical Effects of Anti-predator Behaviour of Adult Prey  in  a  Predator-Prey  Model  with  Ratio-Dependent  Functional  Response,  Asian  Journal  of Mathematics and Physics, Vol 1, pp. 19–32.
Perko, L., 1991, Differential Equation and Dynamical System, New York, SpringerVerlag Berlin Heidelberg.
Verhulst, F.,1990, Nonlinear Differential Equations and Dynamical Systems, Springer – Verlag, New- York, USA.