Poisson modelling of log normal claims for non life insurance

Question

Suppose we have a portfolio of $J = 1000$ policies. Assume further that the number of claims $N$ is Poisson distributed with intensity $µ = 0.01$, and that the claim sizes $Z_i$ follow a log-normal distribution with $\mathbb{E}(Z_i) = 2$ and $\text{sd}(Z_i) = 1.0, 3.0, 5.0$.

First, I want to find the parameters $ξ$ and $σ$ of the log-normal distribution for each of the three values of $\text{sd}(Z_i)$ and plot the probability density function of each of the three distributions. Secondly, I would like to compute the $95\%$ and $99\%$ reserve of the portfolio for each of the three sets of parameters.

For the first part, I have tried to set up the equation for the mean and variance of a log-normal distribution by,

$$\mathbb{E}(Z_i) = e^{\mu+\frac{\sigma^2}{2}}$$ $$\text{sd}(Z_i)=(e^{\sigma^2}-1)e^{\mu+\frac{\sigma^2}{2}}$$

From here, I insert my values for the mean and sd to find $\mu$ and $\sigma$, however for $\text{sd}(Z_i)=1$, there does noe seem to be a feasible solution for $\mu$ and $\sigma$.

For the second case, I use the rules for the possion distribution which says that the number of claims $N'$ over the portfolio is Poisson distirbuted with intensity $\mu=0.01\cdot 1000 = 10$

Now I define the random variable,

$$X=\sum_{i=1}^{N'}Z_i$$

And then I calculate the following probabilites,

$$P(X>q_{\epsilon})=0.05$$ $$P(X>q_{\epsilon})=0.01$$

To find the $95\%$ and $99\%$ reserve respectively.

Finally, If I want to add a deductible of $a=0.5$ and maximum of $b=3.0$ to this policy, how can I then recompute the $95\%$ and $99\%$ reserves for all three cases of $\text{sd}(Z_i)$?

Thanks in advance.

Answer 1

Since you reuse $\mu$ in two contexts, I will instead use $\lambda = 0.01$ for the Poisson intensity and reserve $\mu$ for the lognormal location parameter.

Next, your second equation is not correct for the standard deviation. $$\operatorname{E}[Z_i] = e^{\mu + \sigma^2/2}, \\ \operatorname{Var}[Z_i] = e^{2\mu + \sigma^2} (e^{\sigma^2} - 1)$$ requires that the standard deviation be $$\operatorname{SD}[Z_i] = \operatorname{E}[Z_i] \sqrt{e^{\sigma^2} - 1}.$$ Consequently, for the first pair of parameters, we have $$\log 2 = \mu + \frac{\sigma^2}{2}, \\ \frac{1}{4} = e^{\sigma^2} - 1,$$ hence $$\mu = \frac{1}{2} \log \frac{16}{5} \approx 0.581575, \quad \sigma = \sqrt{ \log \frac{5}{4}} \approx 0.472381.$$

Unfortunately, the sum of lognormal random variables does not have a nice closed form. The product remains lognormal, but not the sum. Therefore, I think the only algebraically tractable approach is to perform an approximation. The basic idea is to use a lognormal approximation that matches the mean and variance of the sum; namely, $$\sigma_X^2 = \log \left( (e^{\sigma^2} - 1) \frac{N e^{2\mu}}{ (N e^\mu )^2} + 1 \right) = \log \left( \frac{e^{\sigma^2} - 1}{N} + 1 \right), \\ \mu_X = \log N e^\mu + \frac{\sigma^2}{2} - \frac{\sigma_X^2}{2} = \mu + \log N + \frac{\sigma^2 - \sigma_X^2}{2}.$$ This clearly runs into problems if $N = 0$ with positive probability, but in such a case, the aggregate claim size is $0$. We can then use a computer to perform the relevant calculations; e.g., $$0.01 = \Pr[X > q_\epsilon] = \sum_{n = 1}^{1000} \Pr[X > q_\epsilon \mid N' = n]\Pr[N' = n] \tag{1}$$ where $N' \sim \operatorname{Poisson}(1000\lambda)$. This still requires us to solve for a sum of Poisson-weighted lognormal quantiles, which is hardly ideal. It is worth noting that when $N > 30$ the Poisson weight is tiny, so we don't need to evaluate the full sum.

If a deductible and policy limit apply, then we are best off performing simulations to compute the reserves. In fact, I would perform simulation for the above case as well, to confirm that our approximation serves as a good model, since there are multiple sources of imprecision.

For example, in Mathematica, I would use the following code:

Quantile[ParallelTable[Total[RandomVariate[LogNormalDistribution[
    Log[16/5]/2, Sqrt[Log[5/4]]], RandomVariate[PoissonDistribution[10]
    ]]], {10^6}], {0.95, 0.99}]

This gives me {32.4733, 38.6364} as the $95\%$ and $99\%$ reserves for $10^6$ simulations where $\mu, \sigma$ are as stated above. This also conveniently furnishes a starting point to evaluate the equation $(1)$:

msn[m_, s_, n_] := {m + Log[n] + (s^2 - #)/2, Sqrt[#]} &[Log[(Exp[s^2] - 1)/n + 1]]
Sum[(1 - CDF[LogNormalDistribution @@ msn[Log[16/5]/2, Sqrt[Log[5/4]], n], 38.6364)
    PDF[PoissonDistribution[10], n], {n, 1, 100}]

The result I obtained is $0.00996836$, noting that we truncated the sum at $n = 100$ which gives us more than enough precision for our purposes. In a previous version of this answer, I had gotten a much less reasonable result which I could not explain, but after restarting my session, something seems to have been fixed. This time, we get a very close result, which suggests the approximation is indeed good. The $99^{\rm th}$ percentile we find through this approximation is $38.62521$, which agrees with my simulations.

Subject to an individual ordinary deductible of $a = 0.5$ and policy limit of $b = 3.0$, we modify the code accordingly:

Quantile[ParallelTable[Total[Min[Max[0, # - 0.5], 3] & /@ RandomVariate[
         LogNormalDistribution[Log[15/6]/2, Sqrt[Log[5/4]]],
         RandomVariate[PoissonDistribution[10]]]], {10^6}, {0.95, 0.99}]

gives me {23.3995, 27.8966}, substantially lower since the heavy lognormal tail implies exposure to extreme aggregate losses, and imposing a policy limit curtails this effect.