Teknik Interpolasi PieceWise Menggunakan RStudio

Interpolasi Piecewise

Interpolasi dengan polinomial sering memberikan hasil yang tidak dapat diterima. Interpolasi polinomial yang dihasilkan dari sejumlah besar data titik biasanya berderajat tinggi. Polinomial berderajat tinggi pada umumnya bersifat osilatif (grafiknya naik turun secara cepat). Akibatnya, perubahan data pada interval kecil dapat menyebabkan fluktuasi besar pada keseluruhan interval. Karena alasan ini, biasanya interpolasi hanya menggunakan polinomial berderajat rendah.

Interpolasi piecewise menawarkan alternatif lain. Pada interpolasi piecewise, pada titik yang berbeda sepanjang kurva, nilai fungsi lebih mungkin lebih baik didekati menggunakan dua atau lebih interpolasi. Pada metode ini kita akan membuat fungsi interpolasi ditiap antara dua titik observasi.

Pada sub-Chapter ini akan dijelaskan 2 buah metode interpolasi piecewise, yaitu: interpolasi linier piecewise dan interpolasi kubik spline. Interpolasi pertama dilakukan menggunakan persamaan linier, sehingga kurva yang terbentuk bukan merupakan kurva kontinu. Interpolasi selanjutnya dilakukan menggunakan persamaan polinomial berderajat tinggi sehingga kurva yang dihasilkan lebih halus (tidak ada sudut siku pada setiap titik).

8.2.1 Interpolasi Linier Piecewise

Interpolasi linier piecewise merupakan interpolasi yang menggunakan pendekatan interpolasi linier. Fungsi linier akan dibentuk pada setiap dua titik observasi. Untuk lebih memahaminya perhatikan kembali Gambar 8.2. Pada gambar tersebut sebelumnya kita telah membentuk persamaan kudratik untuk menghubungkan titik-titik tersebut. Dibanding menggunakan persamaan polinomial seperti kuadratik tersebut, interpolasi piecewise akan menghubungkan tiap dua titik observasi tersebut dengan garis lurus.

Algoritma Interpolasi Linier Piecewise

Tentukan set titik berpasangan (x,y)(x,y) yang telah diurutkan berdasarkan nilai sumbu xx.
Hitung mm pada setiap dua titik berdekatan menggunakan Persamaan (8.4)
Hitung bb pada setiap dua titik berdekatan menggunakan Persamaan (8.5)
Definiskan fungsi linier berdasarkan nilai mm dan bb
Hitung yy dengan cara substitusi nilai xx pada persamaan linier untuk melakukan interpolasi nilai yy yang ingin dicari.
Untuk melakukan ekstrapolasi dengan titik observasi diluar rentang titik diketahui, gunakan persamaan linier yang berada pada bagian ujung terdekat dengan nilai xx yang hendak dicari nilai yy-nya.

Berdasarkan algoritma tersebut, kita dapat menyusun fungsi pada R untuk membentuk persamaan linier piecewise. Berikut adalah sintaks yang digunakan:

pwise_linterp <- function(x, y){
  n <- length(x)-1
  
  y <- y[order(x)]
  x <- x[order(x)]
  
  m_vec <- b_vec <- c()
  
  for(i in 1:n){
    m <- (y[i+1]-y[i]) / (x[i+1]-x[i])
    b <- y[i+1] - m*x[i+1]
    m_vec <- c(m_vec, m)
    b_vec <- c(b_vec, b)
  }
  
  return(list(b = b_vec, m = m_vec))
}

Contoh 8.7 Tentukan persamaan-persamaan linier yang dihasilkan dari titik observasi yang ditampilkan pada Contoh 8.5?

Jawab:

Untuk menentukan persamaan-persamaan linier yang menghubungkan setiap titik, kita akan menggunakan fungsi pwise_linterp() yang telah kita buat. Berikut adalah sintaks yang digunakan:

