Teori Risiko dan Simulasi
Sapto Indratno
August 30, 2018
Simulasi pembentukan banyak klaim dan waktu eksposur pemegang polis yang loyal
Diberikan \(n\) pemegang polis dalam periode waktu \(D1\) sampai \(D2\).
n=983
D1=as.Date("01/01/1993",'%d/%m/%Y')
D2=as.Date("31/12/2013",'%d/%m/%Y')
# Tanggal-tanggal yang ada dalam kurun waktu [D1,D2]
L=D1+0:(D2-D1)
# Tetapkan himpunan random
set.seed(1)
# Ambil sampel sebanyak n dengan pengembalian. Sampel ini menggambarkan waktu kedatangan klaim pertama dalam periode [D1,D2].
arrival=sample(L,size=n,replace=TRUE)
# Pesan vektor exposure dan depature
exposure=N=rep(NA,n)
depature=rep(D2,n)
set.seed(2)
# Pembentukkan kedatangan klaim beserta exposurenya untuk pemegang polis yang loyal
for(i in 1:n){
expo=D2-arrival[i]
w=0
while(max(w)<expo){
w=c(w,max(w)+1+trunc(rexp(1,1/1000)))} # trunc untuk membulatkan waktu antar klaim harapannya 1000 hari
exposure[i]=depature[i]-arrival[i]
N[i]=max(0,length(w)-2) # banyaknya klaim orang ke i, kita kurangi 2 untuk menghilangkan klaim pertama dan klaim akhir yang di luar D2.
}
# Tabel banyaknya Klaim dan eksposur yang dihasilkan
df=data.frame(N=N,E=exposure/365)
head(df)## N E
## 1 6 15.432877
## 2 4 13.191781
## 3 2 8.975342
## 4 0 1.928767
## 5 2 16.775342
## 6 1 2.134247i1=which(df$N>0) # kumpulan pemegang polis yang klaimnya >0
i0=which(df$N==0) # kumpulan pemegang polis yang klaimnya=0
h1=hist(df$E[i1],probability = TRUE, main="Histogram Eksposur dengan Klaim >0")
h0=hist(df$E[i0],probability = TRUE,main="Histogram Eksposur dengan Klaim =0")
plot(h1$mids,h1$density,type='s',lwd=2,col='red')
lines(h0$mids,h0$density,type='s',lwd=2,col='blue')
Simulasi untuk pemegang polis yang tidak loyal
n=983
D1=as.Date("01/01/1993",'%d/%m/%Y')
D2=as.Date("31/12/2013",'%d/%m/%Y')
# Tanggal-tanggal yang ada dalam kurun waktu [D1,D2]
L=D1+0:(D2-D1)
# Tetapkan himpunan random
set.seed(1)
# Ambil sampel sebanyak n dengan pengembalian. Sampel ini menggambarkan waktu kedatangan klaim pertama dalam periode [D1,D2].
arrival=sample(L,size=n,replace=TRUE)
# Pesan vektor exposure dan depature
exposure=N=rep(NA,n)
depature=rep(D2,n)
set.seed(2)
# Pembentukkan kedatangan klaim beserta exposurenya untuk pemegang polis yang loyal
for(i in 1:n){
expo=D2-arrival[i]
w=c(0,0)
while((max(w)<expo)&(max(diff(w))<1500)){
w=c(w,max(w)+1+trunc(rexp(1,1/1000)))} # trunc untuk membulatkan waktu antar klaim harapannya 1000 hari
if(max(diff(w))>1500){
depature[i]=arrival[i]+max(w[-length(w)])+1500
}
exposure[i]=depature[i]-arrival[i]
N[i]=max(0,length(w)-3) # banyaknya klaim orang ke i, kita kurangi 2 untuk menghilangkan klaim pertama dan klaim akhir yang di luar D2.
}
# Tabel banyaknya Klaim dan eksposur yang dihasilkan
df=data.frame(N=N,E=exposure/365)
head(df)## N E
## 1 0 4.109589
## 2 2 5.621918
## 3 3 9.128767
## 4 0 1.928767
## 5 4 11.838356
## 6 0 2.134247i1=which(df$N>0) # kumpulan pemegang polis yang klaimnya >0
i0=which(df$N==0) # kumpulan pemegang polis yang klaimnya=0
h1=hist(df$E[i1],probability = TRUE, main="Histogram Eksposur dengan Klaim >0")
h0=hist(df$E[i0],probability = TRUE,main="Histogram Eksposur dengan Klaim =0")
plot(h1$mids,h1$density,type='s',lwd=2,col='red')
lines(h0$mids,h0$density,type='s',lwd=2,col='blue')
Simulasi pemegang polis yang keluar karena pelayanan
Kadang pemegang polis tidak puas dengan pelayanan pembayaran klaim. Misalkan setelah proses pembayaran klaim pertama terdapat peluang 50% pemegang polis keluar dari asuransi ini. Tetapi bisa juga dikarenakan klaim pertama ini menyebabkan mobilnya rusak sehingga tidak dapat dikendarai lagi, sehingga tidak ada keuntungan membeli asuransi bermotor lagi.
n=983
D1=as.Date("01/01/1993",'%d/%m/%Y')
D2=as.Date("31/12/2013",'%d/%m/%Y')
# Tanggal-tanggal yang ada dalam kurun waktu [D1,D2]
L=D1+0:(D2-D1)
# Tetapkan himpunan random
set.seed(1)
# Ambil sampel sebanyak n dengan pengembalian. Sampel ini menggambarkan waktu kedatangan klaim pertama dalam periode [D1,D2].
arrival=sample(L,size=n,replace=TRUE)
# Pesan vektor exposure dan depature
exposure=N=rep(NA,n)
depature=rep(D2,n)
set.seed(2)
# Pembentukkan kedatangan klaim beserta exposurenya untuk pemegang polis yang loyal
for(i in 1:n){
expo=D2-arrival[i]
w=0
stay=TRUE
while((max(w)<expo)&(stay==TRUE)){
w=c(w,max(w)+1+trunc(rexp(1,1/1000))) # trunc untuk membulatkan waktu antar klaim harapannya 1000 hari
stay=sample(c(TRUE,FALSE),prob=c(0.5,0.5),size=1)}
N[i]=length(w)-2
if(stay==FALSE){
depature[i]=arrival[i]+max(w)
}
exposure[i]=depature[i]-arrival[i]
N[i]=max(0,length(w)-2) # banyaknya klaim orang ke i, kita kurangi 2 untuk menghilangkan klaim pertama dan klaim akhir yang di luar D2.
}
# Tabel banyaknya Klaim dan eksposur yang dihasilkan
df=data.frame(N=N,E=exposure/365)
head(df)## N E
## 1 1 5.515068
## 2 1 3.191781
## 3 0 1.430137
## 4 0 1.928767
## 5 2 7.909589
## 6 0 2.183562i1=which(df$N>0) # kumpulan pemegang polis yang klaimnya >0
i0=which(df$N==0) # kumpulan pemegang polis yang klaimnya=0
h1=hist(df$E[i1],probability = TRUE, main="Histogram Eksposur dengan Klaim >0")
h0=hist(df$E[i0],probability = TRUE,main="Histogram Eksposur dengan Klaim =0")
plot(h1$mids,h1$density,type='s',lwd=2,col='red')
lines(h0$mids,h0$density,type='s',lwd=2,col='blue')
Pemodelan dengan data sebenarnya
reg=glm(N~1+offset(log(E)), data=df, family=poisson)
summary(reg)##
## Call:
## glm(formula = N ~ 1 + offset(log(E)), family = poisson, data = df)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.5003 -0.8623 -0.4449 0.3851 3.4103
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.73600 0.03594 -48.3 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 924.93 on 982 degrees of freedom
## Residual deviance: 924.93 on 982 degrees of freedom
## AIC: 1956.7
##
## Number of Fisher Scoring iterations: 5data1<- read.csv("datapolismotor.csv")
str(data1)## 'data.frame': 413169 obs. of 11 variables:
## $ X : int 1 2 3 4 5 6 7 8 9 10 ...
## $ PolicyID : int 1 2 3 4 5 6 7 8 9 10 ...
## $ ClaimNb : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Exposure : num 0.09 0.84 0.52 0.45 0.15 0.75 0.81 0.05 0.76 0.34 ...
## $ Power : Factor w/ 12 levels "d","e","f","g",..: 4 4 3 3 4 4 1 1 1 6 ...
## $ CarAge : int 0 0 2 2 0 0 1 0 9 0 ...
## $ DriverAge: int 46 46 38 38 41 41 27 27 23 44 ...
## $ Brand : Factor w/ 7 levels "Fiat","Japanese (except Nissan) or Korean",..: 2 2 2 2 2 2 2 2 1 2 ...
## $ Gas : Factor w/ 2 levels "Diesel","Regular": 1 1 2 2 1 1 2 2 2 2 ...
## $ Region : Factor w/ 10 levels "R11","R23","R24",..: 9 9 5 5 6 6 9 9 5 1 ...
## $ Density : int 76 76 3003 3003 60 60 695 695 7887 27000 ...i1=which(data1$ClaimNb>0)
i0=which(data1$ClaimNb==0)
h1=hist(data1$Exposure[i1],probability = TRUE)
h0=hist(data1$Exposure[i0],probability = TRUE)
plot(h1$mids,h1$density,type='s',lwd=2,col='red')
lines(h0$mids,h0$density,type='s',lwd=2,col='blue')