Soal No 1

Diketahui pdf (probability density function) dari distribusi Kumaraswamy adalah sebagai berikut:

\[f\left ( x\mid a,b \right )=abx^{a-1}\left ( 1-x^{\left ( a \right )} \right )^{b-1}\]

dimana \(x\in \left ( 0,1 \right )\). Jika diketahui \(a=5\) dan \(b=3\).

bagian a

Hitunglah \(E\left ( X \right )\) dengan menggunakan metode Trapezoida untuk \(n=4\), Simpson untuk \(n=4\) dan four-point Gauss Quadrature tanpa menggunakan software apapun. (Boleh diketik ataupun ditulis tangan kemudian difoto)

Jawaban bagian a

Metode Trapezoida untuk n=4

Diketahui:

  • \(x\in \left ( 0,1 \right )\)
  • \(a=5\)
  • \(b=3\)
  • \(f\left ( x\mid a,b \right )=abx^{a-1}\left ( 1-x^{\left ( a \right )} \right )^{b-1}\)

Substitusi nilai \(a\) dan \(b\) ke dalam fungsi, \(f\left ( x\mid 5,3 \right )=5 \cdot 3 \cdot x^{5-1}\left ( 1-x^{\left ( 5 \right )} \right )^{3-1}\) \(\ni\)

\[f\left ( x\mid 5,3 \right )=15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\]

Nilai \(E\left ( X \right )\)

\(E\left ( X \right )=x \cdot f\left ( x\mid 5,3 \right )\) \(E\left ( X \right )=x \cdot 15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\) \(\ni\) \[E\left ( X \right )=15 x^{5}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\]

Panjang setiap sub interval dimana \(x\in \left ( 0,1 \right )\) adalah sebagai berikut: \[h=\frac{b-a}{n}=\frac{1-0}{4}=\frac{1}{4}\]

Karena dalam soal diperintahkan untuk menghitung nilai \(E\left ( X \right )\) dengan menggunakan metode trapezoidal maka, \(E\left ( X \right )\) dianggap sebagai fungsi. Sehingga, karena n=4 maka bentuk dari metode trapezoidal adalah: \[T_{4}=\frac{h}{2}\left [ E\left ( X_{0} \right )+2E\left ( X_{1} \right )+2E\left ( X_{2} \right )+2E\left ( X_{3} \right )+E\left ( X_{4} \right ) \right ]\]

Menghitung fungsi \(E\left ( X \right )\) pada setiap \(x_{i}\):

\(E\left ( X_{0} \right )=E\left ( 0 \right )=15 \cdot 0^{5}\left ( 1-0^{\left ( 5 \right )} \right )^{2}=0\) \(E\left ( X_{1} \right )=E\left ( \frac{1}{4} \right )=15 \cdot \frac{1}{4}^{5}\left ( 1-\frac{1}{4}^{\left ( 5 \right )} \right )^{2}=0.014619841\) \(E\left ( X_{2} \right )=E\left ( \frac{2}{4} \right )=15 \cdot \frac{2}{4}^{5}\left ( 1-\frac{2}{4}^{\left ( 5 \right )} \right )^{2}=0.439910889\) \(E\left ( X_{3} \right )=E\left ( \frac{3}{4} \right )=15 \cdot \frac{3}{4}^{5}\left ( 1-\frac{3}{4}^{\left ( 5 \right )} \right )^{2}=2.070616786\) \(E\left ( X_{4} \right )=E\left ( 1 \right )=15 \cdot 1^{5}\left ( 1-1^{\left ( 5 \right )} \right )^{2}=0\)

sehingga nilai

\(T_{4}=\frac{\left ( \frac{1}{4} \right )}{2}\left [ 0 + 2 \cdot \left ( 0.014619841 \right ) + 2 \cdot \left ( 0.439910889 \right ) + 2 \cdot \left ( 2.070616786 \right ) + 0 \right ]\)

\(T_{4}=\frac{1}{8}\left [ 0 + 0.029239682 + 0.879821777 + 4.141233573 + 0 \right ]\)

\(T_{4}=\frac{1}{8}\left [ 5.050295033 \right ]\)

\(T_{4}=0.631286879\)

sehingga: \[\int_{0}^{1}x\cdot15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\approx T_{4}=0.631286879\]

Metode Simpson untuk n=4