x <- c(-1, 1, 0, -2)
y <- c(-2, 2, -1, -1)

pwise_linterp(x, y)

## $b
## [1] -3 -1 -1
## 
## $m
## [1] -1  1  3

Jika persamaan-persamaan yang terbentuk tersebut divisualisasikan akan terlihat seperti pada Gambar 8.3.

Interpolasi linier piecewise empat titik (Sumber:Howard, 2017).

Gambar 8.3: Interpolasi linier piecewise empat titik (Sumber:Howard, 2017).

R juga menyediakan fungsi untuk melakukan interpolasi linier piecewise. Fungsi approxfun() dapat digunakan untuk melakukan interpolasi tersebut. Fungsi approxfun() hanya memerlukan dua input yaitu vektor xx dan vektor xx. Untuk lebih memahaminya, kita akan menggunakan kembali data yang disajikan pada Contoh 8.5 untuk membuat fungsi linier piecewise. Berikut adalah sintaks yang digunakan:

f <- approxfun(x, y)

# tentukan nilai y jika x= 0
f(0)

## [1] -1

# tentukan nilai y jika x = 0.5
f(0.5)

## [1] 0.5

8.2.2 Interpolasi Spline Kubik

Jika menggunakan interpolasi polinomial berderajat satu (sebuah garis) lebih dari beberapa interval merupakan peningkatan dari satu baris interpolasi, dan jika menggunakan polinomial berderajat tinggi juga merupakan peningkatan dari satu garis interpolasi tunggal, maka dapat disimpulkan bahwa penggunaan polinomial berderajat tinggi pada selang beberapa interval juga akan menjadi peningkatan dalam proses interpolasi. Dalam beberapa kasus, hal tersebut benar, tetapi kita masih menghadapi sudut tajam di mana masing-masing kurva interpolasi tergabung (interpolasi linier piecewise). Sudut tajam ini mencegah diferensiasi dan pada prakteknya tidak dapat digunakan untuk memodelkan beberapa fungsi di dunia nyata, seperti roller coaster span.

Interpolasi spline kubik memecahkan masalah ini. Interpolasi ini akan memberikan kurva tergabung yang halus. Hal tersebut juga membuat spline terintegrasi. Karena setiap bagian individu diwakili oleh kurva kubik (polinomial derajat 3), maka masing-masing bagian individu juga dapat dianalisis sebagai kurva kubik. Dengan asumsi ada nn titik data untuk interpolasi, kita akan mendefinisikan SiSi sebagai fungsi polinomial kubik yang mewakili kurva pada domain [xi;xi+1][xi;xi+1]. Kemudian untuk nn titik observasi, ada n−1n−1 interpolasi polinomial kubik.

Bentuk umum seri polinomial dituliskan pada Persamaan (8.9).

Si=di(x−xi)3+ci(x−xi)2+bi(x−xi)+ai(8.9)(8.9)Si=di(x−xi)3+ci(x−xi)2+bi(x−xi)+ai

Ini mengarah ke 4 nilai yang tidak diketahui tidak diketahui, yaitu: aiai, bibi, cici, dan didi pada setiap Persamaan (8.9). Oleh karena itu, ada 4(n−1)=4n−44(n−1)=4n−4 nilai yang tidak diketahui. Karena kita ingin spline membentuk garis kontinu dan dapat didiferensiasi, ada satu set persamaan yang menentukan spline kubik:

Si(xi)=yi, i=1,……,n−1(8.10)(8.10)Si(xi)=yi, i=1,……,n−1

Si(xi+1)=yi+1, i=1,……,n−1(8.11)(8.11)Si(xi+1)=yi+1, i=1,……,n−1

S′i(xi+1)=S′i+1(xi), i=1,……,n−2(8.12)(8.12)Si′(xi+1)=Si+1′(xi), i=1,……,n−2

