A k-Inflated Negative Binomial Mixture Regression Model: Application to Rate--Making Systems

0.32 MB·

by  Amir T. Payandeh Najafabadi

A k-Inflated Negative Binomial Mixture Regression Model: Application to
Rate–Making Systems
January 20, 2017
Amir T. Payandeh Najafabadia,1 & Saeed MohammadPourb
a Department of Mathematical Sciences, Shahid Beheshti University, G.C. Evin, 1983963113,
Tehran, Iran.
7
1
b E.C.O College of Insurance, Allameh Tabatabái University, Tehran, Iran.
0
2

n
a Abstract
J

9
1
This article introduces a k-Inflated Negative Binomial mixture distribution/regression model as a

]
E
more flexible alternative to zero-inflated Poisson distribution/regression model. An EM algorithm
M
. has been employed to estimate the model’s parameters. Then, such new model along with a Pareto
t
a
t
s mixture model have been employed to design an optimal rate–making system. Namely, this article
[

1 employs number/size of reported claims of Iranian third party insurance dataset. Then, it employs
v
2 the k-Inflated Negative Binomial mixture distribution/regression model as well as other well devel-
5
4
oped counting models along with a Pareto mixture model to model frequency/severity of reported
5
0
. claims in Iranian third party insurance dataset. Such numerical illustration shows that: (1) the k-
1
0
7 Inflated Negative Binomial mixture models provide more fair rate/pure premiums for policyholders
1
:
v under a rate–making system; and (2) in the situation that number of reported claims uniformly
i
X
distributed in past experience of a policyholder (for instance k = 1 and k = 1 instead of k = 0
r 1 2 1
a
and k = 2). The rate/pure premium under the k-Inflated Negative Binomial mixture models are
2
more appealing and acceptable.
Keywords: Negative Binomial regression; Poisson regression; Mixture model; Overdispersed be-
havior; Heavy–tail behavior; Inflated model; EM algorithm; Rate–making System.
1
Corresponding author: [email protected]
1
2010 Mathematics Subject Classification: 62Jxx, 91B30, 97M30
1. Introduction
Modeling count data is an interesting topic in variety fields of applied sciences, such as actuarial
sciences, economics, sociology, engineering, etc. In many practical situation the popular classical
Poisson regression model fails to model count data which exhibit overdispersion (i.e., the variance of
the response variable exceeds its mean). Moreover, strict assumptions on Poisson distribution make
it more less applicable in situation that such assumption cannot be strictly verified. The Negative
Binomial distribution/regression has become more and more popular as a more flexible alternative
to Poisson distribution/regression. In a situation that strict requirements for Poisson distribution
cannot be verified, the Negative Binomial distribution is an appropriate choice (Johnson et al.,
2005). Moreover, the Negative Binomial is an appropriate choice for overdispersed count data that
are not necessarily heavy–tailed (Aryuyuen & Bodhisuwan, 2013).
For count data, the overdispersed behavior has been arrived by either observing excess of a single
value more than number of expected under the model or the target population consisting of several
sub-populations. Using k-Inflated and mixture models are two popular statistical approach to deal-
ing with an overdispersed behavior. Simar (1976) and Laird (1978) were two authors who employed
Poisson mixture models to considering an overdispersed behavior. Lambert (1992) considered a
zero-inflated Poisson regression model to take into account an overdispersed behavior. Wedel et
al. (1993) Brännaäs & Rosenqvist (1994), Wang et al. (1996), Alfó & Trovato (2004), Wang et al.
(1998), among others, developed idea of using a finite mixture Poisson regression model to handel
overdispersion.
Greene (1994) and Hall (2000) were pioneer authors who employed zero-inflated Negative Binomial
regression to model overdispersion. The ordinary Negative Binomial distribution can be viewed as
a mixture of Poisson and gamma distributions (Simon, 1961). To handel an overdispersion phe-
2
nomena, several extension of Negative Binomial distribution have been introduced by authors. For
instance Negative Binomial exponential distribution (Panjer & Willmot, 1981), Negative Binomial
Pareto distribution (Meng et al., 1999), Negative Binomial Inverse Gaussian distribution (Gómez-
Déniz et al., 2008), Negative Binomial Lindley distribution (Zamani & Ismail, 2010), Negative
Binomial Beta Exponential distribution (Pudprommarat, 2012), and Negative Binomial General-
ized Exponential distribution (Pudprommarat et al., 2012).
In 2014, Lim et al. considered a k-Inflated Poisson mixture model which simultaneously takes into
account both inflated and mixture approaches to handel an overdispersion phenomena. Moreover
Tzougas et al. (2014)’ introduced a Negative Binomial mixture model to model an overdispersion
phenomena. This article follows Lim et al. (2014)’sand generalized Tzougas et al. (2014)’sfindings.
More precisely, It introduces a k-Inflated Negative Binomial mixture distribution/regression. To
show practical application of our finding, we consider the problem of designing an optimal rate–
making system. Then, premium of such optimal rate–making system has been evaluated using the
result of this article.
This article has been structured as follows. The k-Inflated Negative Binomial mixture model,
some of its properties, and an EM algorithm, to estimated its parameters, have been developed
in Section 2. The Pareto mixture regression model has been given in Section 3. Application of
the k-Inflated Negative Binomial mixture model along with a Pareto mixture model to design an
optimal rate–making system have been given in Section 4. Section 5 employs our model, as well as
other well-known model, to evaluate rate, base, and pure premiums under a rate–making system for
Iranian third party insurance dataset. Base upon three comparison methods, Section 6 shows that
our model provides more accurate (in some sense) results. Section 7 employs our well fitted models
to calculate rate and pure premiums under two different Scenarios. Conclusion and suggestions
have been given in Section 8
3
2. k-Inflated Negative Binomial mixture regression model
Thek-InflatedNegativeBinomialmixture,saykINBM,distributionarrivesbycombiningmweighted
mixture Negative Binomialdistribution witha single mass atpoint k.The probability mass function
for a kINBM distribution has been given by
m y+α −1 τ αj α y
P(Y =y|θ) = p1I(y=k)+ pj yj α +j τ α +jτ IN(y), (1)
j=2 (cid:18) (cid:19)(cid:18) j j(cid:19) (cid:18) j j(cid:19)
X
wherek ∈ Nandθ standsforall3munknownparameters. Moreover, m p = 1andp ,α ,τ ≥ 0,
j=1 j j j j
for all j = 1,··· ,m. By a straightforward calculation, one may showPthat
m
E(Y) = p k + p α2/τ
1 j j j
j=2
X
m τ αj α
M (t) = p ekt + p j ,t ≤ −max log( j ),j = 2,··· ,m
Y 1 j τ +α (1−et) α +τ
j=2 (cid:18) j j (cid:19) (cid:26) j j (cid:27)
X
m
F (r) = p I (r)+1−p + p RIBeta (r +1;α ),
Y 1 [k,∞) 1 j αj/(αj+τj) j
j=2
X
where RIBeta (a;b) = xta−1(1−t)b−1dtΓ(a+b)/(Γ(a)Γ(b)), for a,b ≥ 0,x ∈ [0,1], stands for the
x 0
regularized incomplete bReta function.
It is well-known that a Negative Binomial distribution can be arrived by mixing two Poisson and
gamma distributions (Simon, 1961). The following generalized the above fact to the kINBM distri-
bution.
Corollary 1. Suppose random variable Y, given parameter Λ = λ, has been distributed accord-
ing to a k-Inflated Poisson distribution with probability mass function P(Y = y|λ) = pI (y) +
{k}
qexp(−λ )(λ )y/y!, where p & q ∈ [0,1] and p+q = 1. Moreover, suppose that parameter λ has been
i i
distributed according to a finite mixture gamma distribution f (λ) = m ϕ λαj−1τ αje−τjλ/Γ(α ),
Λ j=1 j j j
where, for all j = 1,··· ,m, ϕ ∈ [0,1] and m ϕ = 1. Then, unconPditional distribution of Y has
j j=1 j
a kINB finite mixture distribution with probPability mass function
m y +αj −1 τj αj 1 y
P (Y = y) = pI (y)+q ϕ ( ) ( ) . (2)
{k} j
y 1+τ 1+τ
j=1 (cid:18) (cid:19) j j
X
4
For practical application, in Equation (2), we set q := m ω and ϕ := ω / m ω .
s=1 s j j s=1 s
P P
th
Now, to formulated a kINBM regression model, suppose that for an i individual, information on
countresponsevariablesY ,··· ,Y alongwithinformationonpcovariatesX ,··· ,X areavailable.
i1 it 1 p
Also suppose that Y given parameter Λ = λ has been distributed according to a k-Inflated Pois-
il il il
son distribution with probability mass function P(Y = y |λ ) = pI (y )+qexp(−λ )(λ )y/y !,
il il il {k} il il il il il
where p & q ∈ [0,1] and p+q = 1. Moreover, suppose that parameter λ can be evaluated by the
il
following regression model
p
log(λ ) = β + β x +ǫ ,
il 0il kil ki i
k=1
X
where β ,··· ,β are regression coefficients and u = exp(ǫ ) has been distributed according to a
0i pi i i
finite mixture gamma distribution with density function
m u αj−1α αje−αjui
i j
f (u ) = ϕ , (3)
Ui i j Γ(α )
j
j=2
X
where m ϕ = 1 and α ≥ 0. To have E(ǫ ) = 0, we set both parameters of all gamma
j=2 j j i
distribuPtions, in the finite mixture gamma distribution, to be equal.
Using the law of total probability and setting d := β + p β x , one may show that
il 0il k=1 kil ki
P
∞
P(Y = y |θ) = P(Y = y |θ,u )f (u )du
il il il il i Ui i i
Z0
m dyilααj ∞
= pI (y )+ qϕ ijl j e−(dijl+αj)uiuyil+αj−1du
{k} il jy !Γ(α ) i i
j=2 il j Z0
X
m Γ(y +α ) dyilααj
il j ijl j
= pI (y )+ qϕ ,
{k} il j yil!Γ(αj) (dijl +αj)yil+αj
j=2
X
where θ stands for all unknown parameters. Now by setting p = 1/(1+ m eωs), qϕ = eωj/(1+
s=2 j
m eωs),forj = 2,...,m,and yil+αj−1 := Γ(yil+αj),thekINBMregresPsionmodelcanberestated
s=2 yil yil!Γ(αj)
(cid:16) (cid:17)
P
5
as
1
P (Y = y ) = I (y ) (4)
il il m {k} il
1+ eωj
j=2
X
m eωj y +α −1
+ il j t αj(1−t )yil,,
m y ilj ilj
Xj=2 1+ eωs (cid:18) il (cid:19)
s=2
X
where t := α /(α +exp{X B }) and X B := β + p β x .
ilj j j i il i il 0il k=1 kil ik
P
Parameters estimation
All unknown parameters of the kINBM regression (4) can be represented as θ := (ω,α,B). Now
to provide a Maximum likelihood estimator, say MLE, for θ, one may employ an EM algorithm.
In statistical literature, the EM algorithm is a well-known and practical method to obtain the
Maximum likelihood estimators for parameters in an arbitrary finite mixture model (McLachlan &
Krishnan, 1997). Now suppose that number of components, m, is given, and v = (v ,··· ,v )
i i1 im
stands for the latent vector of component indicator variables, where for i = 1,··· ,n and j =
th
1,··· ,m, v = 1 whenever observation i comes from j component and v = 0, otherwise.
ij ij
Therefore, we assume that each observation has been arrived from one of the m components, but
the component it belongs to is unobservable and therefore considered to be the missing data.
Now using the Multinomial distribution for the unobservable vector v , the complete data log-
i
likelihood function, for, the kINB regression model, can be written as the following, see Rigby &
Stasinopoulos (2009) for an update information.
n
1
l (θ|y ,v ,X ) = v log I (5)
c i i i i1 1+ m exp(ω ) (yi=k)
i=1 (cid:18) l=2 l (cid:19)
X
n m
P exp(ω ) y +α −1
+ vijlog 1+ m exjp(ω ) i yj tjαj(1−tj)yi ,
i=1j=2 (cid:18) l=2 l (cid:18) i (cid:19) (cid:19)
XX
P
where θ := (ω,α,B) stands for all unknown parameters, t := α /(α +exp{X B }), and X B :=
j j j i j i j
β + p β x .
0j k=1 kj ik
P
6
The EM algorithm employs the following two steps to maximize the above loglikelihood function.
E-step: Inthisstep,usinggivendataalongwithcurrentestimatesθ(r) := (ω(r),α(r),B(r))obtained
th th
from the r iteration, the probability vˆ estimates. This probability, (r+1) iteration, can
ij
b b b b
be stated as
vˆ(r+1) = E v |y ,X ,θ(r) =1×P v =1|y ,X ,θ(r) +0×P v =0|y ,X ,θ(r)
ij ij i i ij i i ij i i
(cid:20) (cid:21) (cid:18) (cid:19) (cid:18) (cid:19)
b (r) (r)b b
f y |v =1,X ,θ P v =1|X ,θ
i in i in i
= (cid:18) (cid:19) (cid:18) (cid:19)
f by |X ,θ(r) b
i i
(cid:18) (cid:19)
= exp ωˆj(r) yi+αby(jir)−b1 tj(r) αb(jr) 1−tj(r) yi ,
1+ ms=2(cid:16)exp(cid:17)ωˆ(cid:16)s(r) yi+αby(cid:17)(sir)(cid:16)−b1 (cid:17)ts(r)(cid:16)αb(sr)b1−(cid:17)ts(r) yi
(cid:16) (cid:17)(cid:16) (cid:17)(cid:16) (cid:17) (cid:16) (cid:17)
P
where t(r) := α(r)/(α(r) +exp{X B(r)}), and X Bb(r) := β(r)b+ p β(r)x .
j j j i j i j 0j k=1 kj ik
P th
M-step: Gibven thebprobbability vˆij, thibs step maximibzes, in bthe (r + 1) b iteration, the following
loglikelihood Q(·) with respect to θ := (ω,α,B).
(r)
Q = E l |y ,X ,θ
c i i
(cid:20) (cid:21)
n b n m exp ω(r)
= vˆ(r+1)log 1 I + vˆ(r+1)log l
Xi=1 i1 (cid:18)1+ ml=2exp(ωl)(cid:19) (yi=k) Xi=1Xl=2 il 1+ ms=2(cid:16)ebxp (cid:17)ωs(r)  
n m P y +α(r)−1 α(r) αb(lr) Pexp(X B(cid:16)b(r))(cid:17) yi
+ vˆ(r+1)log i l l i l
Xi=1Xl=2 il  byi !αl(r)+ebxp(XiBl(r)) αl(r)+exp(bXiBl(r)) 
=: Q +Q .      
1 2
b b b b
Updated parameters θ(r+1) := (ω(r+1),α(r+1),B(r+1)) have been arrived by solving the follow-
ing equation.
b b b b
n n
∂Q ∂ 1 ∂ exp(ω )
1 = vˆ(r+1)log I + vˆ(r+1)log j =0;
∂ω ∂ω i1 1+ m exp(ω ) (yi=k) ∂ω ij 1+ m exp(ω )
j i=1 j (cid:18) l=2 l (cid:19) i=1 j (cid:18) l=2 l (cid:19)
X X
∂Q n ∂ y P+α −1 α αj exp(X B )P yi
2 = vˆ(r+1)log i j j i j =0;
∂B ∂B in y α +exp(X B ) α +exp(X B )
j i=1 j (cid:18)(cid:18) i (cid:19)(cid:18) j i j (cid:19) (cid:18) j i j (cid:19) (cid:19)
X
∂Q n ∂ y +α −1 α αj exp(X B ) yi
2 = vˆ(r+1)log i j j i j =0.
∂α ∂α in y α +exp(X B ) α +exp(X B )
j i=1 j (cid:18)(cid:18) i (cid:19)(cid:18) j i j (cid:19) (cid:18) j i j (cid:19) (cid:19)
X
Since the above three equations cannot solve explicitly, such updated parameters have been
7
obtained using the following Iteratively Reweighted Least Squares, say IRLS, method.
−∂2Q −1 ∂Q
(r+1) (r) 1 1
ωˆ = ωˆ + E . ;
j j ∂ω 2 ∂ω
(cid:18) (cid:18) j (cid:19)(cid:19) j
−∂2Q −1 ∂Q
Bˆ(r+1) = Bˆ(r) + E 2 . 2;
j j ∂B 2 ∂B
(cid:18) (cid:18) j (cid:19)(cid:19) j
−∂2Q −1 ∂Q
(r+1) (r) 2 2
α = α + E . .
j j ∂α 2 ∂α
(cid:18) (cid:18) j (cid:19)(cid:19) j
In IRLS method, E(−∂2b(loglikelihoobd )/∂parameter2) can be viewed as the Fisher informa-
tion matrix and ∂(loglikelihood )/∂parameter as score function.
After updated Parameter estimates θ(r+1) := (ω(r+1),α(r+1),B(r+1)), the complete data loglikeli-
th
hood for (r +1) iteration, arrives by
b b b b
n
1
l(r+1) = vˆ(r+1)log I
c i1 1+ m exp ωˆ(r+1)  (yi=k)
Xi=1 l=2 l
+ n m vˆ(r+1)log P exp(cid:16) ωˆj(r+1(cid:17))  yi+α(jr+1)−1 t(r+1) αb(jr+1) 1−t(r+1) yi ,
ij 1+ m (cid:16)exp ωˆ(cid:17)(r+1) yi ! j j 
Xi=1Xj=2 l=2 l b (cid:16) (cid:17) (cid:16) (cid:17)
 (cid:16) (cid:17) b b 
where t(r+1) := α(r+1)/(α(r+1) +expP{X B(r+1)}), and X B(r+1) := β(r+1) + p β(r+1)x . Now, in
j j j i j i j 0j k=1 kj ik
the E-step v −s have been estimated. This loop has been repeated until thePdifference l(r+1) −l(r)
b ij b b b b b b c c
(cid:12) (cid:12)
has been converged, in some sense. (cid:12) (cid:12)
(cid:12) (cid:12)
It is worthwhile to mention that, since regression coefficients B = (β ,··· ,β )′ have been estimated
0 p
using the MLE methods. therefore, number of mixture component impact on such estimators.
3. Pareto mixture regression model
The Pareto mixture distribution arrives by combining m weighted mixture Pareto distributions.
The density function for a Pareto mixture distribution has been given by
m γαj
f (z|ϑ) = ρ α j , (6)
Z j j(z+γj)αj+1
j=1
X
where ϑ = (ρ,α,γ) stands for all 3m unknown parameters. Moreover, m ρ = 1 and ρ ,α ≥ 0,
j=1 j j j
for all j = 1,··· ,m. More details on this distribution can be found in TPzougas et al. (2014).
8
Tzougas et al. (2014) showed that a Pareto mixture distribution can be arrived by mixing two
exponential and inverse gamma distributions.
th
Now, to formulated a Pareto mixture regression model, suppose that for an i individual, informa-
tiononcontinuousresponsevariablesZ ,··· ,Z alongwithinformationonpcovariatesW ,··· ,W
i1 it 1 p
are available. Also suppose that Z given parameter Θ = θ has been distributed according to an
il il il
exponential distributionwithdensityfunctionf (z ) = exp{−z /θ }/θ .Moreover, suppose
Zil|Θil=θil il il il il
that parameter θ can be evaluated by the following regression model
il
p
log(θ ) = d + d w +ǫ ,
il 0il kil ki i
k=1
X
where d ,··· ,d are regression coefficients and u = exp(ǫ ) has been distributed according to a
0i pi i i
finite mixture Inverse gamma distribution with density function
m (α −1)αju−αj−1
f (u ) = ρ j i e−(αj−1)/ui, (7)
Ui i j Γ(α )
j
j=1
X
where m ρ = 1 and ρ , α ≥ 0. To have E(ǫ ) = 0 in Equation (6) we set γ = α − 1, for
j=1 j j j i j j
j = 1,·P·· ,m.
Using the law of total probability and setting b := d + p d w , one may show that
il 0il k=1 kil ki
P
∞
f (z ) = f (z )f (u )du
Zil|ϑ il Zil||ϑ,ui il Ui i i
Z0
∞ m (α −1)αju−αj−1
= e−zilexp{−bil}/uiexp{−b }u−1 ρ j i e−(αj−1)/uidu
il i j Γ(α ) i
Z0 j=1 j
X
m (α −1)αj
j
= ρ α .
j j(z +α −1)αj+1
j
j=1
X
Similar to the kINB regression/distribution the maximum liklihood estimator for parameters of a
Pareto mixture regression/distribution canbeobtainedusing theEMalgorithm. Fortunately, Rigby
& Stasinopoulos (2001) developed a R package, named ’GAMLSS‘, for such propose, see Rigby &
Stasinopoulos (2001, 2009) for more details.
9
4. Application to posteriori rate–making system
The rate–making system is a non-life actuarial system which rates policyholders based upon their
last t years record (Payandeh Najafabadi et al., 2015). A rate–making system based upon policy-
holders’ characteristics assigns a priori premium for each policyholder. Then, it employs the last
t years claims experience of each insured to update such priori premium and provides posteriori
premium (Boucher & Inoussa, 2014). The Bonus–Malus system is a commercial and practical ver-
sion of the rate–making system which takes into account current year policyholders’ experience to
determine their next year premium.
There is a considerable attention from authors to study rate–making systems (or Bonus–Malus
systems). For instance: Several mathematical tools for pricing a rate–making system has been pro-
vided by Lange (1969). Dionne & Vanasse (1989, 1992) employed available asymmetric information
under Poisson and Negative Binomial regression models to determine premium of a rate–making
system. In 1995, Lemaire designed an optimal Bonus–Malus system based on Negative Binomial
distribution. Pinquet (1997) considered Poisson and Lognormal distributions to design an optimal
Bonus–Malus system. Walhin & Paris (1999) considered a Hofmann’s distribution along with a
finite mixture Poisson distribution to evaluate elements of a Bonus–Malus system. The relatively
premium of a rate–making system under the exponential loss function has been evaluated by Denuit
& Dhaene (2001). In 2001, Frangos & Vrontos designed an optimal Bonus–Malus system using both
ParetoandNegativeBinomialdistributions. UsingthebivariatePoissonregressionmodelBermúdez
& Morata (2009) studied priori rate–making procedure for an automobile insurance database which
has two different types of claims. In 2011, Bermúdez & Karlis employed a Bayesian multivariate
Poisson model to determine premium of a rate–making system which has a non-ignorable correla-
tion between types of its claims. Boucher & Inoussa (2014) introduced a new model to determine
premium of a rate–making system whenever panel or longitudinal data are available. The Sichel
distribution along with a Negative Binomial distribution have been considered by Tzougas & Fran-
10