Diketahui:

  • \(x\in \left ( 0,1 \right )\)
  • \(a=5\)
  • \(b=3\)
  • \(f\left ( x\mid a,b \right )=abx^{a-1}\left ( 1-x^{\left ( a \right )} \right )^{b-1}\)

Substitusi nilai \(a\) dan \(b\), \(f\left ( x\mid 5,3 \right )=5 \cdot 3 \cdot x^{5-1}\left ( 1-x^{\left ( 5 \right )} \right )^{3-1}\) \(\ni\)

\[f\left ( x\mid 5,3 \right )=15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\]

Nilai \(E\left ( X \right )\)

\(E\left ( X \right )=x \cdot f\left ( x\mid 5,3 \right )\) \(E\left ( X \right )=x \cdot 15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\) \(\ni\) \[E\left ( X \right )=15 x^{5}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\]

Panjang setiap sub interval dimana \(x\in \left ( 0,1 \right )\) adalah sebagai berikut: \[h=\frac{b-a}{n}=\frac{1-0}{4}=\frac{1}{4}\]

Karena dalam soal diperintahkan untuk menghitung nilai \(E\left ( X \right )\) dengan menggunakan metode simpson maka, \(E\left ( X \right )\) dianggap sebagai fungsi. Sehingga, karena n=4 maka bentuk dari metode simpson adalah: \[S_{4}=\frac{h}{3}\left [ E\left ( X_{0} \right )+4E\left ( X_{1} \right )+2E\left ( X_{2} \right )+4E\left ( X_{3} \right )+E\left ( X_{4} \right ) \right ]\]

Menghitung \(E\left ( X \right )\) pada setiap \(x_{i}\):

\(E\left ( X_{0} \right )=E\left ( 0 \right )=15 \cdot 0^{5}\left ( 1-0^{\left ( 5 \right )} \right )^{2}=0\) \(E\left ( X_{1} \right )=E\left ( \frac{1}{4} \right )=15 \cdot \frac{1}{4}^{5}\left ( 1-\frac{1}{4}^{\left ( 5 \right )} \right )^{2}=0.014619841\) \(E\left ( X_{2} \right )=E\left ( \frac{2}{4} \right )=15 \cdot \frac{2}{4}^{5}\left ( 1-\frac{2}{4}^{\left ( 5 \right )} \right )^{2}=0.439910889\) \(E\left ( X_{3} \right )=E\left ( \frac{3}{4} \right )=15 \cdot \frac{3}{4}^{5}\left ( 1-\frac{3}{4}^{\left ( 5 \right )} \right )^{2}=2.070616786\) \(E\left ( X_{4} \right )=E\left ( 1 \right )=15 \cdot 1^{5}\left ( 1-1^{\left ( 5 \right )} \right )^{2}=0\)

sehingga nilai

\(S_{4}=\frac{\left ( \frac{1}{4} \right )}{3}\left [ 0 + 4 \cdot \left ( 0.014619841 \right ) + 2 \cdot \left ( 0.439910889 \right ) + 4 \cdot \left ( 2.070616786 \right ) + 0 \right ]\)

\(S_{4}=\frac{1}{12}\left [ 0 + 0.058479365 + 0.879821777 + 8.282467145 + 0 \right ]\)

\(S_{4}=\frac{1}{12}\left [ 9.220768288 \right ]\)

\(S_{4}=0.768397357\)

sehingga: \[\int_{0}^{1}15 x^{5}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\approx S_{4}=0.768397357\]

Metode four-point Gauss Quadrature

Metode gaussian quadratur memiliki algoritma:

  1. Menentukan fungsi \(f\left ( x \right )\) dan selang integrasinya \(\left [ a,b \right ]\).
  2. Melakukan transformasi fungsi tersebut hingga diperoleh fungsi dengan selang \(\left [ -1,1 \right ]\) dengan menggunakan persamaan: \[\int_{a}^{b}f\left ( x \right )dx=\frac{b-a}{2}\int_{-1}^{1}f\left ( \frac{t(b-a)+a+b}{2} \right )dt\]
  3. Menentukan orde polinomial \(n\) yang akan digunakan.
  4. Lakukan proses integrasi dengan mengalikan faktor bobot \(c\) dengan \(f\left ( x_{i} \right )\) dengan menggunakan persamaan : \[\int_{-1}^{1}f\left ( x \right )dt\approx C_{1}f\left ( x_{1} \right )+C_{2}f\left ( x_{2} \right )+...+C_{n}f\left ( x_{n} \right )\]

