1 Program Studi Teknik Elektro (Mekatronika), Universitas Katolik Parahyangan,  Jl. Ciumbuleuit no. 94, Bandung 40141, INDONESIA  (tlp: 022-203-2655; fax: 022-203-1110; e-mail: ttamba@unpar.ac.id)

## Estimasi Parameter Model Nonlinear Menggunakan Analisis Sensitivitas dan Pengoptimalan Berbasis Turunan
Tua A. Tamba 1
Abstract—Model estimation based on the observation data of a system’s states is an important subject in the study of dynamical systems. Maximum likelihood (ML) estimation is a stochastic estimation method which can be used to obtain an optimal set of parameter based on noisy measurements. This paper describes the method and implementation of the ML estimator  to identify an optimal parameter set in a discrete-time nonlinear state space model. In particular, the optimal parameter set is defined as the value that minimizes the error between the actual and estimated model outputs of the system. This paper discusses a gradient- based optimization that is equipped with sensitivity analysis method for searching such a parameter set. Simulation results which describe an implementation of the proposed estimation method in a nonlinear system model are also discussed.
Intisari—Estimasi model berdasarkan data hasil pengukuran variabel sistem merupakan salah satu topik penting dalam studi sistem dinamik. Estimasi maximum likelihood  (ML) merupakan suatu teknik estimasi stokastik yang dapat digunakan untuk menentukan  nilai parameter optimal pada  suatu  model  sistem dinamik  berdasarkan data observasi  yang bercampur dengan sinyal gangguan acak. Makalah  ini memaparkan implementasi estimator ML untuk mengidentifikasi parameter sistem dinamik yang dimodelkan dengan  persamaan ruang keadaan waktu diskrit.  Secara khusus,  parameter optimal didefinisikan  sebagai nilai yang meminimalkan selisih antara data pengukuran aktual dan  sinyal keluaran yang dihasilkan oleh model sistem dengan nilai  parameter  hasil  estimasi.  Penentuan parameter optimal dilakukan  dengan metode  pengoptimalan  berbasis gradien yang dikombinasikan dengan analisis sensitivitas. Hasil simulasi numerik yang mengilustrasikan  implementasi metode  estimasi yang diusulkan pada suatu model chemostat juga didiskusikan.
Kata Kunci—Estimasi parameter,  metode  stokastik,  estimator maximum likelihood, pengoptimalan, analisis sensitivitas.

