Simulasi Monte Carlo

Monte Carlo untuk menghitung integral

N=10000
set.seed(1)
x2<-rep(NA,N)
for(i in 1:N){
  x<- runif(N)
  x2[i]<-sum(x)/N
}

plot(x2)

N=10000
set.seed(1)
x2<-rep(NA,N)
for(i in 1:N){
  x<- runif(N)
  x2[i]<-sum(x^2)/N
}

plot(x2)

Simuasi menghitung peluang banyaknya klaim

Model berhitung Poisson

set.seed(1)
t=2
N1= rpois(1000,t*13)

hist(N1, main="Histogram distribusi Poisson(13t)", xlab="Banyaknya klaim", probability = TRUE)

lines(density(N1))

Misalkan \(N_i\) banyaknya klaim yang dilaporkan oleh pemegang polis ke-\(i\) dengan karakter yang terobservasi \(x_i\) dan eksposur risikonya \(t_i\). Kita asumsikan diberikan \(\theta_i\), \[N_i | \theta \sim Poisson(t_i\lambda_i \theta_i),\] dengan \(\lambda_i=e^{x_i^T \beta}\).

Kita memiliki \[P(N_i =k)=\int_0^\infty P(N_i=k|\theta_i)h(\theta_i)d\theta_i.\] Misalkan kita asumsikan \(\theta_i \sim Gamma(1/2,1/2)\). Maka kita dapat menggunakan Simulasi Monte Carlo untuk menghitung \(P(N_i=k)\) di atas.

set.seed(2)
N=1000
theta<- rgamma(N,1/2,1/2)
lambda=13
t=2
Po=rep(NA,21)
pm<- t*lambda*theta
for(k in 0:20){
  Po[k]<- sum( dpois(k,pm))/N
}

plot(Po)

Prediksi N

Diberikan \(n_{i,1}, n_{i,2}, ..., n_{i,m}\). Kita ingin menentukan \(P(N_{i,m+1}=k|N_{i,1}=n_{i,1}, ...,N_{i,m}=n_{i,m})\).

Kita dapat tuliskan \[P(N_{i,m+1}=k|N_{i,1}=n_{i,1}, ...,N_{i,m}=n_{i,m})=\int P(N_{i,m+1}=k,\theta_i|N_{i,1}=n_{i,1}, ...,N_{i,m}=n_{i,m})d\theta_i\] \[=\int P(N_{i,m+1}=n_{i,m+1}|\theta_i, n_{i,j},j=1,...,m)P(\theta_i|n_{i,j},j=1,2,...,m)d\theta_i\] \[= \int P(N_{i,m+1}=n_{i,m+1}|\theta_i, n_{i,j},j=1,...,m)\frac{P(n_{i,j},j=1,2,...,m|\theta_i)f_H(\theta_i)}{P(N_{i,j}=n_{i,j}, j=1,2,...,m)}d\theta_i\] \[=\int \prod_{k=1}^{m+1}\frac{P(N_{i,k}=n_{i,k}|\theta_i)f_H(\theta_i)}{P(N_{i,j}=n_{i,j}, j=1,2,...,m)} d\theta_i.\] ## PR lakukan simulasi untuk distribusi di atas