Diketahui:

  • \(x\in \left ( 0,1 \right )\)
  • \(a=5\)
  • \(b=3\)
  • \(f\left ( x\mid a,b \right )=abx^{a-1}\left ( 1-x^{\left ( a \right )} \right )^{b-1}\)

Substitusi nilai \(a\) dan \(b\), \(f\left ( x\mid 5,3 \right )=5 \cdot 3 \cdot x^{5-1}\left ( 1-x^{\left ( 5 \right )} \right )^{3-1}\) \(\ni\)

\[f\left ( x\mid 5,3 \right )=15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\]

Nilai \(E\left ( X \right )\)

\(E\left ( X \right )=x \cdot f\left ( x\mid 5,3 \right )\) \(E\left ( X \right )=x \cdot 15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\) \(\ni\) \[E\left ( X \right )=15 x^{5}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\]

Karena dalam soal diperintahkan untuk menghitung nilai \(E\left ( X \right )\) dengan menggunakan metode trapezoidal maka, \(E\left ( X \right )\) dianggap sebagai fungsi. Maka pada soal dengan metode four-point gaussian quadratur, berarti nilai \(n=4\). Berdasarkan tabel gaussian quadratur nilai-nilai koefisien dan titik gauss adalah sebagai berikut:

n Koefisien Titik Gauss
4 \(C_{1} = 0.3478548\) \(C_{2} = 0.6521452\) \(C_{3} = 0.6521452\) \(C_{4} = 0.3478548\) \(x_{1} = -0.861136311\) \(x_{2} = -0.339981043\) \(x_{3} = 0.339981043\) \(x_{4} = 0.861136311\)

langkah 1 Fungsi \(f\left ( x \right )\) dimana dari soal diminta menggunakan \(E\left ( X \right )\) dan selang integrasinya \(\left ( a,b \right )\) telah diketahui yaitu:

\(E\left ( X \right )=x \cdot 15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}\) dan \(x \in \left [ 0,1 \right ]\)

langkah 2 Sebelum kita menggunakan nilai-nilai koefisien dan titik gauss, perlu dilakukan transformasi sehingga domain integral menjadi \(\left [ -1,1 \right ]\).

\(x=\left ( \frac{t(b-a)+a+b}{2} \right )\) \(x=\left ( \frac{t(1-0)+1+0}{2} \right )\) \(x=\left ( \frac{t+1}{2} \right )\)

dan

\(dx=\frac{b-a}{2}dt\) \(dx=\frac{1-0}{2}dt\) \(dx=\frac{1}{2}dt\)

sehingga, integralnya menjadi: \[\int_{0}^{1}15 x^{5}\left ( 1-x^{\left ( 5 \right )} \right )^{2}dx=\int_{-1}^{1}15 \left ( \frac{t+1}{2} \right )^{5}\left ( 1-\left ( \frac{t+1}{2} \right )^{\left ( 5 \right )} \right )^{2}\frac{1}{2}dt\] langkah 3 Orde \(n\) yang diminta yaitu \(4\) atau \(n=4\) langkah 4 Kemudian, kita menggunakan nilai-nilai koefisien dan titik gauss saat \(n=4\) \[\int_{-1}^{1}15 \left ( \frac{t+1}{2} \right )^{5}\left ( 1-\left ( \frac{t+1}{2} \right )^{\left ( 5 \right )} \right )^{2}\frac{1}{2}dt\approx I\] \[I=C_{1}E\left ( t_{1} \right )+C_{2}E\left ( t_{2} \right )+C_{3}E\left ( t_{3} \right )+C_{4}E\left ( t_{4} \right )\] \(I\)=\(0.3478548E\left ( -0.861136311 \right )+\)\(0.6521452E\left ( -0.339981043 \right )+\)\(0.6521452E\left ( 0.339981043 \right )+\)\(0.3478548E\left ( 0.861136311 \right )\)

Menghitung \(E\left ( X \right )\) pada setiap \(X_{i}\)

