Month: June 2017

Mathematical models for insurance payments – part 2 – ordinary deductible

Posted on Updated on

The post is the second post on models for insurance payments. The focus here is on the insurance payment when the loss is adjusted by a deductible. The first post is on insurance policy with a limit.

Ordinary Deductible

Suppose that the random variable X is the size of the loss (if a loss occurs). Under an insurance policy with a deductible d, the first d dollars of a loss is responsible by the insured and the amount of the loss in excess of d is paid by the insurer. Under this policy provision, what is the expected amount of insurance payment for a given loss? The deductible described here is called an ordinary deductible. The insurance payment (per loss) under the presence of an ordinary deductible d is denoted by (X-d)_+. The idea for this notation is that the insurance payment is X-d as long as X-d is not negative. The following describes the payment rule more explicitly.

    \displaystyle  (X-d)_+=\left\{ \begin{array}{ll}                     \displaystyle  0 &\ X \le d \\           \text{ } & \text{ } \\           \displaystyle  X-d &\ X > d           \end{array} \right.

The goal is to describe the probability model of the insurance payment (X-d)_+. First and foremost, we need to know how much it will cost the insurer on average to pay the insured, which is the mean E[(X-d)_+]. The variance of (X-d)_+ is also important as it is a measure of risk to the insurer. The mean E[(X-d)_+] and the variance Var[(X-d)_+] can be calculated by using the density function of X. At times it will be useful to know the the probability density function (PDF) and cumulative distribution function (CDF) of (X-d)_+ and other distributional quantities.

Moments

First, let’s focus on the moments of (X-d)_+. One way to obtain moments of (X-d)_+ is through the distribution of the unmodified loss X. Let f(x), F(x) and S(x)=1-F(x) be the probability density function (PDF), cumulative distribution function (CDF) and the survival function of the random loss, respectively. Then the moments of (X-d)_+ can be expressed using these distributional quantities of X.

    \displaystyle E[(X-d)_+]=\int_d^\infty (x-d) \ f(x) \ dx

    \displaystyle E[(X-d)_+^k]=\int_d^\infty (x-d)^k \ f(x) \ dx

    \displaystyle Var[(X-d)_+]=E[(X-d)_+^2]-E[(X-d)_+]^2

The above calculation is based on first principle. When the loss X is less than or equal to d, the insurance payment is zero. Note that the zero payment does not have to be explicitly stated in the above integrals. When the loss X exceeds the deductible d, the payment is X-d, the amount of the loss in excess of d. Example 1 and Example 2 below demonstrate how this calculation is performed.

Another way to obtain the moments of (X-d)_+ is obtained by first finding the PDF of (X-d)_+. This is discussed in a section below.

Basic Example

Example 1
Suppose that the loss distribution has PDF f(x)=1/5 (1- x/10) \ \ \ 0<x<10. Evaluate E(X) and E[(X-2)_+]. What is the expected amount of the loss eliminated as a result of imposing the ordinary deductible of 2 (eliminated from the insurer’s perspective)?

    \displaystyle \begin{aligned} E[(X-2)_+]&=\int_2^{10} (x-2) \ \frac{1}{5} (1- x/10) \ dx\\&=\int_2^{10} \frac{1}{5} (\frac{6}{5} x -x^2/10 -2) \ dx\\&=\frac{1}{5} \biggl[\frac{3 x^2}{5}-\frac{x^3}{30}-2x \biggr]_2^{10}=\frac{128}{75}=1.7067  \end{aligned}

From Example 1 in the previous post, E(X)=10/3. Thus when a loss occurs, the expected amount of the loss that is responsible by the insured is

    \displaystyle E(X)-E[(X-2)_+]=\frac{10}{3}-\frac{128}{75}=\frac{122}{75}=1.6267

,which is also the expected amount of the loss eliminated (from the insurer’s perspective) as a result of imposing an ordinary deductible.

Example 2
Continue with Example 1. Calculate the variance of the insurance payment (X-2)_+.

    \displaystyle \begin{aligned} E[(X-2)_+^2]&=\int_2^{10} (x-2)^2 \ \frac{1}{5} (1- x/10) \ dx \\ &=\int_2^{10} \frac{1}{5} \ (x^2-4x+4) \  (1- x/10) \ dx \\&=\int_2^{10} \frac{1}{5} (\frac{7}{5} x^2 -x^3/10-\frac{22}{5} +4) \ dx\\&=\frac{1}{5} \biggl[\frac{7 x^3}{15}-\frac{x^4}{40}-\frac{11x^2}{5} +4x \biggr]_2^{10}=\frac{512}{75}  \end{aligned}

    \displaystyle \begin{aligned} Var[(X-2)_+]&=E[(X-2)_+^2]-E[(X-2)_+]^2 \\ &=\frac{512}{75}-\biggl(\frac{128}{75} \biggr)^2 \\&=\frac{22016}{5625}=3.9140  \end{aligned}

From Example 2 in the previous post, Var(X)=50/9=5.5556. This shows that imposing a deductible that is significant enough has a variance reducing effect on the insurance payment. Had the deductible not present, there would be a greater fluctuation in the payment.

Connection with Policy Limit

The calculation in Example 1 and Example 2 are based on first principle. For many parametric families of distributions, the calculation for E[(X-d)_+] can be done using formulas. To that end, consider summing the payment (X-d)_+ and the payment X \wedge d.

    \displaystyle (X-d)_+ + (X \wedge d)=\displaystyle  \left\{ \begin{array}{ll}                     \displaystyle  0 &\ X \le d \\           \text{ } & \text{ } \\           \displaystyle  X-d &\ X > d           \end{array} \right. +\displaystyle  \left\{ \begin{array}{ll}                     \displaystyle  X &\ X \le d \\           \text{ } & \text{ } \\           \displaystyle  d &\ X > d           \end{array} \right. =X

Summing up each case produces the original loss X. Thus buying a policy with an ordinary deductible d along with a policy with a limit of d equals to full coverage. From a calculation standpoint, E[(X-d)_+] can be computed as follows:

    E[(X-d)_+]=E(X)-E(X \wedge d)

The advantage of using the above idea is that the calculation of E(X \wedge d) is routine for a large number of parametric families of distributions. For examples, most of the distributions in the table in this link have formulas for E(X \wedge x) and E[(X \wedge x)^k]. Then subtracting E(X \wedge x) into E(X) would give E[(X-d)_+]. The following table gives the limited expectation E(X \wedge x) for three distributions, taken from the table in the given link.

Limited Expectation
  • Exponential: E[X \wedge x]=\theta (1-e^{-x/\theta})
  • Pareto: \displaystyle E[X \wedge x]=\frac{\theta}{\alpha-1} \biggl[1-\biggl(\frac{\theta}{x+\theta} \biggr)^{\alpha-1} \biggr], \alpha \ne 1
  • Lognormal: \displaystyle E[(X \wedge x)^k]=\exp(k \mu+k^2 \sigma^2/2) \Phi \biggl(\frac{\log(x)-\mu-k \sigma^2}{\sigma} \biggr)+x^k \biggl[1-\Phi \biggl(\frac{\log(x)-\mu}{\sigma} \biggr) \biggr]

Examples are given to demonstrate how these formulas are used.

Example 3
You are given the following:

  • Losses follow a Pareto distribution with parameters \alpha=3 and \theta=500.
  • The coverage has an ordinary deductible of 100.

For the next loss that will occur, determine the expected amount that will be paid by the insurer.

The following shows the calculation.

    \displaystyle E(X)=\frac{\theta}{\alpha-1}=\frac{500}{2}=250

    \displaystyle \begin{aligned} E[X \wedge 1000]&=\frac{\theta}{\alpha-1} \biggl[1-\biggl(\frac{\theta}{100+\theta} \biggr)^{\alpha-1} \biggr] \\ &=\frac{500}{2} \biggl[1-\biggl(\frac{500}{600} \biggr)^{2} \biggr]=\frac{2750}{36}=76.389  \end{aligned}

    \displaystyle E[(X-1000)_+]=E(X)-E[X \wedge 1000]=250-\frac{2750}{36}=\frac{3125}{18}=173.611

As a result of imposing a deductible of 1000, the payment made by the insurer per loss goes from 250 to 173.611. The reduction of payment of E[X \wedge 100] = 76.389.

Example 4
You are given the following:

  • Losses follow an exponential distribution with mean 2500.
  • The coverage has an ordinary deductible of 1000.

Determine the expected amount per loss that is paid by the insurer.

The following shows the calculation.

    \displaystyle E(X)=2500

    \displaystyle E(X \wedge 1000)=2500 (1-e^{-1000/2500})=2500 (1-e^{-0.4})=824.20

    \displaystyle \begin{aligned} E[(X-1000)_+]&=E(X)-E[X \wedge 1000] \\ &=2500-2500 (1-e^{-0.4})  =2500 e^{-0.4}=1675.80 \end{aligned}

Note that the reduction of payment by the insurer per loss is \displaystyle E(X \wedge 1000)=824.20, which is the amount per loss that has to be met by the insured.

One observation about exponential loss. When the loss X is an exponential distribution, E[(X-d)_+] can actually be directly calculated using first principle since the exponential distribution is very tractable mathematically.

Example 5
You are given the following:

  • Losses follow a lognormal distribution with parameters \mu=6.5 and \sigma=1.75.
  • The coverage has an ordinary deductible of 1000.

Determine the expected amount per loss that is paid by the insurer.

The following shows the calculation.

    \displaystyle E(X)=e^{\mu+(1/2) \sigma^2}=e^{6.5+0.5(1.75)^2}=e^{8.03125}=3075.583751

    \displaystyle \begin{aligned} E[X \wedge 1000]&=e^{6.5+0.5(1.75)^2} \ \Phi \biggl(\frac{\log(1000)-6.5-1.75^2}{1.75} \biggr)+1000 \biggl[1-\Phi \biggl(\frac{\log(1000)-6.5}{1.75} \biggr) \biggr] \\ &=e^{8.03125} \ \Phi (-1.52)+1000 [1-\Phi (0.23) ]  \\ &=e^{8.03125} \ (1-0.9357)+1000 [1-0.5910 ]=606.7600352 \end{aligned}

    \displaystyle \begin{aligned} E[(X-1000)_+]&=E(X)-E[X \wedge 1000] \\ &=e^{8.03125}-606.7600352=2468.82 \end{aligned}

PDF and CDF of Insurance Payment

The preceding section shows that E[(X-d)_+]=E(X)-E(X \wedge d) is a great way to evaluate mean claim cost per loss when the loss distribution is from this list of distributions. For the distributions in this list, the limited expectation E(X \wedge d) is readily available. To calculate the higher moments of (X-d)_+, it is helpful to know more about its distribution. To that end, we derive its PDF and CDF and the survival function.

Before deriving the PDF and CDF, note that the payment variable (X-d)_+ is a censored and shifted random variable. It is censored from below and is shifted to the left by the amount of the deductible d. As far as payment is concerned, any loss that is below the deductible is considered zero. In effect, any loss in the interval (0,d) is recorded as zero. The insurance policy only pays for losses that are in excess of the deductible d. Hence the positive payments are derived by shifting losses in the interval (d,\infty) to the left by d. The resulting distribution for the payment (X-d)_+ is a mixed random variable. It has a point mass at y=0 with probability F_X(d) where F_X(x) is the CDF of the loss X. On the interval (0,\infty) or some appropriate interval (0,M), the density curve of (X-d)_+ is continuous and is the resulting of shifting the density curve of X to the left by the amount d. Let Y=(X-d)_+. Let f_X(x) and F_X(x) be the PDF and CDF of X, respectively. The following is the PDF of Y.

    \displaystyle  f_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  F_X(d) &\ \ \ \ y=0 \\           \text{ } & \text{ } \\           \displaystyle  f_X(y+d) &\ \ \ \ y > 0           \end{array} \right.

Note the point mass at y=0, which is the probability that the loss is less than or equal to the deductible. The curve f_X(y+d) is the density curve f_X(y) shifted to the left by the amount d. The following is the CDF and the survival function of Y=(X-d)_+.

    \displaystyle  F_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  0 &\ \ \ \ y<0 \\           \text{ } & \text{ } \\           \displaystyle  F_X(y+d) &\ \ \ \ y \ge 0                      \end{array} \right.

    \displaystyle  S_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  1 &\ \ \ \  y<0 \\           \text{ } & \text{ } \\           \displaystyle  S_X(y+d) &\ \ \ \ y \ge 0                      \end{array} \right.

Now that we have a description of the model for Y=(X-d)_+, moments and other distributional quantities can be derived. Consider the following two examples.

Example 6
You are given the following:

  • Losses follow a uniform distribution on the interval (0, 10).
  • The coverage has an ordinary deductible of 4.

Determine the mean and variance of Y=(X-4)_+.

The following shows the PDF, CDF and survival function of the loss X.

    \displaystyle f_X(x)=\frac{1}{10} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 0<x<10

    \displaystyle F_X(x)=\frac{x}{10} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 0<x<10

    \displaystyle S_X(x)=1-\frac{x}{10} \ \ \ \ \ \ \ \ \ \ \ \ \ 0<x<10

The point mass is F_X(4)=0.4. The following is the PDF f_Y(y).

    \displaystyle  f_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  0.4 &\ \ \ \ y=0 \\           \text{ } & \text{ } \\           \displaystyle  \frac{1}{10} &\ \ \ \ 0<y < 6           \end{array} \right.

Note that the PDF f_X(x) over the interval (4,10) is shifted to the interval (0, 6). The following gives the CDF and survival function of Y.

    \displaystyle  F_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  0 & \ \ \ \ y < 0 \\           \text{ } & \text{ } \\           \displaystyle  \frac{y+4}{10} & \ \ \ \ 0 \le y < 6 \\            \text{ } & \text{ } \\           \displaystyle  1 & \ \ \ \ y \ge 6            \end{array} \right.

    \displaystyle  S_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  1 & \ \ \ \ y < 0 \\           \text{ } & \text{ } \\           \displaystyle  1-\frac{y+4}{10}=\frac{6-y}{10} & \ \ \ \ 0 \le y < 6 \\            \text{ } & \text{ } \\           \displaystyle  0 & \ \ \ \ y \ge 6            \end{array} \right.

The moments of Y=(X-4)_+ can be can be obtained by integrating x^k with the PDF f_Y(y).

    \displaystyle  E(Y)=\int_0^6 \frac{1}{10} \ y \ dy=\frac{9}{5}=1.8

    \displaystyle  E(Y^2)=\int_0^6 \frac{1}{10} \ y^2 \ dy=\frac{36}{5}=7.2

    Var(Y)=7.2-1.8^2=3.96

Note that the variance of the unmodified loss is Var(X)=100/12=8.33. In this case, imposing a deductible of 4 reduces the variance by a little more than half.

Example 7
Consider the loss distribution and the coverage in Example 3. Determine the PDF, CDF and the survival function of Y=(X-100)_+. Then find the mean and variance of Y=(X-100)_+.

The loss X has a Pareto distribution with parameters \alpha=3 and \theta=500. The following shows its PDF, CDF and survival function.

    \displaystyle f_X(x)=\frac{3 \cdot 500^3}{(x+500)^4} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ x>0

    \displaystyle F_X(x)=1-\biggl(\frac{500}{x+500} \biggr)^3 \ \ \ \ \ \ \ \ x>0

    \displaystyle S_X(x)=\biggl(\frac{500}{x+500} \biggr)^3 \ \ \ \ \ \ \ \ \ \ \ \ \ x>0

The point mass is F_X(100)=91/216. The following is the PDF f_Y(y).

    \displaystyle  f_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  \frac{91}{216} &\ y=0 \\           \text{ } & \text{ } \\           \displaystyle  \frac{3 \cdot 500^3}{(y+600)^4} &\ y > 0           \end{array} \right.=  \left\{ \begin{array}{ll}                     \displaystyle  \frac{91}{216} &\ y=0 \\           \text{ } & \text{ } \\           \displaystyle  \biggl(\frac{500}{600} \biggr)^3 \ \frac{3 \cdot 600^3}{(y+600)^4} &\ y > 0           \end{array} \right.

Note that the part for y > 0 is S_X(100) times the Pareto density function for parameters \alpha=3 and \theta=600. Thus the moments of Y=(X-1000)_+ can be obtained by multiplying S_X(1000) with the moments of this new Pareto distribution. Since the point mass is for y=0, we can ignore it.

    \displaystyle E(Y)=E[(X-100)_+]=\biggl(\frac{500}{600} \biggr)^3 \ \frac{600}{3-1}=\frac{3125}{18}=173.611

    \displaystyle E(Y^2)=E[(X-1000)_+^2]=\biggl(\frac{125}{216} \biggr) \ \frac{600^2 \ 2}{(3-1)(3-2)}=\frac{625000}{3}

    \displaystyle Var(Y)=E(Y^2)-E(Y)^2=\frac{625000}{3}-\biggl(\frac{3125}{18}\biggr)^2=178192.5154

    \sigma_Y=\sqrt{178192.5154}=422.129

To complete the example, the following gives the CDF and survival function.

    \displaystyle  F_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  0 & \ \ \ \ y < 0 \\           \text{ } & \text{ } \\           \displaystyle  1-\biggl( \frac{500}{600+y} \biggr)^3 & \ \ \ \ 0 \le y < \infty                      \end{array} \right.

    \displaystyle  S_Y(y)=\left\{ \begin{array}{ll}                     \displaystyle  1 & \ \ \ \ y < 0 \\           \text{ } & \text{ } \\           \displaystyle  \biggl( \frac{500}{600+y} \biggr)^3 & \ \ \ \ 0 \le y < \infty                      \end{array} \right.

Insurance Payment Per Loss versus Per Payment

Note that the E[(X-d)_+] is calculated over all losses (whether they are below or above the deductible). The expected value E[(X-d)_+] reflects the probability P(X \le d) that is assigned to the payment of zero (the point mass), and the probabilities that are applied to the payment x-d over all x in the interval (d,\infty) (the integral). Thus E[(X-d)_+] is the average payment per loss (or over all possible losses). Another average payment to consider is the average payment over all payments, i.e. over all losses larger than the deductible. This topic is continued in the next post.

\text{ }

\text{ }

\text{ }

\copyright 2017 – Dan Ma

Advertisements

Mathematical models for insurance payments – part 1 – policy limit

Posted on Updated on

Suppose an individual (or an entity) faces a random loss and the loss is modeled by the random variable X. The individual purchases an insurance coverage to cover this loss. If the policy pays the loss in full (if a loss occurs), then the expected insurance payment is mean loss E(X). Usually the expected insurance payment is less than E(X) due to the presence of policy provisions such as a deductible and/or limit. This post is the first post in a series of posts in discussing the probability models of insurance payments.

Policy Limit

In this post, we assume that the random loss X is a continuous random variable taking on the positive real numbers or numbers from an interval such as (0,M). The methodology can be adjusted to handle discrete loss variable (mostly replacing integrals with summation).

Suppose that the random variable X is the size of the loss (if a loss occurs). Under an insurance policy with a limit u, the insurer pays the loss in full if the loss X is less than the limit. Furthermore if the loss is u or greater, the insurance payment is capped at u. Under this policy provision, what is the expected amount of insurance payment if there is a loss? The insurance payment in the presence of a limit u is denoted by X \wedge u. This is called the limited loss random variable.

    \displaystyle  X \wedge u=\left\{ \begin{array}{ll}                     \displaystyle  X &\ X \le u \\           \text{ } & \text{ } \\           \displaystyle  u &\ X > u           \end{array} \right.

First, let’s consider the mean E(X \wedge u) (called limited expectation). Let f(x), F(x) and S(x)=1-F(x) be the probability density function (PDF), cumulative distribution function (CDF) and the survival function of the random loss, respectively.

    \displaystyle \begin{aligned} E(X \wedge u)&=\int_0^u x \ f(x) \ dx+ u \ S(u) \\&=\int_0^u S(x) \ dx  \end{aligned}

The first integral is based on the definition of the limited loss variable, summing up the payments in the interval (0,u) and in the interval (0, \infty). The second integral is an alternative way to evaluate the first integral (see here for more explanation). The following gives the second moment, which can be used in evaluating the variance of insurance payment.

    \displaystyle \begin{aligned} E[(X \wedge u)^2]&=\int_0^u x^2 \ f(x) \ dx+ u^2 \ S(u) \\&=\int_0^u 2x \ S(x) \ dx  \end{aligned}

    Var(X \wedge u)=E[(X \wedge u)^2]-E(X \wedge u)^2

Once again, the second moment is also expressed using the alternative calculation. The following two examples demonstrate the evaluation of the first two moments of X \wedge u.

\text{ }

Basic Example

Example 1
Suppose that the loss distribution has PDF f(x)=1/5 (1- x/10) \ \ \ 0<x<10. Evaluate E(X) and E(X \wedge 4)

The CDF of the loss X is F(x)=1/5 (x -x^2/20).

    \displaystyle E(X)=\int_0^{10} x \ 1/5 (1- x/10) \ dx=\frac{1}{5} \biggl[\frac{x^2}{2}-\frac{x^3}{30} \biggr]_0^{10}=\frac{10}{3}=3.33

    \displaystyle \begin{aligned} E(X \wedge 4)&=\int_0^4 x \ 1/5 (1- x/10) \ dx +4 \ S(4)\\&=\int_0^4 1/5 (x -x^2/10) \ dx+4 \ S(4)\\&=\frac{1}{5} \biggl[\frac{x^2}{2}-\frac{x^3}{30} \biggr]_0^4+4 \ \frac{9}{25} \\&=\frac{196}{75}=2.6133  \end{aligned}

    \displaystyle \begin{aligned} E(X \wedge 4)&=\int_0^4 \biggl[1-\frac{1}{5} \biggl(x-\frac{x^2}{20} \biggr) \biggr] \ dx \\&=\int_0^4 \biggl[1-\frac{x}{5}+\frac{x^2}{100} \biggr] \ dx \\&=x-\frac{x^2}{10}+\frac{x^3}{300} \biggr|_0^4 \\&=\frac{196}{75}=2.6133  \end{aligned}

The limited expected value is calculated in two ways to demonstrate that the two approaches are equivalent. Without the policy limit, the expected loss is 3.3333. With the benefit payment capped at 4, the expected amount paid by the insurance policy is 2.6133 (per loss). The difference is 3.3333 – 2.6133 = 0.72, which is the expected loss responsible by the insured.

Example 2
Continue with Example 1. Calculate the variance of the insurance payment X \wedge 4.

    \displaystyle \begin{aligned} E[(X \wedge 4)^2]&=\int_0^4 x^2 \ 1/5 (1- x/10) \ dx +16 \ S(4)\\&=\int_0^4 1/5 (x^2 -x^3/10) \ dx+16 \ \frac{9}{25} \\&=\frac{1}{5} \biggl[\frac{x^3}{3}-\frac{x^4}{40} \biggr]_0^4+16 \ \frac{9}{25} \\&=\frac{656}{75}  \end{aligned}

    \displaystyle Var(X \wedge 4)=\frac{656}{75}- \biggl(\frac{196}{75} \biggr)^2=\frac{10784}{5625}=1.9172.

To compare, the variance of X is \displaystyle Var(X)=\frac{50}{9}=5.5556 as calculated below.

    \displaystyle \begin{aligned} E[X^2]&=\int_0^{10} x^2 \ 1/5 (1- x/10) \ dx \\&=\int_0^{10} 1/5 (x^2 -x^3/10) \ dx \\&=\frac{1}{5} \biggl[\frac{x^3}{3}-\frac{x^4}{40} \biggr]_0^{10}=\frac{50}{3}  \end{aligned}

    \displaystyle Var(X)=\frac{50}{3}- \biggl(\frac{10}{3} \biggr)^2=\frac{50}{9}=5.5556

The example demonstrates that policy provisions such as limit has a variance reducing effect. Had the insurer been liable to pay for the entire loss amount, there would be a greater fluctuation in the payment.

A Censored Variable

Mathematically, the random variable X \wedge u is an upper censored random variable. Any realized value of X above the limit is known only known as the limit u (as far as payment is concerned). On the other hand, applying a limit on the loss X turns the continuous random variable X into a mixed random variable. In the interval (0,u), X \wedge u is continuous while there is a point mass at the point X=u with probability mass P[(X \wedge u)=u]=S(u).

The policy provisions such as deductible and limit have the effect of turning the unmodified loss variable X into censored or truncated variables. As a result, these variables are mixed random variables. The subsequent posts will further demonstrate this point.

Parametric Distributions

There is a vast inventory of parametric distributions that are potential candidates for models of random losses. For example, the table in this link is an inventory of distributions that may be useful in modeling insurance losses. Included in the table are the limited expectations for various distributions. The following table shows three of the distributions. They are listed here primarily because the calculation is tractable and is better suited for demonstration of the concept.

Limited Expectation
  • Exponential: E[X \wedge x]=\theta (1-e^{-x/\theta})
  • Pareto: \displaystyle E[X \wedge x]=\frac{\theta}{\alpha-1} \biggl[1-\biggl(\frac{\theta}{x+\theta} \biggr)^{\alpha-1} \biggr], \alpha \ne 1
  • Lognormal: \displaystyle E[(X \wedge x)^k]=\exp(k \mu+k^2 \sigma^2/2) \Phi \biggl(\frac{\log(x)-\mu-k \sigma^2}{\sigma} \biggr)+x^k \biggl[1-\Phi \biggl(\frac{\log(x)-\mu}{\sigma} \biggr) \biggr]

For formulas for the higher moments, refer to the table. This table will be used throughout the discussion (in subsequent blog posts) on estimating insurance payments. Here’s some blog posts on these three distributions – exponential, Pareto and lognormal.

The following example shows how these formulas are used.

Example 3
Evaluate E(X \wedge u) for the following loss distributions.

  • Exponential: mean 1000, u = 2000.
  • Pareto: mean 5, variance 75, u = 10.
  • Lognormal: mean 177.682811, variance 13680.72152, u = 250.

The exponential distribution with mean 1000 has parameter \theta=1000 (the scale parameter). The following gives the limited expected value.

    E(X \wedge 2000)=1000 (1-e^{-2000/1000})=864.6647

\text{ }

The Pareto distribution with mean 5 and variance 75 translates to the parameters \alpha=3 and \theta=10. The following gives the limited expected value.

    \displaystyle E(X \wedge 10)=\frac{10}{3-1} \biggl[1-\biggl(\frac{10}{10+10} \biggr)^{3-1} \biggr]=3.75.

\text{ }

The lognormal distribution with mean 177.682811 and variance 13680.72152 corresponds to \mu=5 and \sigma=0.6. The following gives the limited expected value.

    \displaystyle \begin{aligned} E[X \wedge 250]&=e^{5+0.6^2/2} \ \Phi \bigg(\frac{\log(250)-5-0.6^2}{0.6} \biggr)+250 \ \biggl[1-\Phi \bigg(\frac{\log(250)-5}{0.6} \biggr) \biggr] \\&=e^{5.18} \ \Phi (0.269)+250 \ [1-\Phi (0.869 )]\\&=e^{5.18} \ 0.6064+250(1-0.8072)=155.7969 \end{aligned}

Comments

The discussion here is just the beginning. The limited expected value E(X \wedge u) is for a simple policy provision, having just one modification of loss. It is a building block for other payments under more complicated insurance payment rules.

\text{ }

\text{ }

\text{ }

\copyright 2017 – Dan Ma

Mathematical models for insurance payments – part 0 – introduction

Posted on

Insurance losses depend on two random variables. The first one is the number of losses that will occur in a specified period in connection to an insured (or a group of insureds). This is commonly referred to as the loss frequency or claim frequency. Its probability distribution is called the loss (claim) frequency distribution. The second random variable is the amount of the loss (if a loss has occurred). The amount of loss is usually referred to as the severity and its probability distribution is called the severity distribution. Then putting these two distributions together will lead to the total loss distribution.

In the several posts that follow, the focus is on the severity distribution, i.e. we will focus on the size of the loss or size of the claim. Once the methodology for the severity is discussed in a fair amount of details, we will add claim frequency.

Given that a loss has occurred and that the amount is X, the insurance company will make a payment to the insured to cover the loss X. The insurance payment is a random quantity since X is random. Due to the presence of policy provisions such as deductible and limit, the insurance payment is likely less than X. The focus is on determining the distributional quantities of the probability model of the insurance payment to the insured. First and foremost, we would like to calculate the mean and variance of the distribution of payment as well as probability density function and cumulative distribution function and other distributional quantities.

The subsequent post will show, from a mathematical standpoint, the random variable of the insurance payment is a truncated and/or censored variable of the loss X due to policy provisions such as deductible and policy limit.

In practice, the mean and variance of the insurance payments are often estimated from historical observations rather than using a hypothesized distribution of losses. However, the parametric distribution approach is a good starting point of the discussion of estimation of insurance losses.

\text{ }

\text{ }

\text{ }

\text{ }

\copyright 2017 – Dan Ma

Integrating survival function to calculate the mean

Posted on Updated on

For a continuous random variable X that are non-negative in values, the mean E(X) is obtained by integrating x f(x) from 0 to \infty where f(x) is the probability density function of X. In some cases, this integral is hard (or impossible) to do. In these cases, it may be possible to find the mean and higher moments by integrating the survival functions.

Here’s the usual way to calculate moments.
\displaystyle (1) \ \ \ \ E(X)=\int_0^\infty x f(x) \ dx
\text{ }
\displaystyle (2) \ \ \ \ E(X^k)=\int_0^\infty x^k f(x) \ dx

Moments can be calculated by integrating the survival function.
\displaystyle (3) \ \ \ \ E(X)=\int_0^\infty S(x) \ dx
\text{ }
\displaystyle (4) \ \ \ \ E(X^k)=\int_0^\infty k x^{k-1} S(x) \ dx

In fact, the mean and moments of many of the distributions in the Exam C table are calculated using the survival function method, e.g. the Pareto distribution and the exponential distribution.

The integrals in (3) and (4) are derived from (1) and (2) by performing integration by parts. The derivation is not very hard. You will find out that some assumptions are needed to make the integration by parts work. Specifically, the assumption is that the following limit is 0.

    \displaystyle \lim_{x \rightarrow \infty} x^k \ S(x)= 0

The integrals in (3) and (4) works only if this limit converges to zero. Essentially this limit is zero if the expectation E(X^k) exists.

The same idea can be applied to a discrete distribution.

\displaystyle (5) \ \ \ \ E(X)=\sum_{x=0}^\infty x P(X=x)
\text{ }
\displaystyle (6) \ \ \ \ E(X^k)=\sum_{x=0}^\infty x^k P(X=x)

\displaystyle (7) \ \ \ \ E(X)=\sum_{x=0}^\infty P(X>x)
\text{ }
\displaystyle (8) \ \ \ \ E(X^k)=\sum_{x=0}^\infty [(x+1)^k-x^k] P(X>x)

\text{ }

\text{ }

Examples

Example 1
Suppose that the useful life (in months) of a device is modeled by the CDF F(x)=1-(1- x/120)^{1/2} for 0 \le x \le 120. Calculated the expected useful of a brand new device. Calculate the variance of the lifetime of such a device.

Note that the survival function is S(x)=(1- x/120)^{1/2}. The following shows the calculation.

    \displaystyle E(X)=\int_0^{120} (1- x/120)^{1/2} \ dx=80

    \displaystyle \begin{aligned} E(X^2)&=\int_0^{120} 2x \ (1- x/120)^{1/2} \ dx \\&=\int_1^{0} 2 \cdot 120(1-U) \ U^{1/2} (-120) \ dU \\&=\int_0^{1} 28800 \ (U^{1/2}-U^{3/2}) \ dU \\&=28800 \  \biggl[\frac{2}{3} U^{3/2}-\frac{2}{5} U^{5/2} \biggr]_0^1 \\&=7680 \end{aligned}

    Var(X)=E(X^2)-E(X)^2=7680-80^2=1280

Note that the integral for E(X^2) uses the method of substitution.

Example 2
Suppose that the lifetime of a machine follows a distribution with the following CDF.

    \displaystyle F(x)=1+e^{-x/6}-2 e^{-x/12} \ \ \ \ \ x>0

Calculate the mean and variance of the lifetime.

The survival function is \displaystyle S(x)=2 e^{-x/12}-e^{-x/6}. The first is to calculate the first 2 moments.

    \displaystyle \begin{aligned} E(X)&=\int_0^\infty (2 e^{-x/12}-e^{-x/6}) \ dx \\&=\int_0^\infty (2 \cdot 12 \frac{1}{12} e^{-x/12}-6 \ \frac{1}{6} e^{-x/6}) \ dx \\&=24 \int_0^\infty \frac{1}{12} e^{-x/12} \ dx- 6 \int_0^\infty \frac{1}{6} e^{-x/6} \ dx=24-6=18  \end{aligned}

    \displaystyle \begin{aligned} E(X^2)&=\int_0^\infty 2x (2 e^{-x/12}-e^{-x/6}) \ dx \\&=\int_0^\infty (4x e^{-x/12}-2x e^{-x/6}) \ dx  \\&=\int_0^\infty (4 \cdot 12 \ x \frac{1}{12} e^{-x/12}-2 \cdot 6 \ x \ \frac{1}{6} e^{-x/6}) \ dx \\&=48 \int_0^\infty x \frac{1}{12} e^{-x/12} \ dx -12 \int_0^\infty x \ \frac{1}{6} e^{-x/6} \ dx=48(12)-12(6)=504  \end{aligned}

    Var(X)=E(X^2)-E(X)^2=504-18^2=180

Note that the integral for E(X) is manipulated to be in terms of integrals of exponential density functions while the integral for E(X^2) is manipulated to be in terms of integrals of x times exponential density functions (each becoming the mean of that exponential distribution).

Example 3
The two-parameter Pareto distribution (Type II Lomax) is discussed in here in a companion blog. Its survival function is

    \displaystyle S(x)=\biggl(\frac{\theta}{x+\theta} \biggr)^\alpha \ \ \ \ \ x>0

Derive the mean of the Pareto distribution using the survival function approach. What assumption is made on the shape parameter \alpha?

    \displaystyle \begin{aligned} E(X)&=\int_0^\infty \biggl(\frac{\theta}{x+\theta} \biggr)^\alpha \ dx \\&=\int_0^\infty \theta^\alpha (x+\theta)^{- \alpha} \ dx \\&=\frac{\theta^\alpha}{-\alpha+1} \   (x+\theta)^{-\alpha+1} \biggr|_0^\infty \\&=\frac{\theta^\alpha}{-\alpha+1} \ \frac{1}{(x+\theta)^{\alpha-1}} \biggr|_0^\infty =\frac{\theta}{\alpha-1} \end{aligned}

In order for the integral to converge, we need to assume \alpha>1.

The next two examples are left as exercises.

Practice Problems

Practice Problem 1
The survival function for the distribution of a lifetime of a type of electronic devices is \displaystyle S(x)=\frac{1}{40} \ e^{-x/10}+\frac{1}{20} \ e^{-x/15} \ \ \ \ \ x>0.

Calculate the mean and variance of the lifetime of such devices.

Practice Problem 2
The probability that the size of a randomly selected auto collision loss is greater than x is \displaystyle S(x)=\biggl(1-\frac{x}{10} \biggr)^6 \ \ \ \ 0<x<10.

Calculate the mean and variance of the loss size.

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

Answers

Practice Problem 1

    \displaystyle E(X)=1

    \displaystyle Var(X)=\frac{53}{2}=26.5

Practice Problem 2

    \displaystyle E(X)=\frac{10}{7}

    \displaystyle Var(X)=\frac{75}{49}

\text{ }

\text{ }

\text{ }

\copyright 2017 – Dan Ma

Practice Problem Set 3 – basic lognormal problems

Posted on Updated on

This post has several practice problems to go with this previous discussion on lognormal distribution.

Practice Problem 3A
The amount of annual losses from an insured follows a lognormal distribution with parameters \mu and \sigma = 0.6 and with mode = 2.5. Calculate the mean annual loss for a randomly selected insured.

Practice Problem 3B
Claim size for an auto insurance coverage follows a lognormal distribution with mean 149.157 and variance 223.5945. Determine the probability that a randomly selected claim will be greater than 170.

Practice Problem 3C
For x-ray machines produced by a certain manufacturer, the following is known.

  • Lifetime in years follows a lognormal distribution with \mu = 0.9 and \sigma.
  • The expected lifetime of such machines is 15 years.

Calculate the probability that an x-ray machine produced by this manufacturer will last at least 12 years.

Practice Problem 3D
Claim sizes expressed in US dollars follow a lognormal distribution with parameters \mu = 5 and \sigma = 0.25. One Canadian dollar is currently worth $0.75 US dollars. Calculate the 75th percentile of a claim in Canadian dollars.

Practice Problem 3E

For a commercial fire coverage, the size of a loss follows a lognormal distribution with parameters \mu = 2.75 and \sigma = 0.75. Determine y-x where y is the 75th percentile of a loss and x is the 25th percentile of a loss. Note that y-x is known as the interquartile range.

Practice Problem 3F

Claim sizes in the current year follow a lognormal distribution with \mu = 4.75 and \sigma = 0.25. In the next year, all claims are expected to be inflated uniformly by 25%.

One claim is expected in the next year for an insured. Determine y-x where y is the 80th percentile of the size of this claim and x is the 40th percentile of the size of this claim.

Practice Problem 3G
  • In the current year, losses follow a lognormal distribution with \mu = 1.6 and \sigma = 1.35.
  • In the next year, inflation of 20% will impact all losses uniformly.
  • Determine the median of the portion of next year’s loss distribution that is above 10.

Practice Problem 3H

Losses follow a lognormal distribution with mean 17 and variance 219. Determine the skewness of the loss distribution.

Practice Problem 3I
  • Losses from an insurance coverage follow a lognormal distribution with parameters \mu and \sigma = 2.
  • The 80th percentile of the losses is 5884.
  • Determine the probability that a loss is less than 5000.

Practice Problem 3J
  • Losses from an insurance coverage follow a lognormal distribution.
  • The 25th percentile of the losses is 133.62.
  • The 75th percentile of the losses is 997.25.
  • Determine the mean of the losses.

All normal probabilities are obtained by using the normal distribution table found here.

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

Problem Answer
3A 4.29
3B 0.0869
3C 0.2033
3D 233.9675
3E 16.39085
3F 42.5155
3G 21.143268
3H 3.271185
3I 0.7764
3J 1124.394559

______________________________________________________________________________________________________________________________
\copyright 2017 – Dan Ma

Basic properties of lognormal distribution

Posted on Updated on

A detailed discussion of the mathematical properties of lognormal distribution is found in this previous post in a companion blog. This post shows how to work basic calculation problems for lognormal distribution. A summary of lognormal distribution is given and is followed by several examples. Practice problems are in the next post.

______________________________________________________________________________________________________________________________

Basic Properties

The random variable Y is said to follow a lognormal distribution with parameters \mu and \sigma if \log(Y) follows a normal distribution with mean \mu and variance \sigma^2. Here, \log is the natural logarithm in base e = 2.718281828…. It is difficult (if not impossible) to calculate probabilities by integrating the lognormal density function. Since the lognormal distribution is intimately related to the normal distribution, the basic lognormal calculation is performed by calculating the corresponding normal distribution. The following summary shows how.

In the following points, Y has a lognormal distribution with parameters \mu and \sigma and X=\log(Y) is the corresponding normal distribution with mean \mu and variance \sigma^2. The notation \exp(t) means raising e to the number t.

1. Lognormal observations and normal observations
  • Taking natural log of a lognormal observation gives a normal observation.
  • Raising e to a normal observation gives a lognormal observation.
2. Lognormal CDF and normal CDF
  • \Phi(z) is the CDF of the standard normal distribution.
  • \displaystyle F_Y(y)=\Phi \biggl(\frac{\log(y)-\mu}{\sigma} \biggr)
  • In words, lognormal CDF evaluated at y equals to standard normal CDF evaluated at \frac{\log(y)-\mu}{\sigma}.
    • Derivation:
      \displaystyle \begin{aligned} F_Y(y)=P(Y \le y)&=P[\log(Y) \le \log(y)] \\&=P \biggl[\frac{\log(Y)-\mu}{\sigma} \le \frac{\log(y)-\mu}{\sigma} \biggr] \\&=\Phi \biggl(\frac{\log(y)-\mu}{\sigma} \biggr)  \end{aligned}

3. Lognormal density function and normal density function
  • normal density: \displaystyle f_X(x)=\frac{1}{\sqrt{2 \pi} \ \sigma} \exp(-\frac{(x-\mu)^2}{2 \sigma^2})
  • lognormal density: \displaystyle F_Y(y)=\frac{1}{\sqrt{2 \pi} \ \sigma \ y} \exp(-\frac{(\log(y)-\mu)^2}{2 \sigma^2})
  • In words, lognormal density evaluated at y equals to normal density evaluated at \log(y) times \frac{1}{y}.
4. Lognormal moments and normal moment generating function
  • normal mgf: M(t)=e^{\mu \ t+\frac{1}{2} \sigma^2 \ t^2}=\exp(\mu \ t+\frac{1}{2} \sigma^2 \ t^2)
  • lognormal moment: E(x^k)=M(k)=\exp(\mu \ k+\frac{1}{2} \sigma^2 \ k^2)
  • In words, lognormal kth raw moment equals to normal mgf evaluated at k.
5. Examples of lognormal moments
  • E(Y)=E[e^{\log(Y)}]=M(1)=e^{\mu+\frac{1}{2} \sigma^2}
  • E(Y^2)=E[e^{2 \log(Y)}]=M(2)=e^{2 \mu+2 \ \sigma^2}
  • Var(Y)=e^{2 \ \mu+2 \ \sigma^2}-e^{2 \mu+ \sigma^2}=e^{2 \ \mu+\sigma^2} (e^{\sigma^2}-1)
  • skewness: \gamma_1=(e^{\sigma^2}+2) \sqrt{e^{\sigma^2}-1}
  • kurtosis: e^{4 \sigma^2}+2 e^{3 \sigma^2}+3 e^{2 \sigma^2}-3
6. Lognormal percentiles and normal percentiles
  • (100p)th percentile of the normal distribution is \displaystyle z_p.
  • (100p)th percentile of the normal distribution with mean \mu and variance \sigma^2 is \displaystyle \mu+z_p \times \sigma.
  • (100p)th percentile of the lognormal distribution with parameters \mu and \sigma is \displaystyle e^{\mu+z_p \times \sigma}.
  • In words, to find the (100p)th percentile of the lognormal distribution, find the (100p)th percentile of the corresponding normal distribution and then raise e to it.
7. Constant multiple of lognormal distribution
  • Let c>0. If Y has a lognormal distribution with parameters \mu and \sigma, then c Y has a lognormal distribution with parameters \mu+\log(c) and \sigma.
  • The effect of the multiplicative constant is on the parameter \mu in the form of an additive adjustment of \log(c).
8. Mode of lognormal distribution
  • \displaystyle e^{\mu - \sigma^2}

______________________________________________________________________________________________________________________________

Examples

Two examples are given to illustrate the calculation discussed here. The next post has practice problems.

All normal probabilities are obtained by using the normal distribution table found here.

Example 1
Suppose that the random variable Y has a lognormal distribution with parameters \mu = 1 and \sigma = 2. Calculate the following.

  • P(Y \le 75.19) and P(Y > 0.9)
  • The 67th, 95th and 99th percentiles of Y.
  • Let Y_1=1.1Y. Find P(Y_1 \le 75.19) and P(Y_1 > 0.9)

\text{ }

    \displaystyle P(Y \le 75.19)=\Phi \biggl(\frac{\log(75.19)-1}{2} \biggr)=\Phi(1.66)=0.9515

    \displaystyle \begin{aligned}P(Y > 0.9)&=1-\Phi \biggl(\frac{\log(0.9)-1}{2} \biggr) \\&=1-\Phi(-0.55) \\&=1-(1-0.7088) \\&=0.7088  \end{aligned}

_______________________________________________

To find the percentiles, first find the standard normal percentiles, either by using calculator or by looking up a table. Using a standard normal table, we get 0.44 (67th percentile), 1.645 (95th percentile) and 2.33 (99th percentile). The following gives the lognormal percentiles.

    \displaystyle e^{1+0.44(2)}=e^{1.88} = 6.5535 (67th percentile)

    \displaystyle e^{1+1.645 (2)}=e^{4.29} = 72.9665 (95th percentile)

    \displaystyle e^{1+2.33(2)}=e^{5.66} = 287.1486 (99th percentile)

_______________________________________________

The random variable Y_1=1.1Y has a lognormal distribution with parameters \mu=1+\log(1.1) and \sigma = 2.

    \displaystyle P(Y_1 \le 75.19)=\Phi \biggl(\frac{\log(75.19)-1-\log(1.1)}{2} \biggr)=\Phi(1.61)=0.9463

    \displaystyle \begin{aligned}P(Y_1 > 0.9)&=1-\Phi \biggl(\frac{\log(0.9)-1-\log(1.1)}{2} \biggr) \\&=1-\Phi(-0.57) \\&=1-(1-0.7157) \\&=0.7157  \end{aligned}

Note. One interpretation of Y_1=1.1Y is that of inflation, in this case a 10% inflation. For example, let Y be the size of a randomly selected auto insurance collision claim in the current year. If the claims are expected to increase 10% in the following year, Y_1=1.1Y is the the size of a randomly selected claim in the following year.

Example 2
Suppose that the random variable Y has a lognormal distribution with mean 12.18 and variance 255.02. Calculation the following.

  • P(Y > 10)
  • The skewness and kurtosis of Y.

\text{ }

First, determine the parameters \mu and \sigma by setting up the following equations.

    \displaystyle E(Y)=e^{\mu+\frac{1}{2} \sigma^2}=12.18

    \displaystyle Var(Y)=[e^{\sigma^2}-1] \ [  e^{ \mu+\frac{1}{2} \sigma^2} ]^2=255.02

Plug the first equation into the second equation and obtain the equation \displaystyle [e^{\sigma^2}-1] \ [12.18 ]^2=255.02. Solving for \sigma produces \sigma = 1. Plug \sigma = 1 into the first equation produces \mu = 2. The following gives the desired probability.

    \displaystyle P(Y > 10)=1-\Phi \biggl(\frac{\log(10)-2}{1} \biggr)=1-\Phi(0.30)=1-0.6179=0.3821

To find the skewness and kurtosis, one way is to find the first 4 lognormal moments and then calculate the third standardized moment (skewness) and the fourth standardized moment (kurtosis). To see how this is done, see this previous post. Another is to use the formulas given above.

    \gamma_1=(e^{\sigma^2}+2) \sqrt{e^{\sigma^2}-1}=(e^1+2) \sqrt{e^1-1}=6.1849

    \text{Kurtosis}=e^{4 \sigma^2}+2 e^{3 \sigma^2}+3 e^{2 \sigma^2}-3=e^{4}+2 e^{3}+3 e^{2}-3=113.9364

Example 3
Suppose that the lifetime (in years) of a certain type of machines follows the lognormal distribution described in Example 2. Suppose that you purchased such a machine that is 10-year old. What is the probability that it will last another 10 years?

This is a conditional probability since the machine already survived for 10 years already.

    \displaystyle \begin{aligned}P(Y > 20 \lvert Y > 10)&=\frac{1-\Phi (\frac{\log(20)-2}{1} )}{1-\Phi (\frac{\log(10)-2}{1} )} \\&=\frac{1-\Phi(1.0)}{1-\Phi(0.30)} \\&=\frac{1-0.8413}{1-0.6179} \\&=\frac{0.1587}{0.3821}=0.4153  \end{aligned}

______________________________________________________________________________________________________________________________
\copyright 2017 – Dan Ma

Practice Problem Set 2 – finding median losses

Posted on Updated on

This post has two practice problems to find the median of a distribution that models insurance losses.

Practice Problem 2a
Losses have a distribution with the following density function:

    \displaystyle f(x)=\frac{1}{6} e^{-\frac{x}{12}}-\frac{1}{6} e^{-\frac{x}{6}} \ \ \ x>0.

Calculate the median loss amount.

Practice Problem 2b
Losses are modeled by a distribution that is a mixture of two exponential distributions, one with mean 6 and the other with mean 12. The weight of each distribution is 50%. Calculate the median loss amount.

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

\text{ }

Solutions

Problem 2a

    Median \displaystyle =-12 \times \log \biggl(\frac{2-\sqrt{2}}{2} \biggr)=14.7354

    log is logarithm to base e = 2.718281828…

Problem 2b

    Median \displaystyle =-12 \times \log \biggl(\frac{-1+\sqrt{5}}{2} \biggr)=5.7745

    log is logarithm to base e = 2.718281828…

_______________________________________________________________________
\copyright 2017 – Dan Ma