S′′i(xi+1)=S′′i+1(xi), i=1,……,n−2(8.13)(8.13)Si″(xi+1)=Si+1″(xi), i=1,……,n−2

Persamaan (8.10) dan (8.11) sudah cukup jelas. Persyaratan ini memastikan bahwa jika kita mengevaluasi spline di salah satu node internal, hasil yang akan kita peroleh merupakan jawaban yang telah ditentukan, yaitu spline yang dievaluasi pada xixi untuk beberapa ii adalah yiyi, dan setiap komponen spline bergabung dengan rapi. Persamaan (8.12) memastikan bahwa kita memiliki turunan pertama yang berkelanjutan di setiap simpul internal. Ini mencegah terbentuknya sudut tajam pada tiap node. Persamaan (8.13) memastikan turunan kedua juga kontinu, dimana kondisi ini menguntungkan karena itu berarti turunan pertama itu sendiri dapat didiferensiasi juga.

Kondisi ini menyebabkan ada 4n−64n−6 kondisi yang harus kita penuhi. Jika 4n−44n−4 kita yang tidak diketahui dipecahkan sebagai sebuah matriks, dan akhirnya matriks tersebut akan terpecahkan, matriks akan menjadi kurang ditentukan. Kita dapat menyelesaikan kondisi tersebut dengan memasukkan dua ketentuan tambahan. Dengan splines kubik, secara normal adalah menentukan akhir di kedua ujung untuk mencapai dua kondisi tambahan. Untuk contoh ini, dua kondisi yang akan kita tambahkan adalah S′′i(xi)=0Si″(xi)=0 dan S′′i(xn)=0Si″(xn)=0. Kedua kondisi ini memastikan bahwa pada titik akhir, turunan pertamanya linier dan oleh karena itu fungsi spline berlanjut ke arah yang sudah berjalan. Interpolasi spline kubik ini disebut juga sebagai “natural spline”.

Awalnya, kita dapat melihat bahwa setiap polinomial kubik SiSi digeser ke kanan oleh unit xixi atau bergeser ke kiri jika xixi negatif. Pada xixi nilai fungsi adalah aiai yang berarti ai=yiai=yi untuk setiap nilai ii. Menyelesaikan sisa koefisien yang ada akan lebih kompleks, tetapi sekarang 3n3n tidak diketahui dengan 3n3n kondisi. Derivasi penuh tersedia dari berbagai sumber, tetapi secara garis besar dari Persamaan (8.12) kita mendapatkan S′i(xi+1)=S′i+1Si′(xi+1)=Si+1′, kemudian S′i+1−S′i(x)=0Si+1′−Si′(x)=0, dan kita dapat mensubstitusi Persamaan (8.9) ke dalam kedua komponen, pemecahan untuk didi dalam hal ini xixi, yiyi, dan cici. Proses yang sama dapat direplikasi dengan Persamaan (8.13) dan bibi. Hasilnya adalah matriks tridiagonal, seperti yang ada pada Chapter 6.3.3, yang dapat dipecahkan untuk menemukan koefisien. Terdapat sebuah matriks, A, sedemikian rupa sehingga,

AC=V(8.14)(8.14)AC=V

dimana AA merupakan matriks tridiagonal. Pada matriks tridiagonal ini (u1,u2,…,un−1)(u1,u2,…,un−1) merupakan vektor UU, dan vektor MM, LL, dan VV ditentukan dengan cara serupa. Matriks ini sedemikian rupa,

ui=li=xi+1−xi(8.15)(8.15)ui=li=xi+1−xi

kecuali pada u0=lnu0=ln. Lebih jauh, diagonal utama DD ditentukan menggunakan Persamaan (8.16).

mi=2(xi+1−xi+xi−xi−1)(8.16)(8.16)mi=2(xi+1−xi+xi−xi−1)

kecuali pada d0=dn=1d0=dn=1. Akhirnya vektor VV ditentukan dengan Persamaan (8.17).