\(E\left ( -0.861136311 \right )=\frac{15}{2} \left ( \frac{-0.861136311+1}{2} \right )^{5}\left ( 1-\left ( \frac{-0.861136311+1}{2} \right )^{\left ( 5 \right )} \right )^{2}\) \(E\left ( -0.861136311 \right )=1.21019\cdot10^{-5}\)

\(E\left ( -0.339981043 \right )=\frac{15}{2} \left ( \frac{-0.339981043+1}{2} \right )^{5}\left ( 1-\left ( \frac{-0.339981043+1}{2} \right )^{\left ( 5 \right )} \right )^{2}\) \(E\left ( -0.339981043 \right )=0.029126407\)

\(E\left ( 0.339981043 \right )=\frac{15}{2} \left ( \frac{0.339981043+1}{2} \right )^{5}\left ( 1-\left ( \frac{0.339981043+1}{2} \right )^{\left ( 5 \right )} \right )^{2}\) \(E\left ( 0.339981043 \right )=0.757589241\)

\(E\left ( 0.861136311 \right )=\frac{15}{2} \left ( \frac{0.861136311+1}{2} \right )^{5}\left ( 1-\left ( \frac{0.861136311+1}{2} \right )^{\left ( 5 \right )} \right )^{2}\) \(E\left ( 0.861136311 \right )=0.477908858\)

Sehingga

\(I\)=\(0.3478548 \cdot 1.21019 \cdot 10^{-5}+\)\(0.6521452 \cdot 0.029126407+\)\(0.6521452 \cdot 0.757589241+\)\(0.3478548 \cdot 0.477908858\)

\(I = 0.679299934\)

Jadi,

\[\int_{0}^{1}x \cdot 15 x^{4}\left ( 1-x^{\left ( 5 \right )} \right )^{2}dx\approx I=0.679299934\]

Komparasi Metode

komparasi metode digunakan untuk menentukan hasil metode terbaik dengan dilihat error terkecil dari masing-masing metode terhadap nilai exactnya.

Error
exact <- 0
etrap <- 0.7102273-0.631286879
esim <- 0.768397357-0.7102273
egaqu <- 0.7102273-0.679299934
Error <- rbind(exact, etrap,esim,egaqu)
Error
##             [,1]
## exact 0.00000000
## etrap 0.07894042
## esim  0.05817006
## egaqu 0.03092737
Nilai_EX <- rbind(0.7102273,
                   0.631286879,
                   0.768397357,
                   0.679299934
                     )
Nama_Metode <- c("Nilai Exact",
               "Metode Trapezoida untuk n=4",
               "Metode Simpson untuk n=4",
               "Metode four-point Gauss Quadrature"
)

tibble(hasil <- data.frame(Nama_Metode, Nilai_EX, Error))
## # A tibble: 4 x 3
##   Nama_Metode                        Nilai_EX  Error
##   <chr>                                 <dbl>  <dbl>
## 1 Nilai Exact                           0.710 0     
## 2 Metode Trapezoida untuk n=4           0.631 0.0789
## 3 Metode Simpson untuk n=4              0.768 0.0582
## 4 Metode four-point Gauss Quadrature    0.679 0.0309

Kesimpulan

Berdasarkan perhitungan manual dari distribusi kumaraswamy, metode yang terbaik adalah metode four-point Gauss Quadrature karena memiliki nilai error paling kecil dan mendekati exactnya.

bagian b

Jika diketahui nilai exact dari \(E\left ( X \right )=\frac{b\Gamma \left ( 1+\frac{1}{a} \right )\Gamma \left ( b \right )}{\Gamma \left ( 1+\frac{1}{a}+b \right )}\), metode mana yang paling mendekati untuk menghitung nilai \(E\left ( X \right )\) dengan toleransi \(0.0001\)? Gunakan R!

Jawaban bagian b

Nilai exact

a <- 5
b <- 3
E <- (b*gamma(1+(1/a))*gamma(b))/(gamma(1+(1/a)+b))
E
## [1] 0.7102273

Fungsi

EX <- function(x){
  x*a*b*x^(a-1)*(1-x^(a))^(b-1)
} 
# Menggambar kurva f(x)
curve(EX,0,1)
abline(v=seq(0,1,length.out = 7),col=2,lty=2)

Metode trapezoida untuk \(n=4\)

Menghitung integral menggunakan fungsi oleh user-defined function:

trapezoid <- function(ftn, a, b, n = 100) {
     h <- (b-a)/n
     h
     x.vec <- seq(a, b, by = h)
     f.vec <- sapply(x.vec, EX)     # ftn(x.vec)
     Trap <- h*(f.vec[1]/2 + sum(f.vec[2:n]) + f.vec[n+1]/2)
     return(Trap)
}
trapezoid(EX,0,1,n = 4)
## [1] 0.6312869

Menggunakan R fungsi trapzfun dari package pracma

trapzfun(EX,0,1,maxit =  4)
## $value
## [1] 0.7098069
## 
## $iter
## [1] 4
## 
## $rel.err
## [1] 0.005901199

menggunakan R hingga hasil integral mendekati hasil exact-nya dengan toleransi 0.0001

exact_value=0.7102273
tol <-  0.0001
err <- 1
n = 4

while(err>tol){
  res_trap <- trapezoid(EX,0,1,n = n)
  
  err <- abs(res_trap-exact_value)
  
  cat("n=",n,", result=",res_trap,", error=",err,"\n",sep = "")
  
  n=n+1
  if(n==1000){
    break
  }
  
}
## n=4, result=0.6312869, error=0.07894042
## n=5, result=0.6738123, error=0.03641497
## n=6, result=0.6914928, error=0.01873454
## n=7, result=0.6997126, error=0.01051473
## n=8, result=0.7039057, error=0.006321603
## n=9, result=0.7062116, error=0.004015656
## n=10, result=0.7075597, error=0.002667593
## n=11, result=0.7083885, error=0.001838822
## n=12, result=0.7089199, error=0.001307439
## n=13, result=0.7092729, error=0.0009544148
## n=14, result=0.7095147, error=0.000712649
## n=15, result=0.7096846, error=0.0005426739
## n=16, result=0.7098069, error=0.0004204042
## n=17, result=0.7098966, error=0.0003306616
## n=18, result=0.7099637, error=0.0002636072
## n=19, result=0.7100146, error=0.0002127014
## n=20, result=0.7100538, error=0.0001734994
## n=21, result=0.7100844, error=0.0001429188
## n=22, result=0.7101085, error=0.0001187833
## n=23, result=0.7101278, error=9.953088e-05

Metode Simpson untuk \(n=4\)

Menghitung integral menggunakan fungsi oleh user-defined function:

simpson_n <- function(ftn, a, b, n = 100) {
    n <- max(c(2*(n %/% 2), 4))
    h <- (b-a)/n
    x.vec1 <- seq(a+h, b-h, by = 2*h)       # ganjil
    x.vec2 <- seq(a+2*h, b-2*h, by = 2*h)   # genap
    f.vec1 <- sapply(x.vec1, ftn)   # ganjil
    f.vec2 <- sapply(x.vec2, ftn)   # genap
    S <- h/3*(ftn(a) + ftn(b) + 4*sum(f.vec1) + 2*sum(f.vec2))
    return(S)
}
simpson_n(EX,0,1,n = 4)
## [1] 0.7683974

menggunakan R hingga hasil integral mendekati hasil exact-nya dengan toleransi 0.0001

exact_value=0.7102273
tol <-  0.0001
err <- 1
n = 4

while(err>tol){
  res_simp <- simpson_n(EX,0,1,n = n)
  
  err <- abs(res_simp-exact_value)
  
  cat("n=",n,", result=",res_simp,", error=",err,"\n",sep = "")
  
  n=n+1
  if(n==1000){
    break
  }
  
}
## n=4, result=0.7683974, error=0.05817006
## n=5, result=0.7683974, error=0.05817006
## n=6, result=0.7497082, error=0.03948094
## n=7, result=0.7497082, error=0.03948094
## n=8, result=0.728112, error=0.01788467
## n=9, result=0.728112, error=0.01788467
## n=10, result=0.7188088, error=0.008581532
## n=11, result=0.7188088, error=0.008581532
## n=12, result=0.7147289, error=0.004501593
## n=13, result=0.7147289, error=0.004501593
## n=14, result=0.712782, error=0.002554712
## n=15, result=0.712782, error=0.002554712
## n=16, result=0.711774, error=0.001546662
## n=17, result=0.711774, error=0.001546662
## n=18, result=0.7112144, error=0.0009870756
## n=19, result=0.7112144, error=0.0009870756
## n=20, result=0.7108852, error=0.0006578651
## n=21, result=0.7108852, error=0.0006578651
## n=22, result=0.7106819, error=0.0004545628
## n=23, result=0.7106819, error=0.0004545628
## n=24, result=0.7105511, error=0.0003237819
## n=25, result=0.7105511, error=0.0003237819
## n=26, result=0.710464, error=0.0002366803
## n=27, result=0.710464, error=0.0002366803
## n=28, result=0.7104042, error=0.0001769138
## n=29, result=0.7104042, error=0.0001769138
## n=30, result=0.7103621, error=0.0001348299
## n=31, result=0.7103621, error=0.0001348299
## n=32, result=0.7103318, error=0.0001045198
## n=33, result=0.7103318, error=0.0001045198
## n=34, result=0.7103096, error=8.225045e-05