## I. P ENDAHULUAN
Salah satu topik penting dalam studi sistem dinamik adalah terkait penentuan model analitik suatu sistem berdasarkan data observasi  variabel sistem yang tersedia  atau dapat diukur. Aspek terkait nonlinearitas dinamika sistem, data pengukuran yang tercampur  sinyal gangguan acak,  serta  keterbatasan jumlah data yang dapat digunakan dalam analisis merupakan beberapa contoh aspek mendasar yang berkontribusi terhadap kompleksitas proses estimasi model sistem dinamik [1], [2].
Metode  estimasi parameter Bayesian adalah  salah satu teknik  inferensi  berbasis pendekatan stokastik yang sering digunakan untuk menentukan  model sistem dinamik
berdasarkan data pengukuran yang tercampur sinyal gangguan acak [3]-[7]. Teknik estimasi ini pada esensinya menggunakan aturan Bayes  dalam penentuan  rapat probabilitas dari  suatu model sistem dinamik terparameterisasi berdasarkan data hasil pengukuran yang tersedia. Secara khusus,  inferensi terkait statistik  dari  vektor  parameter p yang ditinjau ditentukan menggunakan  aturan Bayes dan diformulasikan sebagai  nilai ekspektasi dari p atau fungsi dari p sebagai berikut [3]-[5]. ( ) , ) ( ) | ( | ) ( Ε dp p y p f y p p ψ ψ ∫ = (1)
dengan  E( . | . ) merupakan  operator ekspektasi bersyarat  dan ) | ( y p f adalah  fungsi  distribusi posterior dari p untuk data observasi y dari sistem  yang tersedia. Isu  mendasar  dalam implementasi  metode  estimasi parameter Bayesian adalah terkait perhitungan nilai integral  pada (1)  yang melibatkan fungsi distribusi parameter yang dianalisis. Secara umum, solusi analitik integral  ini  tidak dapat ditentukan sehingga harus  diperkirakan  dengan  metode  numerik berbasis simulasi Monte Carlo. Sayangnya, selain terbentur  potensi  terjadinya permasalahan  akumulasi kesalahan, metode  simulasi Monte Carlo membutuhkan beban komputasi yang sangat besar jika digunakan pada model  dengan dimensi parameter tinggi  [6]- [7].
Estimator maximum likelihood (ML) merupakan salah satu metode  alternatif untuk proses  estimasi  parameter.  Pada estimator  ML, permasalahan estimasi  diformulasikan secara formal sebagai penentuan nilai parameter optimal  (dengan kriteria tertentu) suatu model terparameterisasi (yang diyakini mendasari  dinamika  sistem yang ditinjau) berdasarkan data pengukuran  variabel keadaan sistem yang ada  [8], [9]. Beberapa keuntungan teknik estimasi ini adalah dimungkinkannya penggunaan informasi awal ( prior information ) yang tersedia terkait model dan sinyal gangguan pada sistem (rentang nilai atau distribusi) serta kemudahan mengasimilasi informasi baru dalam proses estimasi.
Estimator ML menentukan estimasi nilai q buah parameter ) , , ( 1 q p p p  = pada model  sistem dinamik berdasarkan data pengukuran  sistem  dan  keyakinan  awal  ( prior belief )  yang mungkin diketahui tentang  . p Keyakinan awal tersebut secara formal  dinyatakan  sebagai  fungsi distribusi prior ) ( p f yang mengindikasikan  distribusi  nilai probabilitas untuk  setiap kemungkinan nilai p ketika  data pengukuran terbaru belum ada. Ketergantungan nilai pengukuran ) , , ( 1 m y y y  = sistem pada parameter p dinyatakan sebagai fungsi rapat probabilitas bersyarat ). | ( Pr ) | ( p y p y L = Fungsi ) | ( p y L disebut  fungsi likelihood dan  dapat  digunakan untuk memperbarui ) ( p f setiap kali data pengukuran y yang  baru  diperoleh. Proses pembaruan dilakukan menggunakan aturan Bayes berikut [3]- [5]
) | ( ) ( ) | ( ) | ( ) ( ) | ( p y L p f y p f p y L p f y p f ∝ ⇒ ℵ = (2) dengan dp p y L p f p ∫ = ℵ ) | ( ) ( adalah konstanta normalisasi.
Fungsi likelihood dalam hal ini memodifikasi fungsi distribusi prior menjadi posterior ketika data pengukuran baru diperoleh. Dengan demikian, nilai parameter optimal ditentukan sebagai nilai p yang memaksimalkan fungsi likelihood ). | ( p y L
Makalah ini memaparkan metode  dan hasil  implementasi estimator  ML  untuk  mengestimasi  parameter  model sistem nonlinear  dalam bentuk persamaan ruang keadaan waktu diskrit.  Pada  metode  yang dibahas, parameter optimal didefinisikan sebagai nilai  parameter yang meminimalkan selisih  antara  data pengukuran aktual dan sinyal keluaran sistem  yang dihasilkan oleh nilai parameter yang diestimasi. Perhitungan parameter optimal dilakukan dengan metode pengoptimalan berbasis gradien yang dikombinasikan dengan analisis sensitivitas. Contoh implementasi metode  estimasi parameter  model  yang diusulkan  diilustrasikan melalui simulasi numerik penentuan parameter pada model chemostat .

## II. D ESKRIPSI DAN M ODEL S ISTEM
Ditinjau  sistem dinamik dengan Ζ ∈ n variabel keadaan n t x ℜ ∈ ) ( dan Ζ ∈ x  parameter x x p
 ℜ ∈ dengan  masing- masing ℜ dan Ζ adalah himpunan bilangan nyata dan bulat. Untuk waktu 0 ≥ t serta  nilai awal  , 0 x dinamika  sistem diasumsikan memenuhi persamaan diferensial berikut.
, ) 0 ( ), , ), ( ( ) ( 0 x x t p t x f t x x = =  (3)
Makalah ini meninjau metode estimasi parameter di model (3) berdasarkan data pengukuran variabel keadaan x ( t ). Karena  keterbatasan alat  ataupun kesulitan dalam proses pengukuran, data pengukuran seringkali  hanya diperoleh untuk sebagian dari variabel keadaan x ( t )  [9]-[11].  Data pengukuran yang diperoleh juga umumnya  tidak langsung dalam bentuk nilai masing-masing variabel keadaan,  tetapi merupakan nilai kombinasi atau fungsi beberapa variabel keadaan.  Untuk memperhitungkan kondisi ini, pengukuran dari sistem pada (3) diasumsikan secara abstrak sebagai suatu fungsi ) ), ( ( y p t x g yang memetakan variabel keadaan menjadi m pengukuran . ) ( m t y ℜ ∈ Dalam hal ini, model  pengukuran untuk sistem pada  (3) untuk  waktu 0 ≥ t dapat  dituliskan sebagai berikut [10], [11].
), , ), ( ( ) ( t p t x g t y y = (4)
dengan y y p  ℜ ∈ dan Ζ ∈ y  adalah parameter pada fungsi pengukuran.
Data pengukuran pada  (4) yang diperoleh dari sensor umumnya merupakan sinyal diskrit dalam  jumlah terhingga serta  tercampur sinyal gangguan acak. Jika y ij menyatakan data pengukuran ke- j yang diperoleh pada waktu t i ( i =1,..., T ), maka y ij dapat diasumsikan sebagai  gabungan  dari nilai pengukuran yang sebenarnya  ( y j ( t i ))  dan  suatu proses acak tertentu  ( v ij ),  sehingga  model pengukuran dapat dituliskan sebagai (5) [10], [12].
, ) ( ij i j ij v t y y + = (5)
dengan ij v diasumsikan memiliki distribusi Gaussian dengan nilai rata-rata nol dan varians 2 ij σ .
Persamaan  (3)  sampai  (5)  merupakan  dasar untuk  proses estimasi parameter  model  sistem dinamik berdasarkan data pengukuran variabel keadaan  yang  lengkap  ( m = n )  ataupun parsial ). ( n m < Permasalahan estimasi model dapat didefinisikan sebagai proses estimasi nilai  parameter
k y x p x p p p = ≡ ) , , ( 0 untuk K k , , 1  = di  sistem  pada  (3) sampai (4) berdasarkan data pengukuran pada (5).

## III. E STIMASI P ARAMETER M ODEL N ONLINEAR

## A. Estimator Maximum Likelihood (ML)
Salah satu metode  estimasi  parameter model sistem dinamik  berdasarkan data pengukuran yang tercampur sinyal gangguan acak adalah dengan menggunakan estimator ML. Untuk sistem dengan model dinamik  seperti pada (3)  dan (4) serta  pengukuran  pada  (5), estimator ML  menentukan parameter optimal p * sebagai nilai p yang memaksimalkan kecenderungan ( likelihood ) diperolehnya data pengukuran y tersebut.  Secara matematis, p * merupakan    nilai p yang memaksimalkan fungsi likelihood ) | ( p y L terhadap  data pengukuran yang ditinjau [2], [5], [7], [10]-[12].
Jika  sinyal acak pada  (5)  memenuhi  distribusi Gaussian dengan  nilai  rata-rata nol dan varians 2 ij σ ,  maka  fungsi likelihood ) | ( p y L dapat dituliskan sebagai berikut.
( ) ( ) [ ] . ), , ( 2 1 exp 2 1 2 1 1         − − = ∏∏ = = p p t x g y p y L i j ij ij T i m j ij σ π σ (6)
Dari sudut pandang efisiensi komputasi, maksimalisasi fungsi ) | ( p y L di  (6) lebih mudah dilakukan dengan minimalisasi nilai negatif logaritma ) | ( p y L berikut  [2],  [5],
[7], [10]-[12].
( ) [ ]            − −         − = −
∏ ∏ = = T i m j ij ij p R p y L 1 1 2 ) ( 2 1 exp log 2 1 log log π σ ( ) , ) ( 2 1 2 log 1 1 2 ∑∑ = = + = T i m j ij ij p R π σ (7)
dengan ) ( p R ij merupakan fungsi residu sebagai berikut.
( ) [ ] . ), , ( ) ( 1 p p t x g y p R i j ij ij ij − = − σ (8)
Karena komponen ) 2 log( π σ ij pada  (7) tidak bergantung pada nilai p , dapat disimpulkan bahwa nilai parameter optimal p * pada estimator ML memenuhi (9) [2], [5], [7], [10]-[12].
. ) ( 2 1 min arg * 2 1 1       ∑ ∑ ≡ = = p R p ij T i m j p (9)
Persamaan (9)  mengindikasikan bahwa  parameter optimal pada  estimator ML adalah  nilai  yang meminimalkan jumlah kuadrat  kesalahan ( sum of squared errors ) antara data pengukuran sinyal keluaran aktual dan data pengukuran sinyal keluaran yang dihasilkan oleh nilai parameter yang diestimasi.

## B. Pengoptimalan Berbasis Turunan dan Analisis Sensitivitas
Salah satu metode  komputasi nilai p *  pada  (9) adalah dengan teknik pengoptimalan  berbasis gradien. Pada teknik pengoptimalan  ini, nilai solusi optimal  di setiap iterasi diperbarui berdasarkan nilai gradien suatu fungsi residu R ij ( p ). Lebih spesifik, kandidat nilai  solusi  optimal ditentukan sebagai nilai parameter yang menghasilkan perubahan gradien fungsi residu terkecil atau kurang dari  kriteria nilai  toleransi kesalahan yang ditentukan sejak awal pada algoritme pengoptimalan [6]-[7], [11]-[13].
Pada implementasi algoritme  pengoptimalan, perhitungan nilai gradien fungsi residu dapat dilakukan dengan menggunakan berbagai perangkat lunak komputasi,  misalnya MATLAB [14]. Secara khusus, nilai gradien R ij ( p ) terhadap p untuk data pengukuran pada waktu τ = t ditentukan sebagai berikut [5], [7], [11]-[13].
( ) , 1 1 ), , ( 1 ) ( 1 1         ∂ ∂ + ∂ ∂ ∂ ∂ − =         ∂ ∂ + ∂ ∂ ∂ ∂ − = ∂ ∂ − = ∂ ∂ ∑ ∑ = = n i k j k i i j ij n i j n n j ij k i j ij k ij p g p x x g p g p x x g p p p t x g p p R τ τ τ τ τ τ σ σ σ (10) dengan m j , , 1  = dan . , , 1 K k  = Komponen ) ( i j dx dg dan ) ( k j dp dg pada  (8)  merupakan  matriks Jacobian
n m g x J × ℜ ∈ dan
K m g p J × ℜ ∈ dari  fungsi ) , ( p x g terhadap vektor variabel keadaan x dan parameter p berikut.
, 1 1 1 1               ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ = ∂ ∂ ≡ n m m n i j g x x g x g x g x g x g J      (11a)
. 1 1 1 1       
       ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ = ∂ ∂ ≡ K m m K k j g p p g p g p g p g p g J      (11b)
Sementara itu,
) ( k i dp dx dalam (10) adalah komponen dari matriks sensitivitas ) , ( k i S sebagai berikut.
. 1 1 1 1 K n K n n K k i x p p x p x p x p x p x S × ℜ ∈               ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ = ∂ ∂ ≡      (12)
Matriks Jacobian pada  (11) dapat dihitung  menggunakan komputasi simbolik (misalnya pada MATLAB) karena fungsi
) , , ( t p x g y pada (4)  telah  ditentukan/diketahui  sejak  awal.  Di sisi lain, matriks sensitivitas pada (12) umumnya hanya dapat ditentukan dengan aproksimasi numerik. Dalam hal ini, aproksimasi dilakukan dengan terlebih dahulu mengintegralkan (3) secara numerik dan hasilnya  kemudian digunakan untuk menghitung nilai  matriks  sensitivitas (12). Kelemahan teknik aproksimasi numerik ini adalah kemungkinan timbulnya permasalahan akumulasi  kesalahan numerik serta  beban komputasi yang pengaruhnya akan menjadi sangat signifikan  ketika  digunakan pada model dengan banyak parameter. Untuk mengatasi permasalahan ini, matriks sensitivitas dapat ditentukan  dengan  lebih efisien menggunakan persamaan diferensial berikut [2], [5], [7], [11]- [13].
      ∂ ∂ ∂ ∂ = k i x p p x t dt dS       ∂ ∂ ∂ ∂ = t x p i k k x p t p t x f ∂ ∂ = ) , ), ( ( k i k i n i i x i p f p x x t p t x f ∂ ∂ + ∂ ∂ ∂ ∂ = ∑ = 1 ) , ), ( (
(13)
dengan  nilai  awal ). 0 )( ( ) 0 ( k i x p dp dx S = Komponen  matriks
) ) ( ( i i x f ∂ ⋅ ∂ dan ) ) ( ( k i p f ∂ ⋅ ∂ dapat ditentukan secara analitik berdasarkan fungsi eksplisit ) ( ⋅ i f yang ditinjau. Dalam hal ini, matriks sensitivitas dapat ditentukan dengan mengintegralkan (3) dan (12) secara paralel terhadap variabel waktu sehingga beban komputasi tidak dipengaruhi oleh  jumlah parameter pada model yang ditinjau.

## C. Algoritme Estimasi Parameter Model
Metode estimasi parameter dengan pengoptimalan berbasis gradien dan analisis sensitivitas yang diusulkan  diringkaskan pada  Algoritme  I.  Pada  Algoritme  I, proses estimasi dimulai dengan penentuan nilai awal parameter 0 p dan  nilai batas ε sebagai kriteria konvergensi untuk terminasi proses pengoptimalan. Berdasarkan nilai parameter terkini, solusi (3) dan (13) dihitung  secara numerik  dan  kemudian  digunakan untuk menentukan nilai residu (8) dan turunan (10). Pengoptimalan  kemudian dapat  dilakukan  untuk menentukan nilai parameter optimal (9).  Pengoptimalan  ini  dilakukan
secara iteratif sampai  kriteria konvergensi ε terpenuhi. Dalam hal ini, nilai estimasi optimal didefinisikan sebagai nilai parameter dengan batas kesalahan yang memenuhi . ε
Algoritme I [7]-[8], [11]-[12] Estimasi parameter optimal pada (9)
Masukan: Model sistem (3) dan (5), nilai parameter awal 0 p dan nilai batas ε untuk kriteria konvergensi pengoptimalan
Keluaran: Estimasi parameter optimal * p
1. Dengan menggunakan nilai p dari iterasi sebelumnya, tentukan nilai integral dari persamaan (3) dan (13) secara numerik.
2. Tentukan nilai residu (8) dan turunan (10).
3. Tentukan parameter optimal (9) dengan pengoptimalan berbasis gradien (misalnya fungsi “ lsqnonlin” pada MATLAB).
4. Perbaharui nilai estimasi parameter p dan kembali ke langkah 1. Lakukan iterasi 1-4 hingga kriteria konvergensi pengoptimalan ε tercapai.

## IV. S IMULASI N UMERIK
Implementasi  metode  estimasi parameter dengan estimator ML yang diusulkan pada makalah ini  diilustrasikan pada proses penentuan model dinamik untuk interaksi suatu sistem mikroorganisme yang dibiakkan  di dalam chemostat , seperti ditunjukkan pada Gbr. 1.
Gbr. 1 Skematika sistem chemostat .

## A. Deskripsi dan Model Chemostat
Chemostat merupakan  suatu instrumen laboratorium  yang dapat digunakan untuk membiakkan dan mengamati dinamika perkembangan  bakteri atau mikroorganisme ketika diberi bahan makanan atau asupan nutrisi  terbatas.  Sebagaimana ditunjukkan pada Gbr. 1, suatu chemostat umumnya memiliki tiga komponen utama sebagai berikut [15]-[16].
1. Bejana pemasok  ( feed  vessel ),  yaitu  bejana  berisi larutan yang mengandung  nutrisi yang dibutuhkan oleh mikroorganisme yang dibiakkan  untuk berkembang. Salah satu dari nutrisi ini dapat diatur berjumlah terbatas (disebut nutrisi pembatas/kontrol) untuk melihat pengaruhnya pada  perkembangan  mikroorganisme yang dibiakkan.
2. Bejana pembiakan ( media vessel ),  yaitu  bejana berisi larutan nutrisi tempat pembiakan mikroorganisme.
3. Bejana penampungan ( collection vessel ),  yaitu  bejana penampung pembuangan dari bejana pembiakan. Larutan pada  bejana ini mengandung organisme yang jumlahnya dapat diukur dan kemudian digunakan untuk mengestimasi jumlah mikroorganisme pada  bejana pembiakan.
Selain tiga bejana utama  tersebut, chemostat juga dilengkapi dengan beberapa instrumen seperti lampu, pipa saluran udara steril,  serta pengaduk untuk mengatur kondisi lingkungan pendukung perkembangan organisme pada bejana pembiakan.
Salah satu contoh penggunaan  chemostat  adalah  untuk pemodelan  interaksi  mikroorganisme  pemangsa dan mangsa ( predator & prey ) antara alga hijau chlorella dan zooplankton brachionus . Chlorella merupakan  alga hijau bersel tunggal yang  berkembang biak melalui proses fotosintesis dengan bantuan 2 CO , air, cahaya, dan nutrisi. Brachionus merupakan zooplankton air tawar yang pertumbuhannya bergantung pada konsumsi  alga hijau seperti chlorella .  Pembiakan chlorella dan brachinous secara bersamaan pada chemostat dengan demikian dapat memicu interaksi jenis pemangsa ( brachionus ) dan mangsa ( chlorella ).  Lebih jauh, pertumbuhan chlorella dapat dibatasi dengan mengontrol asupan nutrisi tertentu pada bejana pemasok, sedangkan pertumbuhan brachionu s dibatasi oleh jumlah chlorella yang tersedia untuk dikonsumsi.
Makalah ini meninjau model dinamik untuk interaksi antara chlorella dan brachionus pada chemostat sebagai berikut [15].
( ) ( ) C N f N N t N C in ) ( − − = δ  (14a) ( ) C B C f C N f t C ) ( ) ( 2 1 δ ρ − − =  (14b)
( ) R m R C f t R ) ( ) ( 2 λ δ + + − =  (14c)
( ) R m R C f t B ) ( ) ( 2 + − = δ  (14d)
Pada (14), N in adalah  konsentrasi asupan nitrogen N yang dialirkan dari bejana pemasok  secara konstan, δ adalah  laju aliran larutan yang  keluar dari bejana pembiakan, sementara masing-masing m dan λ adalah  laju mortalitas brachionus akibat kematian mendadak dan penuaan. Variabel C menyatakan  konsentrasi chlorella sedangkan R dan B masing-masing  menyatakan  konsentrasi brachionus yang dapat dan tidak dapat  bereproduksi.  Fungsi f 1 ( N )  dan f 2 ( C ) merupakan fungsi logistik sebagai berikut.
, ) ( 1 N K N b N f C C + = . ) ( 2 C K C b C f B B + = (15)
Dengan  mendefinisikan ] ; ; ; [ 4 3 2 1 x x x x x = dengan , 1 N x = ,
2 C x = , 3 R x = dan B x = 4 sebagai  vektor  variabel keadaan dan ] ; ; ; ; ; ; ; ; [ 9 8 7 6 5 4 3 2 1 p p p p p p p p p p x = dengan , 1 δ = p , 2 in N p = , 3 B b p = , 4 C b p = , 5 B K p = , 6 C K p = , 7 m p =
, 8 λ = p ρ = 9 p sebagai vektor parameter sistem, maka model pada  (14)  dan  (15) dapat ditulis  dalam bentuk model (3) sebagai berikut.
( ) ( ) ( ) ( ) ( )             + − + + − − − − − = = ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )) ( ( ) ( ) ; ); ( ( ) ( 3 7 1 3 2 2 3 8 7 1 3 2 2 2 1 9 4 2 2 2 1 1 2 1 1 2 1 t x p p t x t x f t x p p p t x t x f t x p p t x t x f t x t x f t x f t x p p t p t x f t x x  (16) dengan kondisi awal )]. 0 ( ); 0 ( ); 0 ( ); 0 ( [ ) 0 ( 4 3 2 1 0 x x x x x x ≡ =

## B. Estimasi Parameter Model dengan Data Hasil Simulasi
Untuk mengilustrasikan metode  estimasi parameter  yang diusulkan di makalah ini, model (14)-(16) disimulasikan pada MATLAB untuk memperoleh data yang diasumsikan sebagai hasil pengukuran  dari chemostat [14].  Pada simulasi ini diasumsikan hanya ada satu kloning di keseluruhan populasi chlorella .  Data pengukuran hasil simulasi  tersebut kemudian digunakan oleh estimator ML untuk estimasi parameter model.
Model (14) disimulasikan pada MATLAB dengan kondisi awal ] 1 , 0 ; 1 , 0 ; 45 ; 80 [ 0 = x dan nilai parameter  pada  Tabel  I.
Data pengukuran dari chemostat diasumsikan hanya tiga, yaitu konsentrasi nitrogen di dalam bejana pembiakan dan total jumlah sel chlorella dan brachionus .  Dalam hal ini, model pengukuran (4) untuk sistem ini memenuhi persamaan berikut.
( ) ). ( 1 1 0 0 0 0 1 0 0 0 0 1 , ), ( ) ( t x t p t x g t y y           = = (17)

## Gbr. 2 Data pengukuran hasil simulasi.
Data pengukuran dari simulasi kemudian ditambah dengan sinyal acak yang memenuhi distribusi Gaussian dengan nilai rata-rata nol dan varians 0,1. Estimator ML disimulasikan
pada  enam puluh  data yang diasumsikan sebagai sampel pengukuran harian selama enam puluh  hari. Selain itu, diasumsikan pula bahwa parameter 1 p dan 2 p diketahui nilainya berdasarkan konfigurasi awal pengoperasian chemostat . Dengan demikian, estimator ML hanya digunakan untuk mengestimasi nilai tujuh parameter tersisa yang tidak diketahui.
T ABEL I N ILAI P ARAMETER UNTUK S IMULASI M ODEL (16)
δ in N B b C b B K C K m λ ρ
0,04 80 2,25 3,3 15 4,3 0,055 0,4 0,25 T ABEL II N ILAI P ARAMETER A WAL DAN O PTIMAL DALAM E STIMASI B b C b B K C K m λ ρ p 0 4,08 7,31 41,52 5,25 0,07 0,49 0,62 p * 2,116 3,29 13,74  4,849  0,056 0,41 0,252
Gbr.  3  Perbandingan model yang diidentifikasi estimator ML dengan data pengukuran serta batas confidence interval .
Gbr. 2 menunjukkan data pengukuran yang diperoleh dari hasil simulasi. Pada gambar ini, solusi model (16) ditunjukkan oleh kurva dengan garis terputus,  sedangkan  sampel  data pengukuran yang tercampur sinyal gangguan acak ditunjukkan sebagai lingkaran kecil. Data pengukuran diskrit yang tercampur sinyal gangguan acak ini kemudian dianalisis menggunakan estimator ML.
Algoritme  I diimplementasikan pada model (16)-(17)  dan data pengukuran di Gbr. 2 untuk mengestimasi nilai parameter
3 p hingga . 9 p Nilai awal parameter p 0 pada simulasi ditunjukkan di Tabel  II yang diperoleh melalui penambahan sinyal acak pada nilai p di Tabel I. Nilai parameter optimal p * hasil  estimasi  Algoritme  I juga ditunjukkan pada Tabel II. Prediksi rentang ketidakpastian pada respons model yang diestimasi  akibat ketidakpastian nilai parameter optimal yang dihasilkan  menggunakan  estimator ML digambarkan sebagai sabuk kurva abu-abu di Gbr. 3. Kurva rentang ketidakpastian
0 10 20 30 40 50 60 1 1.5 0 10 20 30 40 50 60 4 5 6 0 10 20 30 40 50 60 1.2 1.6 2 0 10 20 30 40 50 60 0.8 1 1.2 1.4 1.6 0 10 20 30 40 50 60 4 5 6 0 10 20 30 40 50 60 1.2 1.6 2
ini ditentukan berdasarkan propagasi nilai error berikut [11]- [13].
( ) ( ) , ) ( ) ( ) ( T x p x p t S t S t y Cov Ψ ≡ (18)
dengan ) ( t S x p adalah matriks  sensitivitas  pada  (12) dan 1 2 − Γ = Ψ merupakan suatu parameter yang berkaitan dengan Fisher Information Matrix (FIM) Γ berikut.
. 2 k ij p R ∂ ∂ = Γ (19)
Perlu diperhatikan bahwa FIM 1 2 − Γ = Ψ pada (18) merupakan salah satu ukuran jumlah informasi terkait p yang terkandung dalam data pengukuran y ( t ). Grafik pada Gbr. 3 ini juga menunjukkan batas 95% confidence interval (CI(95%)) hasil estimasi parameter pada  setiap data yang digunakan. Nilai interval ini ditentukan sebagai berikut [11]-[13].
. ) ( 95 , 1
CI(95%) * Ψ ± = diag p (20)
Terlihat  pada  Gbr. 3  bahwa  respons  dari model dengan parameter optimal yang diperoleh oleh estimator ML dapat menggambarkan dinamika sistem  sebagaimana ditunjukkan oleh  data pengukuran  aktual dari sistem.  Hasil simulasi ini menunjukkan bahwa estimator ML yang diusulkan merupakan salah satu metode yang efektif untuk mengestimasi parameter model  sistem  dinamik  berdasarkan data observasi yang tercampur dengan sinyal gangguan acak.

## V. K ESIMPULAN
Makalah ini telah  memaparkan  metode  dan implementasi estimator ML dalam proses estimasi  nilai  parameter  optimal pada model sistem nonlinear  berbentuk persamaan ruang keadaan waktu diskrit. Secara spesifik, parameter optimal didefinisikan sebagai nilai yang meminimalkan selisih antara data  keluaran  aktual  dari sistem dan  data  keluaran  model ketika menggunakan nilai parameter optimal yang diestimasi. Teknik  pengoptimalan  berbasis gradien yang dikombinasikan dengan  teknik  analisis sensitivitas  diusulkan sebagai metode komputasi untuk penentuan nilai parameter optimal berdasarkan data pengukuran dari sistem yang telah tercampur dengan sinyal gangguan acak.  Hasil simulasi numerik yang mengimplementasikan estimator ML yang diusulkan untuk
penentuan  parameter  optimal pada model chemostat mengilustrasikan potensi metode yang diajukan.
R EFERENSI
[1] J.R. Raol, G. Girija, dan J. Singh, Modelling & Parameter Estimation of Dynamic Systems , IET Control Engineering Series, Vol. 65, United Kingdom: The Institution of Engineering and Technology, 2004.
[2] J.V. Beck dan K.J. Arnold, Parameter Estimation in Engineering and Science , Wiley Series in Probability &  Mathematical Statistics,  New York, USA: Wiley, 1977.
[3] S.D. Silvey, Statistical Inference ,  Baltimore, USA:  Penguin Books, 1970.
[4] D.M. Malakoff, “Bayes Offers ‘New’ Way to Make Sense of Numbers,” Science , Vol. 286, hal. 1460–1464, 1999.
[5] R. Aster  dan  C. Thurber, Parameter Estimation & Inverse Problems , Cambridge, USA: Academic Press, Vol. 90, 2011.
[6] J.D.  Hamilton, Time Series Analysis ( Vol . 2),  Princeton, USA: Princeton University Press, 1994.
[7] W.H.  Press, Numerical Recipes: The Art of  Scientific Computing , 3rd Ed., New York, USA: Cambridge University Press, 2007.
[8] M.A. Tanner, Tools for Statistical Inference: Methods for  Exploration of Posterior Distributions and Likelihood Functions , 2nd Ed.,  New York, USA: Springer-Verlag, 1993.
[9] A.  Raue, C. Kreutz, T. Maiwald, J. Bachmann, M. Schilling, U. Klingmüller,  dan  J. Timmer,  “Structural  & Practical Identifiability Analysis  of Partially Observed Dynamical Models by Exploiting the Profile Likelihood,” Bioinformatics ,  Vol. 25,  No. 15, hal. 1923–1929, 2009.
[10] S. Hengl,  C. Kreutz,  J. Timmer,  dan  T. Maiwald,  “Data-Based Identifiability Analysis of Nonlinear Dynamical Models,”
Bioinformatics , Vol. 23, No. 19, hal. 2612–2618, 2007.
[11] F. Geier,  G. Fengos, F. Felizzi, dan  D. Iber., “Analyzing  and Constraining Signaling Networks: Parameter Estimation for the User,” Computational Modeling of Signaling Networks ,  Totowa, NJ,  USA: Humana Press, hal. 23–39 , 2012.
[12] K.  Jaqaman  dan  G. Danuser, “Linking  Data  to Models: Data Regression,” Nature Review Molecular Cellular Biology ,  Vol. 7, No. 11, hal. 813–819, 2006.
[13] B. Efron dan  R.  Tibshirani,  “Bootstrap  Methods  for Standard Errors, Confidence Intervals and  Other Measures  of Statistical Accuracy,” Statistical Science , Vol. 1, No. 1, hal. 54–75 , 1986.
[14] MATLAB  Optimization Toolbox Version 7.2 Release R2015a , The MathWorks, Inc., 2015.
[15] G.F.  Fussmann, S.P. Ellner, K.W. Shertzer,  dan  N.G. Hairston, “Crossing the Hopf Bifurcation  in a Live  Predator-Prey  System,” Science , Vol. 290, hal. 1358 – 1360, 2000.
[16] K.W. Shertzer, S.P. Ellner, G.F. Fussmann,  dan  N.G. Hairston, “Predator-Prey  Cycles in an Aquatic  Microcosm:  Testing  Hypotheses of  Mechanism,” Journal of Animal Ecology ,  Vol.  71,  hal.  802–815, 2002.