Pricing equity-linked life insurance contracts with multiple risk factors by neural networks
PPricing equity-linked life insurance contracts with multiple risk factors by neural networks
Karim Barigou ∗ and (cid:32)Lukasz Delong † SGH Warsaw School of Economics, Collegium of Economic Analysis, Institute of Econometrics,Niepodleg(cid:32)lo´sci 162, Warsaw 02-554, Poland
Version: July 20, 2020
Abstract
This paper considers the pricing of equity-linked life insurance contracts with deathand survival benefits in a general model with multiple stochastic risk factors: interestrate, equity, volatility, unsystematic and systematic mortality. We price the equity-linked contracts by assuming that the insurer hedges the risks to reduce the local vari-ance of the net asset value process and requires a compensation for the non-hedgeablepart of the liability in the form of an instantaneous standard deviation risk margin. Theprice can then be expressed as the solution of a system of non-linear partial differentialequations. We reformulate the problem as a backward stochastic differential equationwith jumps and solve it numerically by the use of efficient neural networks. Sensitivityanalysis is performed with respect to initial parameters and an analysis of the accuracyof the approximation of the true price with our neural networks is provided.
Keywords : Equity-linked contracts; Neural networks; BSDEs with jumps; Stochas-tic mortality; Heston stochastic volatility; Hull-White stochastic interest rates. ∗ [email protected] (Corresponding author) † [email protected] a r X i v : . [ q -f i n . P R ] J u l Introduction
An equity-linked life insurance contract is a contract between a policyholder and an insurancecompany which presents insurance and investment features (Palmer (2006), Bauer et al.(2008)). Typically, the policyholder pays either a single premium or a stream of periodicpremiums during an accumulation phase. In return, the insurer guarantees a stream ofperiodic payments starting either immediately or at a future date. For investors, equity-linked life insurance contracts allow to participate in the equity market and offer higherreturns than fixed annuity contracts. Moreover, they offer protection against downside shocksin the financial market by the presence of diverse guarantees. In the context of increasinglife expectancies and the population facing the need for sustainable income, the demand forguaranteed income after retirement offered by equity-linked contracts is important (Haefeli(2013)).Pricing of equity-linked life insurance contracts is a classical theme in the actuarial lit-erature. Typically, this involves the choice of a financial and actuarial model to capturethe financial and insurance risks. Based on no-arbitrage arguments and using the standardassumption that the financial and insurance risks are independent, these two components aremerged together by a product of risk-neutral expectations (see Fung et al. (2014), Da Fon-seca & Ziveyi (2015), Ignatieva et al. (2016) among others). However, this approach does notaddress the fundamental issue of incompleteness arising from the insurance and the financialmarket. Since the risks cannot be perfectly hedged, the valuation is non-linear and we shouldproceed in two-steps by disentangling the hedgeable and non-hedgeable part of the liability(Pelsser & Stadje (2014), Dhaene et al. (2017) and Barigou et al. (2019)).We follow the pricing approach from Delong et al. (2019 b ) where the authors derived asystem of partial differential equations (PDEs) for the fair valuation of a stream of paymentscontingent on financial and insurance risk factors. In Delong et al. (2019 b ), three risk factorswere considered and numerical methods for solving the PDEs were not discussed. In practice,we would like to apply a valuation operator which takes into account multiple correlatedrisk factors. In this paper we deal with six stochastic financial and insurance risk factors:interest rate (two factors), equity, volatility, unsystematic and systematic mortality. Theinterest rate is modelled with two-factor Hull-White model, the equity is modelled withHeston model where the volatility is modelled with Cox-Ingersoll-Ross model, the lifetimesof policyholders are modelled, conditionally on the force of mortality, with independentexponential distributions and the force of mortality is modelled with Feller process. Weassume that equity-linked contracts include two types of minimum guarantees: death benefitsin case of the deaths of the policyholder and survival benefits in case the policyholder is stillalive at the contract termination. In such a framework, we can derive a system of PDEsfor the fair valuation of a portfolio of equity-linked contracts. Yet, solving such a system of2DEs seems to be challenging due to a large number of state variables and a large numberof PDEs which reflects the number of policies in the insurance portfolio.In recent years, some papers have considered the pricing of equity-linked contracts byneural networks. Doyle & Groendyke (2019) explored the use of neural networks to priceand hedge equity-linked contracts. They showed that neural networks offer an importantcomputational gain compared to crude Monte-Carle methods, however they did not considerstochastic mortality neither stochastic volatility in their framework. In Hejazi & Jackson(2016), the authors proposed a neural network approach to price equity-linked contractsfor large insurance portfolios. Their framework prices some representative contracts viaMonte-Carlo and uses a spatial interpolation method to price the whole portfolio.Another stream of literature investigated the use of neural networks to accelerate PDEsolvers. Han et al. (2018), Weinan et al. (2017), Beck et al. (2019) take advantage of re-inforcement learning to speed up solving high-dimensional partial differential equations byusing a so-called Deep BSDE method. Moreover, Chan-Wai-Nam et al. (2019) proposedan improved version of the Deep BSDE algorithm using different network architectures andparameterizations.We see three contributions of this paper. First, we consider the pricing of equity-linkedcontracts in a general incomplete financial-actuarial market with stochastic interest rate,equity, volatility and mortality. In the spirit of Delong et al. (2019 b ), we derive a systemof multi-dimensional non-linear PDEs. Second, using the representation of the PDEs as aBSDE with jumps, we show how to solve numerically the pricing problem by neural networks.In particular, we show how the jump component in the BSDE, which reflects the number ofPDEs in the system of equations and the number of policyholders in the portfolio, can beefficiently handled. Finally, we discuss the sensitivity of the price to model parameters andstudy the accuracy of our approach to standard pricing methods.The paper is structured as follows. In Section 2, we present in detail the financial andinsurance model. Section 3 derives the pricing PDEs, discusses their representation as aBSDE with jumps and the numerical implementation by neural networks. In Section 4, wepresent numerical results and Section 5 provides concluding remarks. Throughout the paper, we work on a probability space (Ω , F , P ) equipped with a filtration F = ( F t ) ≤ t ≤ T and a finite time horizon T < ∞ . We consider a general financial marketconsisting of a bank account, a bond and an equity. We use a two-factor interest rate modelto describe the short-rate process and Heston model to describe the dynamics of the equityprice and its volatility. The insurance risk is modeled by a jump process with stochastic3orce of mortality. We assume that the dynamics of the risk-free rate under an equivalent martingale measure Q is given by the G2++ two-additive-factor gaussian model: r ( t ) = ψ ( t ) + x ( t ) + y ( t ) , r (0) = r ,dx ( t ) = − ax ( t ) dt + σ x dW Q ( t ) , x (0) = 0 ,dy ( t ) = − by ( t ) dt + σ y dW Q ( t ) , y (0) = 0 , where ψ is a deterministic function such that the initial term structure predicted by themodel is fitted to the market term structure. The G2++ model is equivalent to two-factorHull-White model under appropriate re-parametrization (see e.g. Section 4.2.5 in Brigo &Mercurio (2007)). For hedging in incomplete markets, models need to be defined underthe real-world measure P . We follow the approach from Ait-Sahalia & Kimmel (2010) andintroduce the dynamics r ( t ) = ψ ( t ) + x ( t ) + y ( t ) , r (0) = r ,dx ( t ) = (cid:0) δ x σ x − ax ( t ) (cid:1) dt + σ x dW ( t ) , x (0) = 0 ,dy ( t ) = (cid:0) δ y σ y − by ( t ) (cid:1) dt + σ y dW ( t ) , y (0) = 0 , where we use constant risk premiums δ x , δ y to switch from the martingale measure to thereal-world measure.The dynamics of the risk-free bank account is given by dR ( t ) = R ( t ) r ( t, x ( t ) , y ( t )) dt, R (0) = 1 , where we emphasize that the risk-free rate r depends on ( t, x ( t ) , y ( t )). In our G2++ model,the dynamics of the bond price P ( t ) = E Q [ e − (cid:82) T ∗ t r ( s ) ds |F t ], which matures at time T ∗ ≥ T , isgiven by dP ( t ) = P ( t ) (cid:16)(cid:0) r ( t, x ( t ) , y ( t )) + ζ ( t ) (cid:1) dt + A ( t ) dW ( t ) + B ( t ) dW ( t ) (cid:17) , for some given deterministic functions A, B (see e.g. Section 4.2.3 in Brigo & Mercurio(2007)): A ( t ) = − σ x − e − a ( T ∗ − t ) a , B ( t ) = − σ y − e − b ( T ∗ − t ) b , and the risk premium ζ : ζ ( t ) = A ( t ) δ x + B ( t ) δ y . ζ should be positive, hence the risk premiums δ x , δ y should be negative(as always, we assume that σ x , σ y are positive).We assume that the equity follows the Heston stochastic volatility model (Heston (1993)): dS ( t ) = S ( t ) (cid:16)(cid:0) r ( t, x ( t ) , y ( t )) + γ (cid:112) v ( t ) (cid:1) dt + (cid:112) v ( t ) dW ( t ) (cid:17) , S (0) = 1 ,dv ( t ) = k ( η − v ( t )) dt + σ v (cid:112) v ( t ) dW ( t ) , v (0) = v , where we also assume a constant risk premium γ for the equity. We assume that v > k > kη ≥ σ v so that v is strictly positive.The policyholder pays a single premium for an equity-linked life insurance policy. Theinitial premium is invested in the policyholder’s account value and the policyholder’s fundsare invested in the bond and the equity. The policyholder’s account dynamics is then givenby dF ( t ) = uF ( t ) (cid:16)(cid:0) r ( t, x ( t ) , y ( t )) + ζ ( t ) (cid:1) dt + A ( t ) dW ( t ) + B ( t ) dW ( t ) (cid:17) +(1 − u ) F ( t ) (cid:16)(cid:0) r ( t, x ( t ) , y ( t )) + γ (cid:112) v ( t ) (cid:1) dt + (cid:112) v ( t ) dW ( t ) (cid:17) − cF ( t ) dt,F (0) = F , where F is the initial premium paid by the policyholder, u is a fixed percentage of theaccount invested in the bond, and c is a constant fee deducted by the insurer to cover thefinancial guarantees and insurance protection embedded in the contract. The portfolio ofequity-linked contracts, which we introduce in the next section, is homogeneous in the sensethat all policyholders invest the same single premium, follow the same investment strategyand have the same guarantees and insurance protection, i.e. F , u and c are equal andfixed for all policyholders. We point out that the value of the fee c results from pricing thecontract. In this paper we assume that the fee c has already been set by a pricing actuaryand we calculate the price of the contracts given c . However, we can also find c such thatthe price of the portfolio of equity-linked contracts is zero and we can call such a fee - a fairfee.Finally, all Brownian motions involved in the interest rate and the equity models arecorrelated with dW i ( t ) dW j ( t ) = ρ ij dt for i, j = 1 , . . . , The portfolio consists of n policies. All policyholders have the same age and are entitled totwo types of benefits: a death benefit D paid at the moment that the insured dies (providedhe dies in [0 , T ]) and a survival benefit S paid at terminal time T if the insured survives thattime. The benefits are contingent on the values of the policyholder’s account F , hence wewill use the notations D ( t, F ( t )) and S ( F ( T )) for death and survival benefits, respectively.5e assume that the lifetimes of the policyholders ( τ k ) k =1 ,...,n at policy inception are,conditional on the mortality intensity, independent and identically distributed with mortalityintensity λ , i.e. we assume that P (cid:0) τ k > t | ( λ ( s ) , ≤ s ≤ t ) (cid:1) = e − (cid:82) t λ ( s ) ds . The spot mortality intensity at calendar time t of a head aged x at time 0 is denoted by λ ( t ). Following Luciano et al. (2008) and Luciano et al. (2012), we assume that λ ( t ) followsthe Feller process without mean reversion: dλ ( t ) = qλ ( t ) dt + σ λ (cid:112) λ ( t ) dW ( t ) , λ (0) = λ , (2.1)where λ > q >
0. The mortality intensity can attain zero value but Luciano & Vigna(2008) showed that the probability of such an event is negligible for calibrated mortality data.The choice of this model is motivated by its parsimony-few parameters to be calibrated andits ability to fit cohort life tables due to the lack of mean reversion. Moreover, under thedynamics (2.1), the survival probability can be obtained in the closed-form: P (cid:0) τ k > t (cid:1) = E (cid:104) e − (cid:82) t λ ( s ) ds (cid:105) = e β ( t ) λ (0) , (2.2)with β ( t ) = − e bt b + b e bt , b = − (cid:112) q + 2 σ λ , b = b + q , b = b − q .The lives of the policyholders are dependent since they are all affected by the path of λ ( t )(systematic mortality risk). However, conditionally to the sample path followed by λ ( t ), thelives are independent (unsystematic mortality risk). Similar as in credit risk, we can computethe probability that two policyholders survive (see e.g. Section 9.5 in McNeil et al. (2005)).Using the conditional independence, we obtain: P ( τ i > T, τ j > T ) = E (cid:2) P (cid:0) τ i > T, τ j > T | ( λ ( s ) , ≤ s ≤ T ) (cid:1)(cid:3) = E (cid:2) P (cid:0) τ i > T | ( λ ( s ) , ≤ s ≤ T ) (cid:1) × P (cid:0) τ j > T | ( λ ( s ) , ≤ s ≤ T ) (cid:1)(cid:3) = E (cid:104) e − (cid:82) T λ ( s ) ds (cid:105) = e ˜ β ( T ) λ (0) , i (cid:54) = j, where ˜ β ( T ) = ( − e ˜ bT ) ˜ b +˜ b e ˜ bT , ˜ b = − (cid:112) q + 4 σ λ , ˜ b = ˜ b + q , ˜ b = ˜ b − q , which is an increasing functionof σ λ . Therefore, the higher the volatility of the mortality intensity, the higher the systematicrisk that both policyholders survive.We introduce the processes N and J which count the number of deaths and in-forcepolicies in the insurance portfolio, respectively: N ( t ) = n (cid:88) k =1 { τ k ≤ t } , N (0) = 0 ,J ( t ) = n − N ( t ) , J (0) = n.
6n addition, the compensated counting process:˜ N ( t ) = N ( t ) − (cid:90) t ( n − N ( s − )) λ ( s ) ds, (2.3)is a martingale, and it will be used in further derivations. Finally, in the remainder of thepaper, we assume that the processes which drive the insurance risk ( N, λ ) are independentof the financial market.The financial and insurance models considered in this paper are very popular in actuarialapplications. Let us remark that different models can also be used and the methods presentedbelow still hold.
Even though the insurer can invest in a bank account, a bond or an equity, the insuranceliabilities are not perfectly hedgeable due to the presence of the stochastic volatility risk andthe independent actuarial risks. Therefore, there is no unique pricing technique for insurancecontracts and one needs to disentangle the hedgeable and non-hedgeable part of the liability.For instance, Bayraktar et al. (2009) used the instantaneous Sharpe ratio and Liang & Lu(2017) considered the principle of equivalent utility. In this paper, we follow the approachof Delong et al. (2019 b ) and use the fair valuation operator with the instantaneous standarddeviation risk margin. Our method for pricing can be described in the following steps: We set up a hedging portfolio composed of a bond, an equity and a bank account. Tofind the optimal investments, we minimize the local variance (quadratic variation) ofthe net asset value process, which is defined as the excess value of the hedging portfolioover the fair price of the liability. In case of a complete market, the minimized localvariance is zero. However, the incompleteness of the market leads to a residual risk, asmeasured by the remaining local variance. The insurer asks a compensation for this residual non-hedgeable risk and we assumethat it is priced via an instanteneous standard deviation risk margin. The price can then be expressed as the solution to a system of non-linear PDEs. Thesystem of PDEs describes the prices given the number of in-force policies and thenumber of equations is equal to the number of policies in the portfolio, In order to solve these PDEs, we reformulate the PDEs as a forward-backward stochas-tic differential equation with jumps and approximate the price using deep neural net-works. 7 .1 Fair pricing with the instantaneous standard deviation riskmargin
We assume that the insurer can invest in the bank account, the bond and the equity in orderto hedge the financial guarantees and the benefits embedded in the contracts. Let θ denotethe amount of money invested in the bond, and θ denote the amount of money investedin the equity. The remaining amount is invested in the bank account. We assume that θ and θ are Markov control strategies and are functions of ( t, x ( t ) , y ( t ) , F ( t ) , v ( t ) , λ ( t ) , J ( t )).Let V θ := ( V θ ( t ) , ≤ t ≤ T ) denote the self-financing hedging portfolio under a strategy θ = ( θ , θ ). The dynamics of the insurer’s hedging portfolio is given by dV θ ( t ) = θ ( t ) (cid:16)(cid:0) r ( t, x ( t ) , y ( t )) + ζ ( t ) (cid:1) dt + A ( t ) dW ( t ) + B ( t ) dW ( t ) (cid:1) + θ ( t ) (cid:16)(cid:0) r ( t, x ( t ) , y ( t )) + γ (cid:112) v ( t ) (cid:1) dt + (cid:112) v ( t ) dW ( t ) (cid:1) +( V θ ( t ) − θ ( t ) − θ ( t )) r ( t, x ( t ) , y ( t )) dt + J ( t − ) cF ( t ) dt − D ( t, F ( t )) dN ( t ) ,V θ (0) = nF (0) . Moreover, let us remark that the survival benefits J ( T ) S ( F ( T )) are paid at the terminaltime T and therefore, do not affect the insurer’s hedging portfolio.In order to simplify the derivations, we introduce the 5-dimensional Itˆo process Z ( t ) =( x ( t ) , y ( t ) , F ( t ) , v ( t ) , λ ( t )) which models the five risk factors: dZ ( t ) = µ ( t, Z ( t )) dt + σ ( t, Z ( t )) dW ( t ) , (3.1)with µ ( t, Z ( t )) = δ x σ x − ax ( t ) δ y σ y − by ( t ) (cid:16) r ( t, x ( t ) , y ( t )) − c + uζ ( t ) + (1 − u ) γ (cid:112) v ( t ) (cid:17) F ( t ) k ( η − v ( t )) qλ ( t ) , (3.2)and σ ( t, Z ( t )) = σ x . . . . . . σ y . . . uF ( t ) A ( t ) uF ( t ) B ( t ) (1 − u ) F ( t ) (cid:112) v ( t ) 0 00 . . . σ v (cid:112) v ( t ) 00 . . . . . . σ λ (cid:112) λ ( t ) . (3.3)In the expression (3.1), W ( t ) is a 5-dimensional correlated Brownian motion with correlationmatrix Q , i.e. dW i ( t ) dW j ( t ) = ρ ij dt where ρ ij is the ( i, j )-component of Q (for i, j =1 , . . . , ϕ := ( ϕ ( t ) , ≤ t ≤ T ) denote the continuous-time fair valuation operator represent-ing the fair price of the liabilities. We assume that ϕ is a function of ( t, x ( t ) , y ( t ) , F ( t ) , v ( t ) , λ ( t ) , J ( t ))and is denoted by ϕ J ( t ) ( t, Z ( t )). We introduce the net asset value process X ( t ) = ϕ J ( t ) ( t, Z ( t )) − V θ ( t ) representing the difference between the price of the liabilities and the value of the in-surer’s hedging portfolio. Let Θ( s, Z ( s )) = ( θ ( s ) A ( s ) , θ ( s ) B ( s ) , θ ( s ) (cid:112) v ( s ) , ,
0) denote a5-dimensional vector related to the hedging strategy. By application of the multidimensionalItˆo’s lemma, we find the dynamics: dX ( s ) = (cid:110) X ( s ) r ( s, x ( s ) , y ( s )) + ϕ J ( s − ) t ( s, Z ( s )) + ∇ ϕ J ( s − ) ( s, Z ( s )) · µ ( s, Z ( s )) − cJ ( s − ) F ( s ) − ϕ J ( s − ) ( s, Z ( s )) r ( s, x ( s ) , y ( s )) + (cid:16) ϕ J ( s − ) − ( s, Z ( s )) + D ( s, F ( s )) − ϕ J ( s − ) ( s, Z ( s )) (cid:17) J ( s − ) λ ( s ) − θ ( s ) ζ ( s ) − θ ( s ) γ (cid:112) v ( s ) (cid:111) ds + 12 Tr (cid:16) σ ( s, Z ( s )) Q σ (cid:124) ( s, Z ( s )) Hess z ϕ J ( s − ) ( s, Z ( s )) (cid:17) ds + (cid:16) ϕ J ( s − ) − ( s, Z ( s )) + D ( s, F ( s )) − ϕ J ( s − ) ( s, Z ( s )) (cid:17) d ˜ N ( s )+ (cid:16) ∇ ϕ J ( s − ) ( s, Z ( s )) σ ( s, Z ( s )) − Θ( s, Z ( s )) (cid:17) dW ( s ) . Here, ∇ ϕ and Hess z ϕ denote the gradient and the Hessian of the function ϕ with respectto Z , Tr denotes the trace of a matrix.We are interested in finding the valuation operator ϕ which satisfies the property:lim h → E (cid:104) X ( t + h ) e − (cid:82) t + ht r ( s ) ds − X ( t ) (cid:12)(cid:12)(cid:12) Z ( t ) = z, J ( t ) = k (cid:105) h + RM (cid:104) X ( t + h ) e − (cid:82) t + ht r ( s ) ds − X ( t ) (cid:12)(cid:12)(cid:12) Z ( t ) = z, J ( t ) = k (cid:105) h = 0 . (3.4)As the one-period risk margin, we consider the standard deviation: RM (cid:104) X ( t + h ) e − (cid:82) t + ht r ( s ) ds − X ( t ) (cid:12)(cid:12)(cid:12) Z ( t ) = z, J ( t ) = k (cid:105) = 12 α (cid:114) V ar (cid:104) X ( t + h ) e − (cid:82) t + ht r ( s ) ds − X ( t ) (cid:12)(cid:12)(cid:12) Z ( t ) = z, J ( t ) = k (cid:105) √ h. The parameter α can be interpreted as the insurer’s risk aversion coefficient towards thenon-hedgeable risks.We can choose the hedging strategy θ so that the insurer does not charge a risk marginfor the hedgeable financial risk. The optimal investments in the bond and the equity are thendetermined by minimizing the instantaneous (local) variance of the net asset value process.Therefore, we need to compute the quadratic variation [ X ] s and find the bivariate process9 θ ( s ) , θ ( s )) which minimizes this quantity. The quadratic variation is given by[ X ] T = (cid:90) T (cid:16) ∇ ϕ ( s, Z ( s )) σ ( s, Z ( s )) − Θ( s, Z ( s )) (cid:17) Q (cid:16) ∇ ϕ ( s, Z ( s )) σ ( s, Z ( s )) − Θ( s, Z ( s )) (cid:17) (cid:124) (cid:124) (cid:123)(cid:122) (cid:125) g ( s,θ ( s ) ,θ ( s )) ds + (cid:16) ϕ J ( s − ) − ( s, Z ( s )) + D ( s, F ( s )) − ϕ J ( s − ) ( s, Z ( s )) (cid:17) dN ( s ) . (3.5)The optimal hedging strategy is then obtained by solving the first-order conditions: ∂g ( s, θ , θ ) ∂θ = − (cid:88) j =1 (cid:16) ρ j A ( s ) + ρ j B ( s ) (cid:17)(cid:16) ∇ ϕ ( s, Z ( s )) σ ( s, Z ( s )) − Θ( s, Z ( s )) (cid:17) j = 0 ,∂g ( s, θ , θ ) ∂θ = − (cid:88) j =1 ρ j (cid:112) v ( s ) (cid:16) ∇ ϕ ( s, Z ( s )) σ ( s, Z ( s )) − Θ( s, Z ( s )) (cid:17) j = 0 . We remark that the index j = 5 does not appear in the derivatives since the mortalityintensity is not correlated with the financial risk factors. The first-order conditions yield asystem of two equations with two unknowns for which the solutions are given by θ ∗ ( t, Z ( t )) = uϕ f ( t, Z ( t )) F ( t )+ ϕ x ( t, Z ( t )) σ x ρ x ( t ) + ϕ y ( t, Z ( t )) σ y ρ y ( t ) + ϕ v ( t, Z ( t )) σ v (cid:112) v ( t ) ρ v ( t ) A ( t ) ρ x ( t ) + B ( t ) ρ y ( t ) ,θ ∗ ( t, Z ( t )) = (1 − u ) ϕ f ( t, Z ( t )) F ( t )+ ϕ x ( t, Z ( t )) σ x (cid:112) v ( t ) (cid:18) ρ − ( ρ A ( t ) + ρ B ( t )) ρ x ( t ) A ( t ) ρ x ( t ) + B ( t ) ρ y ( t ) (cid:19) + ϕ y ( t, Z ( t )) σ y (cid:112) v ( t ) (cid:18) ρ − ( ρ A ( t ) + ρ B ( t )) ρ y ( t ) A ( t ) ρ x ( t ) + B ( t ) ρ y ( t ) (cid:19) + ϕ v ( t, Z ( t )) σ v (cid:18) ρ − ( ρ A ( t ) + ρ B ( t )) ρ v ( t ) A ( t ) ρ x ( t ) + B ( t ) ρ y ( t ) (cid:19) , with some deterministic functions, which take into account the correlation between thedifferent risk factors: ρ x ( t ) = A ( t ) + ρ B ( t ) − ρ ρ A ( t ) − ρ ρ B ( t ) ,ρ y ( t ) = B ( t ) + ρ A ( t ) − ρ ρ A ( t ) − ρ ρ B ( t ) ,ρ v ( t ) = ρ A ( t ) + ρ B ( t ) − ρ ρ A ( t ) − ρ ρ B ( t ) . As it is usually the case, it turns out that the optimal hedging strategy appears as a delta-hedging strategy involving the sensitivities of the liabilities to the underlying financial riskfactors and adjusted for the dependence structure.10lugging the optimal parameters into the quadratic variation (3.5) leads to[ X ] ∗ T = (cid:90) T (cid:16) ∇ ϕ ( s, Z ( s )) σ ( s, Z ( s )) − Θ ∗ ( s, Z ( s )) (cid:17) Q (cid:16) ∇ ϕ ( s, Z ( s )) σ ( s, Z ( s )) − Θ ∗ ( s, Z ( s )) (cid:17) (cid:124) ds + (cid:16) ϕ J ( s − ) − ( s, Z ( s )) + D ( s, F ( s )) − ϕ J ( s − ) ( s, Z ( s )) (cid:17) dN ( s ) , which represents the remaining quadratic variation. We remark that the vector (cid:16) ∇ ϕσ − Θ ∗ (cid:17) can be written as (cid:16) ∇ ϕσ ∗ (cid:17) for an approriate matrix σ ∗ . Similarly, the risk-adjusted drift µ ∗ can be obtained from ∇ ϕ J ( s − ) ( s, Z ( s )) · µ ∗ ( s, Z ( s )) = ∇ ϕ J ( s − ) ( s, Z ( s )) · µ ( s, Z ( s )) − θ ∗ ( s ) ζ ( s, x ( s ) , y ( s )) − θ ∗ ( s ) γ (cid:112) v ( s ) . If we consider the gradient component related to f , we find that ϕ J ( s − ) f ( s, Z ( s )) µ ∗ ,f ( s, Z ( s )) = ϕ J ( s − ) f ( s, Z ( s )) (cid:104)(cid:16) r ( s, x ( s ) , y ( s )) − c + uζ ( s ) + (1 − u ) γ (cid:112) v ( s ) (cid:17) F ( s ) (cid:105) − ϕ J ( s − ) f ( s, Z ( s )) uF ( s ) ζ ( s ) − ϕ J ( s − ) f ( s, Z ( s ))(1 − u ) F ( s ) γ (cid:112) v ( s )= ϕ J ( s − ) f ( s, Z ( s )) ( r ( s, x ( s ) , y ( s )) − c ) . Hence, the risk-adjusted drift of the fund F ( t ) is the risk-free rate r ( t, x ( t ) , y ( t )) minus thefee rate c , and the equity risk-premium is cancelled.Under the optimal hedging strategy, we find thatlim h → E (cid:2) X ( t + h ) e − r ( t ) h − X ( t ) | Z ( t ) = z, J ( t ) = k (cid:3) h = ϕ kt ( t, z ) + ∇ ϕ k ( t, z ) · µ ∗ ( t, z ) + 12 Tr (cid:16) σ ( t, z ) Q σ (cid:124) ( t, z ) Hess z ϕ k ( t, z ) (cid:17) − ckf + (cid:16) ϕ k − ( t, z ) + D ( t, f ) − ϕ k ( t, z ) (cid:17) kλ − ϕ k ( t, z ) r ( t, x, y ) , (3.6)andlim h → Var (cid:2) X ( t + h ) e − r ( t ) h − X ( t ) | Z ( t ) = z, J ( t ) = k (cid:3) h = (cid:16) ∇ ϕ k ( t, z ) σ ∗ ( t, z ) (cid:17) Q (cid:16) ∇ ϕ k ( t, z ) σ ∗ ( t, z ) (cid:17) (cid:124) + (cid:16) ϕ k − ( t, z ) + D ( t, f ) − ϕ k ( t, z ) (cid:17) kλ. (3.7)The last quantity represents the remaining quadratic variation and quantifies the non-hedgeable risks of the net asset value process X . We can conclude that the continuous-time Given the optimal parameters, it is easy to decompose Θ ∗ as Θ ∗ = ∇ ϕ ˜ σ . Then the appropriate matrixis σ ∗ = σ − ˜ σ since ∇ ϕ ( σ − ˜ σ ) = ∇ ϕσ − Θ ∗ . ϕ solves the following system of non-linear PDEs: ϕ kt ( t, z ) + ∇ ϕ k ( t, z ) · µ ∗ ( t, z ) + 12 Tr (cid:16) σ ( t, z ) Q σ ( t, z ) (cid:124) Hess z ϕ k ( t, z ) (cid:17) − ckf + (cid:16) ϕ k − ( t, z ) + D ( t, f ) − ϕ k ( t, z ) (cid:17) kλ ( t ) − ϕ k ( t, z ) r ( t, x, y )+ α (cid:114)(cid:16) ∇ ϕ k ( t, z ) σ ∗ ( t, z ) (cid:17) Q (cid:16) ∇ ϕ k ( t, z ) σ ∗ ( t, z ) (cid:17) (cid:124) + (cid:16) ϕ k − ( t, z ) + D ( t, f ) − ϕ k ( t, z ) (cid:17) kλ = 0 ,ϕ k ( T, z ) = kS ( f ) , (3.8)for k ∈ { , . . . , n } . For each k , we have to deal with a 5-dimensional PDE. The number ofPDEs in the system (3.8) is equal to n + 1. It is well known that there is a close relation between partial differential equations (PDEs)and backward stochastic differential equations (BSDEs). The solutions ( ϕ k ( t, z )) k =0 ,...,n tothe PDEs (3.8) satisfy the following forward-backward stochastic differential equation (FB-SDE) with jumps (see e.g. Barles et al. (1997) or Chapter 4 in Delong (2013)): ϕ J ( t ) ( t, Z ( t )) = J ( T ) S ( F ( T ))+ (cid:90) Tt Υ J ( s − ) (cid:16) s, Z ( s ) , ϕ J ( s − ) ( s, Z ( s )) , ϕ J ( s − ) − ( s, Z ( s )) (cid:17) ds + (cid:90) Tt D ( s, F ( s )) J ( s − ) λ ( s ) ds − (cid:90) Tt (cid:16) cJ ( s − ) F ( s ) + ϕ J ( s − ) ( s, Z ( s )) r ( s, x ( s ) , y ( s )) (cid:17) ds − (cid:90) Tt ∇ ϕ J ( s − ) ( s, Z ( s )) σ ( s, Z ( s )) dW ( s ) − (cid:90) Tt (cid:16) ϕ J ( s − ) − ( s, Z ( s )) − ϕ J ( s − ) ( s, Z ( s )) (cid:17) d ˜ N ( s ) , (3.9) Z ( t ) = Z (0) + (cid:90) t µ ∗ ( s, Z ( s )) ds + (cid:90) t σ ( s, Z ( s )) dW ( s ) ,dJ ( t ) = − dN ( t ) , (3.10)whereΥ J ( s − ) (cid:16) s, Z ( s ) , ϕ J ( s − ) ( s, Z ( s )) , ϕ J ( s − ) − ( s, Z ( s )) (cid:17) = α (cid:110)(cid:16) ∇ ϕ J ( s − ) ( s, Z ( s )) σ ∗ ( s, Z ( s )) (cid:17) Q (cid:16) ∇ ϕ J ( s − ) ( s, Z ( s )) σ ∗ ( s, Z ( s )) (cid:17) (cid:124) + (cid:16) ϕ J ( s − ) − ( s, Z ( s )) + D ( s, F ( s )) − ϕ J ( s − ) ( s, Z ( s )) (cid:17) J ( s − ) λ ( s ) (cid:111) / . (3.11)12he equation (3.9) is a one-dimensional BSDE with jumps where the generator and theterminal condition depend on a process satisfying a forward stochastic differential equation.In our case, the FSDE is a 6-dimensional equation given by (3.10), where five componentsare diffusions and one component is driven by a step process. The BSDE (3.9) can also bewritten as a FSDE: ϕ J ( t ) ( t, Z ( t )) = ϕ n (0 , Z (0)) − (cid:90) t Υ J ( s − ) (cid:16) s, Z ( s ) , ϕ J ( s − ) ( s, Z ( s )) , ϕ J ( s − ) − ( s, Z ( s )) (cid:17) ds − (cid:90) t D ( s, F ( s )) J ( s − ) λ ( s ) ds + (cid:90) t (cid:16) cJ ( s − ) F ( s ) + ϕ J ( s − ) ( s, Z ( s )) r ( s, x ( s ) , y ( s )) (cid:17) ds + (cid:90) t ∇ ϕ J ( s − ) ( s, Z ( s )) σ ( s, Z ( s )) dW ( s )+ (cid:90) t (cid:16) ϕ J ( s − ) − ( s, Z ( s )) − ϕ J ( s − ) ( s, Z ( s )) (cid:17) d ˜ N ( s ) . (3.12)Starting from an initial price of the liabilities: ϕ n (0 , Z (0)) (i.e. for n policyholders at time t = 0 with initial values of the risk factors: z = Z (0)), the FSDE (3.12) describes the pricedynamics of the insurance liabilities through time.The function Υ in (3.12) represents the price for the non-hedgeable risks: the first termin (3.11) is the residual risk from the diffusion component (interest rate, volatility, equityand systematic mortality risks) and the second term in (3.11) is the residual risk from thejump component (unsystematic mortality risk). We can notice that the risk margin for thenon-hedgeable financial risks takes into account the correlation matrix between these risks.The second term in (3.12) represents the price for the expected death benefits. The firsttwo terms in (3.12) decrease the price as they are released over time. The next two termsincrease the price. The third term stands for the continuous fees paid by the policyholders tocover the guarantees and benefits (since the fees to hedge the future liabilities are collectedby the insurer, the price of the future liabilities must increase over time). The fourth termdescribes the effect of interest accrual and the time value of money on the price dynamics.The fifth term shows the impact of the diffusion component. Finally, the last term depictsthe consequence of the death of a policyholder on the price (there is one policyholder less J ( s − ) (cid:55)→ J ( s ) = J ( s − ) − ϕ n (0 , Z (0)),we treat ϕ n (0 , Z (0)) as parameters for the moment and view the FSDE (3.12) as a way ofcomputing the values of ϕ J ( t ) ( t, Z ( t )) knowing the initial price ϕ n (0 , Z (0)) and the dynamicsof ϕ . We consider a temporal discretization for the risk factors and the FSDE (3.12). Fora set of time steps 0 = t < t < . . . < t N = T , we consider the simple Euler scheme for13 = 1 , . . . , N − Z ( t n +1 ) − Z ( t n ) ≈ µ ∗ ( t n , Z ( t n ))∆ t n + σ ( t n , Z ( t n ))∆ W ( t n ) , for the diffusion process (3.1) and˜ N ( t n +1 ) − ˜ N ( t n ) = N ( t n +1 ) − N ( t n ) − J ( t n ) λ ( t n )∆ t n , for the compensated jump process (2.3). For the FSDE, the Euler scheme has the followingform: ϕ J ( t n +1 ) ( t n +1 , Z ( t n +1 )) = ϕ J ( t n ) ( t n , Z ( t n )) − α (cid:110)(cid:16) ∇ ϕ J ( t n ) ( t n , Z ( t n )) σ ∗ ( t n , Z ( t n )) (cid:17) Q (cid:16) ∇ ϕ J ( t n ) ( t n , Z ( t n )) σ ∗ ( t n , Z ( t n )) (cid:17) (cid:124) + (cid:16) ϕ J ( t n ) − ( t n , Z ( t n )) + D ( t n , F ( t n )) − ϕ J ( t n ) ( t n , Z ( t n )) (cid:17) J ( t n ) λ ( t n ) (cid:111) / ∆ t n − D ( t n , F ( t n )) J ( t n ) λ ( t n )∆ t n + (cid:16) cJ ( t n ) F ( t n ) + ϕ J ( t n ) ( t n , Z ( t n )) r ( t n , x ( t n ) , y ( t n )) (cid:17) ∆ t n + ∇ ϕ J ( t n ) ( t n , Z ( t n )) σ ( t n , Z ( t n ))∆ W ( t n )+ (cid:16) ϕ J ( t n ) − ( t n , Z ( t n )) − ϕ J ( t n ) ( t n , Z ( t n )) (cid:17) ∆ ˜ N ( t n ) (3.13)where∆ t n = t n +1 − t n , ∆ W ( t n ) = W ( t n +1 ) − W ( t n ) , ∆ ˜ N ( t n ) = ˜ N ( t n +1 ) − ˜ N ( t n ) . To move from time t n to time t n +1 , there are two unknown functions: the gradient ofthe price, ∇ ϕ J ( t n ) ( t n , Z ( t n )), and the price with one policyholder less, ϕ J ( t n ) − ( t n , Z ( t n )).Following the work of Weinan et al. (2017) and Chan-Wai-Nam et al. (2019) (we remarkthat in both cases, the authors only consider BSDEs without jumps), we approximate bothfunctions by the use of neural networks. Hereafter, we briefly introduce feed-forward neural networks (short neural networks ). Inregression theory, networks can be seen as a broad generalization of generalized linear models(GLMs) where the linear predictor is replaced by a non-linear one.Let us assume that the input is in dimension d (the state variable x ) and the output isin dimension d (the number of value functions to estimate). The network is characterizedby a number of layers L + 1 ∈ N \{ , } with m (cid:96) , (cid:96) = 0 , . . . , L, the number of neurons (unitsor nodes) on each layer: the first layer is the input layer with m = d , the last layer is theoutput layer with m L = d , and the L − m (cid:96) = m, (cid:96) = 1 , . . . , L −
1. A network is a functionfrom R d to R d defined as the composition x ∈ R d (cid:55)→ N ( x ) = A L ◦ (cid:37) ◦ A L − ◦ . . . ◦ (cid:37) ◦ A ( x ) ∈ R d Here, A (cid:96) , (cid:96) = 1 , . . . , L , are affine transformations represented by A (cid:96) ( x ) = W (cid:96) x + β (cid:96) , for a matrix of weights W (cid:96) and a bias term β (cid:96) . The non-linear function (cid:37) : R → R iscalled the activation function and is applied component-wise on the outputs of A (cid:96) , i.e., (cid:37) ( x , . . . , x m ) = ( (cid:37) ( x ) , . . . , (cid:37) ( x m )). Standard examples of activation functions are thesigmoid, the ReLU, the elu and tanh. For a general introduction on neural networks, werefer to Goodfellow et al. (2016).The universal approximation theorem of Hornik et al. (1989) states that networks canapproximate any continuous function on a compact support arbitrarily well if we allow forarbitrarily many neurons q ∈ N in the hidden layer.In the Euler scheme (3.13) of the FBSDE, we approximate the gradient and the pricewith one less policyholder by two respective neural networks. Moreover, in practice, pricesneed to be computed several times for different input parameters. Therefore, as suggestedin Weinan et al. (2017), we determine the price at time 0 in a region of input parametersinstead of a single space-point by including a third neural network. Overall, we have thethree following networks as depicted in Figure 1: N1.
The price at time 0, ϕ J (0) (0 , Z (0)), is approximated by a neural network N : ( z, k ) ∈ R (cid:55)→ ϕ k (0 , z ) ∈ R with parameters φ : ϕ J (0) (0 , Z (0)) = N φ ( Z (0) , J (0)) . N2.
The gradient, ∇ ϕ J ( t n ) ( t n , Z ( t n )), is approximated by a second neural network N :( t, z, k ) ∈ R (cid:55)→ ∇ ϕ k ( t, z ) ∈ R with parameters χ : ∇ ϕ J ( t n ) ( t n , Z ( t n )) = N χ ( t n , Z ( t n ) , J ( t n )) . N3.
The price with one less policyholder, ϕ J ( t n ) − ( t n , Z ( t n )), is approximated by a thirdneural network N : ( t, z, k − ∈ R (cid:55)→ ϕ k − ( t, z ) ∈ R with parameters ψ : ϕ J ( t n ) − ( t n , Z ( t n )) = N ψ ( t n , Z ( t n ) , J ( t n ) − . By the use of the Euler scheme (3.13) and the neural networks, we can obtain a Monte-Carlo approximation (cid:98) ϕ of ϕ J ( T ) ( T, Z ( T )) which depends on the parameters of the neural15igure 1: Illustration of the network architecture for solving the system of non-linear PDEswith N time steps. The network N1 stands for the price at time 0, the network N2 approxi-mates the gradient and the network N3 is related to the jump component and approximatesthe price with one policyholder less.networks ( φ, χ, ψ ). These parameters are then estimated in order to minimize the quadraticloss function: l ( φ, χ, ψ ) = E (cid:20)(cid:12)(cid:12)(cid:12) J ( t N ) S ( F ( t N )) − (cid:98) ϕ J ( t N ) ( t N , Z ( t N )) (cid:12)(cid:12)(cid:12) (cid:21) . (3.14)The parameters can now be optimized by stochastic gradient descent-type (SGD) algorithmas it is usually done to train neural networks. In this paper, we use the Adaptive MomentEstimation (Adam) optimizer proposed by Kingma & Ba (2014) and the algorithm is im-plemented with the Keras R package (Chollet et al. (2017)), which is a user-friendly API toTensorFlow. Remark . • One can notice either from the PDEs (3.8) or the FSDE (3.12) that weface a recursive problem: in order to determine the price of the insurance liabilitiesfor n initial policyholders, we need to know the price for k ∈ { , . . . , n − } . Onepossibility of implementation would have been to solve n recursive optimizations withone neural network for the gradient in each optimization. Instead, we consider onlyone optimization with two neural networks (one for the gradient, one for the jumpcomponent) in order to reduce the computational cost of the problem. • Contrary to Weinan et al. (2017), we consider one single network for all time stepsinstead of one neural network per time step in order to reduce the numbers of param-eters to be estimated. Chan-Wai-Nam et al. (2019) showed that this neural networkarchitecture significantly improves the precision of the results and the stability of thealgorithms. 16
Numerical results
This section presents numerical results of our neural networks to price equity-linked contractswith death and survival benefits. In Section 4.1, we briefly discuss our input parameters andthe hyperparameters of the neural networks. In Section 4.2, we solve the system of PDEs(3.8) with our neural networks and analyse the impact of various model parameters on theprice of equity-linked contracts. In Section 4.3, we assess the accuracy of our approach forsome specific contracts. Let us point out that we only investigate the initial price of thecontracts, i.e. at time t = 0. For the financial and insurance model, we consider a realistic set of model parameters whichare consistent with previous studies: • For the two-factor interest rate model, we consider a = 0 . , σ x = 0 . , b = 0 . , σ y = 0 . , which are based on Russo & Torri (2019). Moreover, we chose: δ x = δ y = − . , ψ ( t ) =0 . . • For the Heston stochastic volatility model, the parameters are κ = 0 . , η = 0 . , σ v = 0 . , γ = 0 . , which are consistent with Li et al. (2016). • For the force of mortality dynamics, we choose the parameters from Luciano & Vigna(2008) which correspond to UK individuals aged 65 at time 0: q = 0 . , σ λ = 0 . . • The correlation matrix for the interest rate factors, the equity, the volatility and themortality intensity is given by Q = − . .
35 0 0 − . .
08 0 00 .
35 0 .
08 1 − . − . (4.1)which is consistent with Grzelak et al. (2011).17dditional input parameters are presented in Table 1. We point out that the pricingneural network should not only deliver one price but a surface of prices which is a functionof the interest rate factors, the premium, the volatility, the number of policyholders and theforce of mortality at time t = 0. Hence, we are interested in determining the price denoted by ϕ J (0) (0 , Z (0)). Here, Z (0) = ( x (0) , y (0) , F (0) , v (0) , λ (0)) stands for the set of initial values forthe interest rate factors, the premium, the volatility and the force of mortality respectively,and J (0) stands for the number of policyholders at the inception. In our numerical examplewe investigate ranges of input parameters by changing the values for v (0) , F (0) , λ (0) , J (0)specified in Table 1 by + / − x (0) ∈ [0 , . , y (0) ∈ [0 , . c x (0) 0 D (cid:63) y (0) 0 S (cid:63) F (0) 1 T v (0) 0 . t λ (0) 0 . α J (0) 100Table 1: The base values of the initial parameters.For the neural network architecture, different hyperparameters need to be chosen. Eachof the neural networks N1 , N2 , N3 defined in Section 3.3 consists of two hidden layers of20 neurons. For the activation function, we considered the ELU (exponential linear unit)which tends to accelerate the learning process (Clevert et al. (2015)). Moreover, neuralnetworks tend to have convergence issues when their inputs are not scaled and centered,and perform best when the inputs are in the range ∼ [ − , g : R (cid:55)→ [ − ,
1] : g ( x ) = 2 × (cid:18) x − min( x, x, − min( x, (cid:19) − . Finally, we chose a batch size and a number of epochs equal to 200. We considered 10000sample paths for each risk factor and the default learning rate of the Adam optimizer. Wetried different hyperparameters and we did not observed significant improvements in trainingof the neural networks. The computation time for the training of the price surface takes about 10 min on a computer with anIntel Core i7-4970MQ processor running at 2.90 GHz with 16 GB of RAM. .2 Pricing in the general case and sensitivity analysis Each equity-linked life insurance contract provides insurance death and survival benefitslinked to the policyholder’s account value. These financial guarantees ensure that the benefitsare not lower than the minimal guarantee levels stipulated at the inception of the contract.We set the benefits: D ( t, F ( t )) = ( D ∗ − F ( t )) + , S ( F ( T )) = ( S ∗ − F ( T )) + . The contract payoffs if the policyholder’s account value goes below D ∗ at the time of deathor below S ∗ in case of survival at the end of the contract. Hence, we are interested in thefair valuation of the financial guarantees embedded in the insurance contracts. First, we consider the convergence of our algorithm for one price with the initial values Z (0) = (0 , , , . , . J (0) = 100 policyholders from Table 1. Figure 2 representsthe convergence of the price and the mean square error against the number of epochs used fortraining the neural networks. Let us recall that our objective in the reinforcement learningis to minimize the mean square error (3.14), hence we present the mean square error lossfunction to validate the convergence. We observe that, with the initialization of our algorithmby default, the price reaches a plateau after around 100 epochs. After 100 epochs, the MSEcontinues to slightly decreases but without significant impact on the price. We arrive at theprice of the portfolio of the equity-linked contracts equal to ϕ J (0) (0 , Z (0)) = 6 . Z (0) and J (0). The sets were constructed by taking the valuesgiven in Table 1 and changing v (0) , F (0) , λ (0) , J (0) by + / −
25% and considering x (0) ∈ [0 , . , y (0) ∈ [0 , . Hereafter, we illustrate the impact of various model parameters on the price of the portfolioof the equity-linked contracts. One of the purpose of this exercise is to validate that our19
246 0 50 100 150 200
Number of epochs P r i ce Convergence of the price by epoch
Number of epochs M SE Accuracy measurement by epoch
Figure 2: Diagnostics of the neural network training. On the left: Convergence of the priceat time 0 per epoch. On the right: Mean Square Error per epoch.PDEs provide a reasonable pricing operator.Figure 4 represents the effect of the initial interest rate factors x (0) and y (0) on the price.As expected, the price decreases as the interest rate factors increase due to the time valueof money. In particular, if both factors increase from 0 to 0 .
01, the price decreases by 16%.Figure 5 shows the effect of varying the initial premium F (0) and the number of policy-holders J (0) on the price. Obviously, the price of the portfolio of the equity-linked contractsincreases with the number of policyholders since more policyholders mean more benefits topay. Moreover, we note that the price of the portfolio of the equity-linked contracts decreaseswith the initial premium. Indeed, the increase of the initial fund value decreases the prob-ability that the fund will go below the death guarantee level D ∗ and the survival guaranteelevel S ∗ , and therefore decreases the price of the contract. We can observe that the pricecan double if the initial premium decreases by 10% as depicted in Figure 5.Figure 6 depicts the impact of varying the initial volatility parameter v (0) and the ini-tial force of mortality λ (0) on the price of the portfolio of the equity-linked contracts. Asexpected, we note that the price is an increasing function of the initial volatility. A highlyvolatile market environment implies that the stock and the fund values are less predictable,which will result in a higher risk charge required by the insurer in order to cover the guar-antees. When the initial force of mortality increases, the value of the death benefit increaseswhile the value of the survival benefit decreases as policyholders are more likely to die. OnFigure 6, the price is a fairly constant function of the initial force of mortality, showing20 Number of epochs P r i ce Convergence of the price by epoch
Number of epochs M SE Accuracy measurement by epoch
Figure 3: On the left: Convergence of the price at time 0 during the neural network trainingfor five randomly chosen initial parameters. On the right: Mean Square Error per epochduring the neural network training of the surface of prices.that the value reduction of the survival benefit compensates the value increase of the deathbenefit, at least for D ∗ = S ∗ = 1 .
02 as in our numerical example.We also investigate the impact of the correlation between the fund and the interest ratefactor x as well as the correlation between the fund and the volatility on the price of theliabilities. The results are presented in Table 2. As expected, the lower the correlationcoefficient, the lower the price for the portfolio of the equity-linked contracts. This propertyreflects diversification effect between the financial risks and is consistent with the findingsreported in Gudkov et al. (2019) and Dai et al. (2015).Change of the correlation Change of the price ρ xf ρ fv +0 . − . ± . α . We expect that the more riskaverse the insurer is, the larger the risk margin he will keep to cover the non-hedgeable risks21 .00000.00250.00500.00750.0100 0.0000 0.0025 0.0050 0.0075 0.0100 Interest rate factor x(0) I n t e r es t r a t e f ac t o r y ( ) -16.0%-12.0%-8.0%-4.0%0.0% Price Change
Price of the equity-linked contract
Figure 4: The price of the portfolio of equity-linked contracts for a range of the initial interestrate factors x (0) ∈ [0; 0 .
01] and y (0) ∈ [0; 0 . α is not linear.The case α = 0 corresponds to the best estimate of the liability, which we discuss in the nextsubsection, when the insurer does not set any risk margin for the non-hedgeable risks. Hereafter, we assess the accuracy of our neural network approach by comparing it withMonte-Carlo pricing for two particular cases.
In insurance, the Best Estimate Liability (BEL) plays a central role in the valuation ofinsurance liabilities. It is defined as the expected present value of the future cashflowsdiscounted by the risk-free rate. For our portfolio of the equity-linked life insurance contracts,the BEL is given by
BEL = E (cid:20) e − (cid:82) Tt r ( τ ) dτ J ( T ) S ( F ( T )) + (cid:90) Tt e − (cid:82) st r ( τ ) dτ D ( s, F ( s )) J ( s ) λ ( s ) ds − (cid:90) Tt e − (cid:82) st r ( τ ) dτ cJ ( s ) F ( s ) (cid:12)(cid:12)(cid:12) J ( t ) , x ( t ) , y ( t ) , F ( t ) , v ( t ) , λ ( t ) (cid:21) (4.2)where the first term is the expected survival benefits, the second term - the expected deathbenefits and the third term - the expected fees. By the Feynman-Kac formula, the BEL is22 .80.91.01.11.2 80 90 100 110 120 Policyholders I n i t i a l p r e m i u m -100%0%100%200%300% Price Change
Price of the equity-linked contract
Figure 5: The price of the portfolio of equity-linked contracts for a range of the initialpremiums F (0) ∈ [0 .
75; 1 .
25] and the numbers of policyholders J (0) ∈ [75; 125].also solution to the following system of PDEs: ϕ kt ( t, z ) + ∇ ϕ k ( t, z ) · µ ∗ ( t, z ) + 12 Tr (cid:16) σ ( t, z ) Q σ ( t, z ) (cid:124) Hess z ϕ k ( t, z ) (cid:17) − ckf + (cid:16) ϕ k − ( t, z ) + D ( t, f ) − ϕ k ( t, z ) (cid:17) kλ ( t ) − ϕ k ( t, z ) r ( t, x, y ) = 0 ,ϕ k ( T, z ) = kS ( f ) , (4.3)for k ∈ { , . . . , n } .The system (4.3) corresponds to our pricing equation (3.8) without the risk margin (with α = 0). To assess the accuracy of our neural network, we compare the price obtained by anapproximation of the BEL (4.2) by 200.000 Monte-Carlo simulations and the price obtainedby our neural network using the same parameters as in Table 1. Figure 8 shows that thereis a good convergence of the neural network to the Monte-Carlo price. After 100 epochs, weobtain a relative error of c.a. 0 . α > α = 0. Ourtests showed that such an approach does not speed up the calibration of our neural networksmainly due to the fact that we still have to initiate the neural network for the gradient froma random starting point. 23 .080.090.100.110.120.011 0.013 0.015 0.017 0.019 Initial force of mortality I n i t i a l v o l a t ili t y -10.0%-5.0%0.0%5.0%10.0% Price Change
Price of the equity-linked contract
Figure 6: The price of the portfolio of equity-linked contracts for a range of the initialvolatilities v (0) ∈ [0 . . λ (0) ∈ [0 . . In the following example, we consider a portfolio of equity-linked contracts in which only asurvival benefit is guaranteed. Such contracts are often called Guaranteed Minimum Ma-turity Benefits (GMMB). We consider a simplified financial-actuarial market, namely theBlack-Scholes financial market model and a jump process J ( t ) with constant mortality λ tomodel the deaths of policyholders. Under these assumptions, the pricing PDEs (3.8) reduceto ϕ kt ( t, f ) + ϕ kf ( t, f ) f r + 12 ϕ kff ( t, f ) f σ f − ϕ k ( t, f ) r + (cid:16) ϕ k − ( t, f ) − ϕ k ( t, f ) (cid:17) kλ (cid:18) − α √ kλ (cid:19) = 0 ,ϕ k ( T, f ) = kS ( f ) , (4.4)for k ∈ { , . . . , n } .We remark that this PDE was considered in (Delong et al. (2019 a ), Example 4.5). Thefirst part of the PDE (4.4) exactly corresponds to the standard Black-Scholes PDE butthe second part adds a risk margin for the uncertainty about the number of survivals whichcannot be completely hedged. Applying the Feynman-Kac formula, the solution of the PDEs(4.4) is given by ϕ n (0 , F (0)) = E ˜ Q (cid:2) e − rT J ( T ) S ( F ( T )) (cid:3) , (4.5)24 .006.056.106.15 0.00 0.05 0.10 0.15 Risk aversion parameter α P r i ce Price of the equity-linked contract
Figure 7: The price of the portfolio of equity-linked contracts in terms of the risk aversionparameter α .where under ˜ Q , the jump process has now a stochastic mortality intensity given by˜ λ ( t ) = λ (cid:32) − α (cid:112) J ( t ) λ (cid:33) . We choose the parameters: r = 0 . , F (0) = 1 , σ f = 0 . , S ∗ = 1 . , T = 1 , n = 100 , α =0 .
1. Under these parameters, ˜ λ ( t ) is strictly positive which allows us to use a step processto simulate the number of deaths under the measure ˜ Q .Figure 9 compares the Monte-Carlo price obtained by approximating (4.5) with 200.000simulations and the price obtained by neural networks by solving the PDEs (4.4). We clearlyobserve a fast convergence after only 25 epochs and a relative error of c.a. 0 .
2% after 100epochs.Summing up, both examples show that our algorithm is accurate. The remaining error isdue to discretization steps applied in reinforcement learning, when we look for the optimalstrategy in a discrete time, and Monte-Carlo simulations of discretized stochastic processes.Moreover, we note that our remaining MSE is of the order [10 − ; 10 − ] which is consistentwith orders obtained in Chan-Wai-Nam et al. (2019). The second example in this sectionalso shows that our approach to handling the jump component performs very well.25
246 0 50 100 150 200
Number of epochs P r i ce NN priceMC price
Convergence of the NN price
Number of epochs M SE Accuracy measurement by epoch
Figure 8: Neural Network price convergence to the Monte-Carlo price for the BEL and theMean Square Error per epoch
In this paper, we have derived a partial differential equation for fair pricing of equity-linkedlife insurance contracts in a general financial-actuarial market with stochastic interest rate,equity, volatility and mortality. We assumed that the insurer can hedge its liabilities witha bank account, a bond and an equity and requires a compensation for non-hedgeable risksin the form of an instantaneous standard deviation risk margin. We derived a system ofPDEs which described the fair valuation operator. We used the connection with BSDEswith jumps and proposed an efficient neural network architecture to solve the PDEs. Thenetwork consists essentially of two subnetworks: one to approximate the gradient componentand one to approximate the jump component. We also included a third subnetwork for theprice at time 0 which allows us to derive the price for a range of initial parameters ratherthan just one price.Numerical results provided convergence diagnostics, price sensitivities to various modelparameters and some comparisons with Monte-Carlo methods. In particular, we observeda fast convergence and stabilization of the price after a few dozens of epochs and a goodapproximation compared to prices obtained by Monte-Carlo simulations. In terms of sensi-tivity, we proved that the calibrated price agreed with intuition.Our paper provides a stable computational algorithm for fair pricing of insurance liabil-ities in incomplete markets which can be used for a wide variety of financial and insurance26
Number of epochs P r i ce NN priceMC price
Convergence of the NN price
Number of epochs M SE Accuracy measurement by epoch
Figure 9: Neural Network price convergence to the Monte-Carlo price for the GMMB andMean Square Error per epoch.products. Moreover, our neural network offers potential other applications in fields wherea system of non-linear PDEs naturally arises, such as quantum mechanics, game theory,dynamic programming and many more.
References
Ait-Sahalia, Y. & Kimmel, R. L. (2010), ‘Estimating affine multifactor term structure modelsusing closed-form likelihood expansions’,
Journal of Financial Economics (1), 113–144.Barigou, K., Chen, Z. & Dhaene, J. (2019), ‘Fair dynamic valuation of insurance liabilities:Merging actuarial judgement with market-and time-consistency’, Insurance: Mathematicsand Economics , 19–29.Barles, G., Buckdahn, R. & Pardoux, E. (1997), ‘Backward stochastic differential equa-tions and integral-partial differential equations’, Stochastics: An International Journal ofProbability and Stochastic Processes (1-2), 57–83.Bauer, D., Kling, A. & Russ, J. (2008), ‘A universal pricing framework for guaranteedminimum benefits in variable annuities’, Astin Bulletin (02), 621–651.Bayraktar, E., Milevsky, M. A., Promislow, S. D. & Young, V. R. (2009), ‘Valuation of27ortality risk via the instantaneous sharpe ratio: applications to life annuities’, Journalof Economic Dynamics and Control (3), 676–691.Beck, C., Weinan, E. & Jentzen, A. (2019), ‘Machine learning approximation algorithms forhigh-dimensional fully nonlinear partial differential equations and second-order backwardstochastic differential equations’, Journal of Nonlinear Science (4), 1563–1619.Brigo, D. & Mercurio, F. (2007), Interest rate models-theory and practice: with smile, infla-tion and credit , Springer Science & Business Media.Chan-Wai-Nam, Q., Mikael, J. & Warin, X. (2019), ‘Machine learning for semi linear pdes’,
Journal of Scientific Computing (3), 1667–1712.Chollet, F., Allaire, J. et al. (2017), ‘R interface to keras’, https://github.com/rstudio/keras .Clevert, D.-A., Unterthiner, T. & Hochreiter, S. (2015), ‘Fast and accurate deep networklearning by exponential linear units (elus)’, arXiv preprint arXiv:1511.07289 .Da Fonseca, J. & Ziveyi, J. (2015), ‘Valuing variable annuity guarantees on multiple assets’, Scandinavian Actuarial Journal pp. 1–22.Dai, T.-S., Yang, S. S. & Liu, L.-C. (2015), ‘Pricing guaranteed minimum/lifetime withdrawalbenefits with various provisions under investment, interest rate and mortality risks’,
In-surance: Mathematics and Economics , 364–379.Delong, (cid:32)L. (2013), Backward stochastic differential equations with jumps and their actuarialand financial applications , Springer.Delong, (cid:32)L., Dhaene, J. & Barigou, K. (2019 a ), ‘Fair valuation of insurance liability cash-flow streams in continuous time: Applications’, ASTIN Bulletin: The Journal of the IAA (2), 299–333.Delong, (cid:32)L., Dhaene, J. & Barigou, K. (2019 b ), ‘Fair valuation of insurance liability cash-flowstreams in continuous time: Theory’, Insurance: Mathematics and Economics .Dhaene, J., Stassen, B., Barigou, K., Linders, D. & Chen, Z. (2017), ‘Fair valuation ofinsurance liabilities: Merging actuarial judgement and market-consistency’,
Insurance:Mathematics and Economics , 14–27.Doyle, D. & Groendyke, C. (2019), ‘Using neural networks to price and hedge variableannuity guarantees’, Risks (1), 1. 28ung, M. C., Ignatieva, K. & Sherris, M. (2014), ‘Systematic mortality risk: An analysisof guaranteed lifetime withdrawal benefits in variable annuities’, Insurance: Mathematicsand Economics , 103–115.Goodfellow, I., Bengio, Y. & Courville, A. (2016), Deep learning , MIT press.Grzelak, L. A., Oosterlee, C. W. & Van Weeren, S. (2011), ‘The affine heston model withcorrelated gaussian interest rates for pricing hybrid derivatives’,
Quantitative Finance (11), 1647–1663.Gudkov, N., Ignatieva, K. & Ziveyi, J. (2019), ‘Pricing of guaranteed minimum with-drawal benefits in variable annuities under stochastic volatility, stochastic interest ratesand stochastic mortality via the componentwise splitting method’, Quantitative Finance (3), 501–518.Haefeli, D. (2013), ‘Variable annuities-an analysis of financial stability’, The Geneva Asso-ciation .Han, J., Jentzen, A. & Weinan, E. (2018), ‘Solving high-dimensional partial differential equa-tions using deep learning’,
Proceedings of the National Academy of Sciences (34), 8505–8510.Hejazi, S. A. & Jackson, K. R. (2016), ‘A neural network approach to efficient valuation oflarge portfolios of variable annuities’,
Insurance: Mathematics and Economics , 169–181.Heston, S. L. (1993), ‘A closed-form solution for options with stochastic volatility withapplications to bond and currency options’, The review of financial studies (2), 327–343.Hornik, K., Stinchcombe, M., White, H. et al. (1989), ‘Multilayer feedforward networks areuniversal approximators.’, Neural networks (5), 359–366.Ignatieva, K., Song, A. & Ziveyi, J. (2016), ‘Pricing and hedging of guaranteed minimumbenefits under regime-switching and stochastic mortality’, Insurance: Mathematics andEconomics , 286–300.Kingma, D. P. & Ba, J. (2014), ‘Adam: A method for stochastic optimization’, arXiv preprintarXiv:1412.6980 .Li, D., Rong, X. & Zhao, H. (2016), ‘Optimal reinsurance and investment problem for aninsurer and a reinsurer with jump-diffusion risk process under the heston model’, Compu-tational and Applied Mathematics (2), 533–557.29iang, X. & Lu, Y. (2017), ‘Indifference pricing of a life insurance portfolio with risky assetdriven by a shot-noise process’, Insurance: Mathematics and Economics , 119–132.Luciano, E., Regis, L. & Vigna, E. (2012), ‘Delta–gamma hedging of mortality and interestrate risk’, Insurance: Mathematics and Economics (3), 402–412.Luciano, E., Spreeuw, J. & Vigna, E. (2008), ‘Modelling stochastic mortality for dependentlives’, Insurance: Mathematics and Economics (2), 234–244.Luciano, E. & Vigna, E. (2008), ‘Mortality risk via affine stochastic intensities: calibrationand empirical relevance’, Belgian Actuarial Bulletin (1), 5–16.McNeil, A. J., Frey, R. & Embrechts, P. (2005), Quantitative Risk Management: Concepts,Techniques and Tools-revised edition , Princeton university press.Palmer, B. (2006), ‘Equity-indexed annuities: Fundamental concepts and issues’,
Work .Pelsser, A. & Stadje, M. (2014), ‘Time-consistent and market-consistent evaluations’,
Math-ematical Finance: An International Journal of Mathematics, Statistics and Financial Eco-nomics (1), 25–65.Russo, V. & Torri, G. (2019), ‘Calibration of one-factor and two-factor hull–white modelsusing swaptions’, Computational Management Science (1-2), 275–295.Weinan, E., Han, J. & Jentzen, A. (2017), ‘Deep learning-based numerical methods for high-dimensional parabolic partial differential equations and backward stochastic differentialequations’, Communications in Mathematics and Statistics5