Metode four-point Gauss Quadrature

Menghitung integral menggunakan fungsi oleh user-defined function tidak dibuat secara lebih umum karena tergantung fungsi dan batas dari integral masing-masing fungsinya.

Hitunglah menggunakan R dengan fungsi gaussLegendre dengan menggunakan domain asal

gL <- gaussLegendre(n = 4,a = 0,1)
Ci <- gL$w # koefisien
Ci
## [1] 0.1739274 0.3260726 0.3260726 0.1739274
xi <- gL$x # gauss point
xi
## [1] 0.06943184 0.33000948 0.66999052 0.93056816
I <- sum(Ci * EX(xi))
I
## [1] 0.6792999

menggunakan R hingga hasil integral mendekati hasil exact-nya dengan toleransi 0.0001

exact_value=0.7102273
tol <-  0.0001
err <- 1
n = 4

while(err>tol){
  
  gL <- gaussLegendre(n = n,a = 0,1)

 Ci <- gL$w # koefisien
 xi <- gL$x # gauss point
  
  res_gl <- sum(Ci * EX(xi))

  err <- abs(res_gl-exact_value)
  
  cat("n=",n,", result=",res_gl,", error=",err,"\n",sep = "")
  
  n=n+1
  if(n==1000){
    break
  }
  
}
## n=4, result=0.6792999, error=0.03092738
## n=5, result=0.707375, error=0.002852279
## n=6, result=0.7101422, error=8.507726e-05

Komparasi Metode

komparasi metode digunakan untuk menentukan hasil metode terbaik dengan dilihat error terkecil dari masing-masing metode terhadap nilai exactnya.

Error
exact <- 0.7102273
etrapUD <- 0.7102273-0.6312869
esimUD <- 0.7683974-0.7102273
egaquUD <- NA
ErrorUserDefine <- rbind(exact, etrapUD,esimUD,egaquUD)
ErrorUserDefine
##              [,1]
## exact   0.7102273
## etrapUD 0.0789404
## esimUD  0.0581701
## egaquUD        NA
exact <- 0.7102273
etrapFR <- 0.7102273-0.7098069
esimFR <- NA
egaquFR <- 0.7102273-0.6792999
ErrorRFunction <- rbind(exact, etrapFR,esimFR,egaquFR)
ErrorRFunction
##              [,1]
## exact   0.7102273
## etrapFR 0.0004204
## esimFR         NA
## egaquFR 0.0309274
exact <- 0.7102273
etrapT <- 0.7102273-0.7101278
esimT <- 0.7103096-0.7102273
egaquT <- 0.7102273-0.7101422
ErrorTolerance <- rbind(exact, etrapT,esimT,egaquT)
ErrorTolerance
##             [,1]
## exact  0.7102273
## etrapT 0.0000995
## esimT  0.0000823
## egaquT 0.0000851
Nama_Metode <- c("Nilai Exact",
               "Metode Trapezoida untuk n=4",
               "Metode Simpson untuk n=4",
               "Metode four-point Gauss Quadrature"
)

tibble(hasil <- data.frame(Nama_Metode, ErrorUserDefine,ErrorRFunction,ErrorTolerance))
## # A tibble: 4 x 4
##   Nama_Metode                        ErrorUserDefine ErrorRFunction ErrorTolerance
##   <chr>                                        <dbl>          <dbl>          <dbl>
## 1 Nilai Exact                                 0.710        0.710         0.710    
## 2 Metode Trapezoida untuk n=4                 0.0789       0.000420      0.0000995
## 3 Metode Simpson untuk n=4                    0.0582      NA             0.0000823
## 4 Metode four-point Gauss Quadrature         NA            0.0309        0.0000851