vi=3(yi+1−yixi+1−xi−yi−yi−1xi−xi−1)(8.17)(8.17)vi=3(yi+1−yixi+1−xi−yi−yi−1xi−xi−1)

kecuali pada v0=vn=0v0=vn=0. Sehingga,

⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣m1u1000l1m2u2000l2⋱⋱000⋱⋱un−1000ln−1mn⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦C=⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣v1v2⋮⋮vn⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦(8.18)(8.18)[m1u1000l1m2u2000l2⋱⋱000⋱⋱un−1000ln−1mn]C=[v1v2⋮⋮vn]

Penyelesaian Persamaan (8.18) akan menghasilkan vektor koefisien cc, sehingga koefisien bb dapat dihitung menggunakan Persamaan (8.19).

bi=yi+1−yixi+1−xi−xi+1−xi3(2ci+ci+1)(8.19)(8.19)bi=yi+1−yixi+1−xi−xi+1−xi3(2ci+ci+1)

Dengan menggunakan vektor CC yang sudah diketahui, koefisien dd dapat dihitung menggunakan Persamaan (8.20).

di=ci+1−ci3(xi+1−xi)(8.20)(8.20)di=ci+1−ci3(xi+1−xi)

Algoritma Interpolasi Spline Kubik

Tentukan set titik berpasangan (x,y)(x,y).
Tentukan koefisien aiai menggunakan Persamaan (8.10), dimana ai=yiai=yi.
Hitung elemen diagonal bawah (lili) dan elemen diagonal atas (uiui) matriks tridiagonal menggunakan Persamaan (8.15), dimana u0=ln=0u0=ln=0.
Hitung elemen diagonal utama (mimi) menggunakan Persamaan (8.16), dimana d0=dn=1d0=dn=1.
Hitung elemen vektor VV menggunakan Persamaan (8.17), dimana v0=vn=0v0=vn=0.
Susunlah elemen lili, uiui, dan mimi menjadi matriks tridiagonal AA.
Deifinisikan persamaan linier seperti pada Persamaan (8.14)
Selesaikan sistem persamaan linier pada Persamaan (8.14) sehingga diperoleh vektor CC yang merupakan kumpulan koefisien cc.
Hitung koefisien bibi menggunakan Persamaan (8.19)
Hitung koefisien didi menggunakan Persamaan (8.20).
Bentuk seri persamaan polinomial menggunakan elemen koefisien aa, bb, cc, dan dd yang telah dihitung.
Untuk melakukan ekstrapolasi dengan titik observasi diluar rentang titik diketahui, gunakan persamaan polinomial yang berada pada bagian ujung terdekat dengan nilai xx yang hendak dicari nilai yy-nya.

Berdasarkan algoritma tersebut, kita dapat menyusun suatu fungsi pada R untuk mencari seri persamaan polinomial derajat tiga. Fungsi tersebut adalah sebagai berikut:

cubic_spline <- function(x, y){
  n <- length(x)
  d_vec <- b_vec <- a_vec <- rep(0, n-1)
  vec <- rep(0, n)
  delta_x <- delta_y <- rep(0, n-1)
  
  ## Menghitung nilai selisih dan vektor A
  for(i in 1:(n-1)){
    a_vec[i] <- y[i]
    delta_x[i] <- x[i+1] - x[i]
    delta_y[i] <- y[i+1] - y[i]
  }
  
  ## Menyusun matriks tridiagona
  Au <- c(0, delta_x[2:(n-1)])
  Am <- c(1, 2*(delta_x[1:(n-2)]+delta_x[2:(n-1)]), 1)
  Al <- c(delta_x[1:(n-2)], 0)
  
  vec[0] <- vec[n] <- 0
  for(i in 2:(n-1))
    vec[i] <- 3 * (delta_y[i]/delta_x[i] - 
                   delta_y[i-1]/delta_x[i-1])
  
  ## penyelesaian tridiagonal matriks
  nm <- length(Am)
  Al <- c(NA , Al)
  
  ### forward sweep
  Au[1] <- Au[1] / Am[1]
  vec[1] <- vec[1] / Am[1]
  for(i in 2:(n - 1)){
      Au[i] <- Au[i] / (Am[i] - Al[i] * Au[i - 1])
      vec[i] <- (vec[i] - Al[i] * vec[i - 1]) /
                (Am[i] - Al[i] * Au[i - 1])
  }
  vec[n] <- (vec[n] - Al[n] * vec[n - 1])/
            (Am[n] - Al[n] * Au[n - 1])
  
  ### backward sweep
  c_vec <- rep.int (0, n)
  c_vec[n] <- vec[n]
  for(i in (n - 1) :1)
      c_vec[i] <- vec[i] - Au[i] * c_vec[i + 1]
  
 ## Hitung vektor B dan D dari vektor C
 for(i in 1:(n-1)){
   b_vec[i] <- (delta_y[i]/delta_x[i])-
               (delta_x[i]/3)*(2*c_vec[i]+c_vec[i+1])
   d_vec[i] <- (c_vec[i+1]-c_vec[i]) / (3*delta_x[i])
 }
  
 return(list(a = a_vec, b = b_vec, c = c_vec[-n], d = d_vec))
}

Contoh 8.8 Tentukan persamaan-persamaan spline kubik yang dihasilkan dari titik observasi yang ditampilkan pada Contoh 8.5?

Jawab:

Untuk menentukan persamaan-persamaan linier yang menghubungkan setiap titik, kita akan menggunakan fungsi cubic_spline() yang telah kita buat. Berikut adalah sintaks yang digunakan:

x <- c(-2, -1, 0, 1)
y <- c(-1, -2, -1, 2)

cubic_spline(x, y)

## $a
## [1] -1 -2 -1
## 
## $b
## [1] -1.4 -0.2  2.2
## 
## $c
## [1] 0.0 1.2 1.2
## 
## $d
## [1]  0.4  0.0 -0.4

Sebagai contoh untuk domain [0,1][0,1], persamaan spline kubiknya adalah −0,4x3+1,2x2+2,2x−1−0,4x3+1,2x2+2,2x−1. Jika persamaan-persamaan yang terbentuk tersebut divisualisasikan akan terlihat seperti pada Gambar 8.4.

Interpolasi spline kubik empat titik (Sumber:Howard, 2017).

Gambar 8.4: Interpolasi spline kubik empat titik (Sumber:Howard, 2017).

Pada R juga terdapat fungsi splinefun() untuk melakukan interpolasi spline. Format fungsi tersebut adalah sebagai berikut:

splinefun(x, y = NULL,
          method = c("fmm", "periodic", "natural", "monoH.FC", "hyman"),
          ties = mean)

Catatan:

x,y : vektor titik yang akan dilakukan interpolasi

method : spesifikasi jenis spline yang digunakan. Nilai yang mungkin antara lain: "fmm", "natural", "periodic", "monoH.FC" dan "hyman".

ties : metode untuk menangani nilai imbang pada titik yang akan diinterpolasi.

Untuk melakukan interpolasi spline kubik, metode yang digunakan adalah "natural". Berikut adalah contoh interpolasi menggunakan kembali data pada Contoh 8.5:

spline <- splinefun(x, y, method="natural")

# uji dengan nilai x yang sama
spline(x)

## [1] -1 -2 -1  2

FUngsi splinefun() menghasilkan nilai yang sama persis dengan proses interpolasi yang telah kita lakukan.

Teknik Interpolasi PieceWise Menggunakan RStudio

Heny Rimadana, Prof. Dr. Suhartono, M.Kom

11/15/2021

Interpolasi Piecewise

8.2.1 Interpolasi Linier Piecewise

8.2.2 Interpolasi Spline Kubik

REFERENSI