Kesimpulan

Berdasarkan perhitungan dengan tiga metode dari distribusi Kumaraswamy, metode yang terbaik adalah metode Simpson untuk n=4 karena memiliki nilai error paling kecil dan mendekati exactnya.

Soal No 2

Diketahui pdf (probability density function) dari distribusi Eksponensial adalah sebagai berikut:

\[ f\left ( x\mid a,b \right )=\lambda e^{-\lambda x}\]

dimana \(x\in \left ( 0,\infty \right )\).

Jika diketahui \(\lambda=2\). Hitunglah CDF dari distribusi Eksponensial tersebut untuk \(x=4\), dengan menggunakan metode Trapezoida, Simpson, Gauss Quadrature dan Monte Carlo dengan R. Metode mana yang paling baik? Silahkan tentukan nilai toleransi, n dan m sendiri.

Jawaban soal no 2

Fungsi distribusi kumulatif (CDF) dari suatu distribusi didefinisikan sebagai:

\[ F\left ( x \right )=\int_{0 }^{\infty }f\left ( t \right )dt \]

Sehingga, fungsi distribusi kumulatif (CDF) distribusi Eksponensial untuk \(x=4\), jika diketahui \(\lambda=2\):

\[ F\left ( x=4\mid \lambda =2 \right )=\int_{0}^{4 }2 e^{-2 t} \]

Nilai exact

a <- 0
b <- 4
F <- function(t){
  2*exp(-2*t)
}
Fint <- integrate(F, 0, 4)
Fint
## 0.9996645 with absolute error < 1.1e-14

Fungsi

F <- function(t){
  2*exp(-2*t)
}

Menggambar kurva F

curve(F,0,4)
abline(v=seq(0,4,length.out = 7),col=2,lty=2)

Metode trapezoida

Menghitung integral untuk \(n=3000\) menggunakan fungsi oleh user-defined function:

trapezoid <- function(ftn, a, b, n = 1000) {
     h <- (b-a)/n
     h
     x.vec <- seq(a, b, by = h)
     f.vec <- sapply(x.vec, F)     # ftn(x.vec)
     Trap <- h*(f.vec[1]/2 + sum(f.vec[2:n]) + f.vec[n+1]/2)
     return(Trap)
}
trapezoid(F,0,4,n = 3000)
## [1] 0.9996651

Menggunakan R fungsi trapzfun dari package pracma

trapzfun(F,0,4,maxit =  3000)
## $value
## [1] 0.9996646
## 
## $iter
## [1] 14
## 
## $rel.err
## [1] 5.958465e-08

Menggunakan R hingga hasil integral mendekati hasil exact-nya dengan toleransi 0.00001

exact_value=0.9996645
tol <-  0.00001
err <- 1
n = 3000

while(err>tol){
  res_trap <- trapezoid(F,0,4,n = n)
  
  err <- abs(res_trap-exact_value)
  
  cat("n=",n,", result=",res_trap,", error=",err,"\n",sep = "")
  
  n=n+1
  if(n==1000){
    break
  }
  
}
## n=3000, result=0.9996651, error=6.297658e-07

Metode Simpson

Menghitung integral untuk \(n=500\) menggunakan fungsi oleh user-defined function:

simpson_n <- function(ftn, a, b, n = 100) {
    n <- max(c(2*(n %/% 2), 4))
    h <- (b-a)/n
    x.vec1 <- seq(a+h, b-h, by = 2*h)       # ganjil
    x.vec2 <- seq(a+2*h, b-2*h, by = 2*h)   # genap
    f.vec1 <- sapply(x.vec1, ftn)   # ganjil
    f.vec2 <- sapply(x.vec2, ftn)   # genap
    S <- h/3*(ftn(a) + ftn(b) + 4*sum(f.vec1) + 2*sum(f.vec2))
    return(S)
}
simpson_n(F,0,4,n = 500)
## [1] 0.9996645

Menggunakan R hingga hasil integral mendekati hasil exact-nya dengan toleransi 0.0001

exact_value=0.9996645
tol <-  0.0001
err <- 1
n = 500

while(err>tol){
  res_simp <- simpson_n(F,0,4,n = n)
  
  err <- abs(res_simp-exact_value)
  
  cat("n=",n,", result=",res_simp,", error=",err,"\n",sep = "")
  
  n=n+1
  if(n==1000){
    break
  }
  
}
## n=500, result=0.9996645, error=3.773605e-08

Metode Gauss Quadrature

Menggunakan fungsi gaussLegendre dari package pracma R dengan \(n=7\)

gL <- gaussLegendre(n = 7,a = 0,4)
Ci <- gL$w # koefisien
xi <- gL$x # gauss point
I <- sum(Ci * F(xi))
I
## [1] 0.9996645

Menggunakan R hingga hasil integral mendekati hasil exact-nya dengan toleransi 0.00001

exact_value=0.9996645
tol <-  0.00001
err <- 1
n = 7

while(err>tol){
  
  gL <- gaussLegendre(n = n,a = 0,4)

 Ci <- gL$w # koefisien
 xi <- gL$x # gauss point
  
  res_gl <- sum(Ci * F(xi))

  err <- abs(res_gl-exact_value)
  
  cat("n=",n,", result=",res_gl,", error=",err,"\n",sep = "")
  
  n=n+1
  if(n==1000){
    break
  }
  
}
## n=7, result=0.9996645, error=1.619311e-08

Metode Monte Carlo

CDF dari suatu distribusi Eksponensial didefinisikan sebagai berikut \[F\left ( x=4\mid \lambda =2 \right )=\int_{0}^{4 }2 e^{-2 t}dt\]

x <- 4
lambda <- 2
F <- function(t){
  lambda * exp(-lambda*t)
}
mc_integral <- function(ftn, a, b,m=1000){
  #Membangkitkan x berdistribusi U(a,b)
  x <- runif(m,a,b)
  # Menghitung rata-rata dari output fungsi
  Gx <- ftn(x)
  Gx_mean <- mean(Gx)
  theta.hat <- (b-a)*Gx_mean
  return(theta.hat)
}
set.seed(1)
int_g <- mc_integral(F,0, 4,m = 1256)
int_g
## [1] 0.99913

Komparasi Metode

komparasi metode digunakan untuk menentukan hasil metode terbaik dengan dilihat error terkecil dari masing-masing metode terhadap nilai exactnya.

Error
exact <- 0.9996645
etrapUD <- 0.9996651-0.9996645
esimUD <- 0.9996645-0.9996645
egaquUD <- NA
emcUD <- 0.9996645-0.99913
ErrorUserDefine <- rbind(exact, etrapUD,esimUD,egaquUD, emcUD)
exact <- 0.9996645
etrapFR <- 0.9996646-0.9996645
esimFR <- NA
egaquFR <- 0.9996645-0.9996645
emcFR <- NA
ErrorRFunction <- rbind(exact, etrapFR,esimFR,egaquFR, emcFR)
exact <- 0.9996645
etrapT <- 0.9996651-0.9996645
esimT <- 0.9996645-0.9996645
egaquT <- 0.9996645-0.9996645
emcT <- NA
ErrorTolerance <- rbind(exact, etrapT,esimT,egaquT, emcT)
Nama_Metode <- c("Nilai Exact",
               "Metode Trapezoida n=3000",
               "Metode Simpson n=500",
               "Metode Gauss Quadrature n=7",
               "Metode Monte Carlo m=1256"
)

tibble(hasil <- data.frame(Nama_Metode, ErrorUserDefine,ErrorRFunction,ErrorTolerance))
## # A tibble: 5 x 4
##   Nama_Metode                 ErrorUserDefine ErrorRFunction ErrorTolerance
##   <chr>                                 <dbl>          <dbl>          <dbl>
## 1 Nilai Exact                     1.00           1.00           1.00       
## 2 Metode Trapezoida n=3000        0.000000600    0.000000100    0.000000600
## 3 Metode Simpson n=500            0             NA              0          
## 4 Metode Gauss Quadrature n=7    NA              0              0          
## 5 Metode Monte Carlo m=1256       0.000534      NA             NA

Kesimpulan

Berdasarkan perhitungan dengan empat metode dari distribusi Kumaraswamy, metode yang terbaik adalah metode Gauss Quadrature n=7, toleransi=0.00001 karena memiliki nilai error paling kecil, nilai \(n\) yang terkecil, dan mendekati exactnya.

TERIMAKASIH