Mortality and Life Expectancy Forecast for (Comparatively) High - - PDF document

Mortality and Life Expectancy Forecast for (Comparatively) High Mortality Countries

Ahbab Mohammad Fazle Rabbi∗ , Stefano Mazzuco† Department of Statistical Sciences, University of Padua, Italy September 27, 2017

Abstract The Lee-Carter method and its later variants are widely accepted extrapolative methods for forecasting mortality and life expectancy in industrial countries due to their simplicity and availability of long time series data with better quality. We compared and contrast the mortality and life expectancy forecast using 7 different variants of Lee-Carter methods for 9 comparatively high mortality countries than industrialized countries along with UN forecast

• method. The models provided under estimated forecast for high mortality countries compare

to comparatively better mortality patterns. Diverging mortality pattern was observed in these countries which is also reflected in future mortality and life expectancy pattern. Better fit of different models in country specific mortality pattern was observed while no models gives better fit uniquely. In the same context, use of probabilistic forecasting technique from Bayesian framework provided better forecast than some of the extrapolative methods. Country specific forecast indicates better fit of certain extrapolative methods may occur in part of the life span rather than that of the whole life span. These findings imply necessity for invention of new forecast method in context of high mortality countries. Keywords: Mortality forecast; Lee-Carter method; Probabilistc forecast; Mortality in Eastern Europe; Bayesian Hirerchical model

1 Introduction

Mortality improvements are observed almost all over the world during the 20th century and it is always considered as a positive change for socio-economical advancement of a country. This change has brought new requirements in support-systems for the elderly, such as health-care and pension provision. Aging became the greatest population problem for many industrialized countries since the 1970s (Brouhns et al, 2002). This has resulted in a surge of interest among government policy makers and planners in accurately modeling and forecasting age-specific mortality rates. Policy makers rely greatly on future mortality rates as an indicator of future population structure also. However, besides the developed countries, population aging are also visible in some other countries of the world for which mortality pattern is not same as the the industrialized countries. This sort of different mortality regime may be visible for several Eastern European countries, Russia and some other countries where causes of deaths may play vital role on age-specific deaths (Vallin & Mesl´ e, 2001). Several parametric (both in Frequentist and Bayesian point of view) and nonparametric methods have been proposed over the years for forecasting age-specific mortality rates and

∗fazlerabbi@stat.unipd.it †mazzuco@stat.unipd.it

1

life expectancy (for instance, Shang, 2012). The simplest way for parametric forecasting is to parametrize the available series of life table and hence to extrapolate each of the parame- ter separately for obtain the forecast from the assumed model (Keyfitz, 1991). Without any doubt, the ground breaking model based on extrapolation was proposed by Lee and Carter (Lee & Carter, 1992). The advantages of the Lee-Carter (LC) method includes its simplicity and robustness in situations where age-specific log mortality rates have linear trends (Booth et al, 2002). LC method attempts to capture the patterns of mortality rates using only one principal component and its scores. However, LC method does not give good results in pres- ence of irregular mortality schedule (Lee & Miller, 2001). The question of whether or how the length of the fitting period affects point forecast accuracy was first identified by Booth et al (2002). These problems made the later studies to impose restriction to start fitting the models from 1950 and following years and other modifications. Also several of studies proposed and compared point and interval forecast accuracy for forecasting age-specific mortality rates and life expectancy at birth (say, Shang, 2012). All of these studies concluded that use of different assumptions leads to different outcomes, and comparing different variants and extensions does not ensure identification of a single best method for all of the countries (Shang, 2012). The reference period used is the main determinant of large differences in outcomes (Pollard 1987), especially when there is considerable non-linearity in the trends or unusual mortality pattern in certain age groups (Stoeldraijer et al, 2013). Thus, the comparison of outcomes from different studies is hampered by differences in the explicit assumptions; such as the choice of the length of the historical period and of the jump-off rates in particular age groups (Stoeldraijer et al, 2013). All these applications and comparisons handled mostly for the industrialized countries. These countries have low mortality regime with high life expectancies; lower adult or early senescence mortality (Shang, 2012). All of the countries studied before for comparing the per- formance of the mortality forecasting models have stable pattern of mortality transition over long period of time; irregular pattern in different age-specific mortality were not the case for any of these countries in the course of demographic transition (except some years off-course). Irregular pattern in age-specific mortality effect forecast of mortality greatly. Distribution of deaths play vital role on modeling mortality, population with lower modal age at deaths sup- posed to have different mortality profile which also have impact on forecast (even for short period). The previous studies also enjoyed the advantages of data quality; most of these coun- tries for illustration have long time series of mortality data with very high quality. For instance, Sweden have series of life tables since 1751 with high quality life tables were available since 1860 (Human Mortality Database, 2016). The available extrapolative methods require certain period

• f mortality data for fitting the model; data quality of certain time periods often restricts the

fitting period. So, performance of these existing mortality models are still subject to analyze in this comparatively high mortality countries. Completely different mortality regime may be observed for several Eastern European coun- tries rather than other developed countries of Europe (Vallin & Mesl´ e, 2001). Several countries

• f these part of the world have higher mortality regime than the neighbor Western European
• r other industrialized countries. Accidental mortality and presence of high frequency of cause-

specific mortality are common in this region (Monostori et al, 2015). Many of these countries showed increasing trend in life expectancies over the decades, however, irregular trend was also

• bserved. Several of these countries showed decreasing pattern of life expectancy in last century

for couple of times (Vallin & Mesl´ e, 2001; Monostori et al, 2015). Another barrier for mod- eling or forecasting mortality for these countries is data quality. Human Mortality Database mentioned country specific caution notes to use the data of some countries during some certain time periods due to data source (Human Mortality Database, 2016). 2

Therefore, we suggest forecasting models should be tested with data characterized by i) higher mortality regime; ii) limited time frame for fitting period of the models and iii) rel- atively lower quality. We encounter these criteria are common for many under-developed or developing countries. In this paper we filled the gap on this line of research on comparison of mortality forecasting models; we compared and contrast the performance of the mortality (and life expectancy) forecasting models for 9 of these (comparatively) high mortality countries; Be- larus, Bulgaria, Estonia, Hungary, Latvia, Lithuania, Russia, Slovakia and Ukraine. All of these countries have some unique characteristic on mortality regime. Diverging regime in longevity is observed in this part of the world whereas other industrialized countries showed converging pattern (Monostori et al, 2015). They also have similar mortality trend over the years with similar pattern of cause-specific mortality. Moreover, the socio-economic background of these countries are also same. Such comparison may help the researchers to understand adequately the further scope of developing mortality forecasting models; also to the policymakers for better policy implication regarding age and cause-specific mortality. There are four parts in this paper. In section 2, we reviewed the forecasting models we used in this study and the source of data. Section 3 has three subsections. In first one, we compared the findings of the models for selected countries with high mortality regime than that of other industrialized countries considered in previous studies. In second subsection we forecast the life expectancy using Bayesian Hierarchical Model and in last subsection we illustrated the country specific results for Hungary to have more insight of the fitted models. The concluding remarks are placed in section 5.

2 Data and Methodology

The data used in this study came from Human Mortality Database (Human Mortality Database, 2016). Data from 9 countries of Eastern Europe, Baltic states and Russia are utilized in this

• study. This part of Europe (and Russia) has higher mortality compare to Western Europe.

Another interesting characteristics about these countries is that mortality levels are diverging for these countries whereas other developed countries in Western Europe are converging. Also, former Socialist countries have some similarities in social inequalities (Mackenbach, 2013). We excluded Poland and Slovenia from this analysis as they have better mortality scenario than these countries. Details of the fitting period of the models are given in Results section (Table- 1). For the UN forecast, the data of life expectancy at birth are taken from Human Mortality Database on 5 years basis rather than complete life tables used for the extrapolative models. We also contrast this method using the data of World Population Prospect 2012 (United Nations, 2013). The models used in this study are reviewed below in brief.

2.1 Lee-Carter Method (1992)

Since its development, Lee-Carter (LC) method is one of the most applicable methods till now. Use of principal components for mortality forecasting came to practice through the work of Lee and Carter (1992). In contrast to parametric approaches which postulate the functional form

• f the age pattern of mortality beforehand; principal components approaches estimate the age

pattern from the data (life tables). Forecasting is then done by extrapolation of the time related parameters using time series methods. The two-factor LC model is given below, ln mx,t = ax + bxkt + ǫx,t (1) Here, mx,t is the central mortality rate at age x for year t; ax represents the average of log- mortality at age x over time; bx is the first principal component capturing relative change in the log mortality rate at each age x; ktrepresents the overall level of mortality in year t; and 3

ex,t is the model residual. The LC model adjusts kt by refitting to the total number of deaths. The adjusted kt are hence extrapolated by a random walk model (with drift), from which point forecasts are obtained by (1) with the fixed ax and bx. The restrictions of the model are;

p

• x=0

bx = 1 and,

n

• t=1

kt = 0 The main advantage of LC method was it’s simplicity to explain and due to less number of parameters, it does not suffer the complicated scenario of over-parametrization (Keyfitz, 1991). A shortcoming of LC model is that it assumes that the ratio of the rates of mortality change at different ages remains constant over time, which is not always the case in real life situation (Booth et al, 2002). Girosi and King (2006) mentioned another problem of LC model in real life problem that the forecast rates lack across-age smoothness and become increasingly jagged

• ver time.

2.2 Lee-Miller method (2001)

An important development on LC model was done by Lee and Miller in the beginning of 2000s. Lee and Miller (2001) proposed three modifications on the basic LC model; (i) the fitting period was restricted to the latter half of the twentieth century to reduce structural shifts, (ii) adjustment of kt was done by matching life expectancy, and (iii) ‘jump-off error’was eliminated by forecasting forward from observed (rather than fitted) rates. The Lee-Miller (LM) variant has now been widely adopted as the standard Lee-Carter method. All the other variants proposed later are also based on this restriction criterion; due to severe heterogeneity in age specific mortality on earlier and later ages. However, still one can run the basic LC model with the old assumptions as well.

2.3 Booth, Maindonald and Smith (2002) method

Lee-Carter model can be extended to include higher order terms also. Higher order terms were modeled by Booth, Maindonald and Smith (2002) and forecast was later developed by using univariate ARIMA processes(Renshaw & Haberman, 2003). The extended version may be written as, ln mx,t = ax + b1

xk1 t + b2 xk2 t + . . . + bp xkp t ǫx,t

(2) The key modifications of this Booth, Maindonald and Smith (BMS) model are; (i) instead of Lee-Miller’s late 1950 referral period, the fitting period is determined by a statistical ‘goodness

• f fit’ criterion, under the assumption that kt is linear; and (ii) the adjustment of kt involves

fitting to the age distribution of deaths rather than to the total number of deaths in the basic LC model. This model is a significant development in the research of forecasting mortality as it slightly eliminates a shortcoming of LC Model. LC model assumes invariant bx whereas evidence of substantial age-time interaction is common (Shang, 2012). Clearly, this model partially addresses this problem by choosing a fitting period that optimally satisfies the model assumptions (Booth et al, 2002).

2.4 Lee-Carter method with Poisson Regression

Although several estimation procedure is later suggested by later studies, this approach was the first one to assume underlying distribution of deaths. Brouhns et al (2002) proposed the following for fitting the Lee-Carter method with ordinary least square (OLS) estimate of the model parameters, Dx,t ∼ Poisson {Ex,tmx,t} with mx,t = exp(ax + bxkt) 4

Here Ex,t are the population exposed for death at age X in time t. The constraints of the basic LC model also holds for this method. One of the main advantages of using this Poisson Regression approach is, this method allows to have maximum likelihood estimation of the model parameters instead of OLS or Gauss-Newton algorithm (Brouhns at al, 2002). This also shows some further development scopes to using Bayesian approach on LC methods later. In this paper, we will mainly analyze the classical approaches instead of Bayesian framework.

2.5 Nonparametric Approaces: Hyndman-Ullah (2007) method

Like older forecasting methods, LC method lacks across-age smoothness and heterogeneity of deaths over long time period (Girosi & King, 2006), which is non-intuitive and may be prob- lematic in many practical applications. Another weakness of the basic LC variants was that they attempt to capture the patterns of age-specific mortality rates using only the first principal component and its associated scores (Booth et al, 2002). To address this problem, Hyndman and Ullah (2007) proposed a functional data model that utilizes second and higher order prin- cipal components to capture additional variation in mortality rates. The method proposed by Hyndman and Ullah (2007) combines the nonparametric penalized regression Spline with func- tional principal component analysis for forecasting mortality rates. Main advantage of using Hyndman and Ullah (HU) method is it allows to model higher order terms and strengthens the estimation basis to allow outliers. There are two more benefits of this approach, now more than

• ne principal component is used and it gives flexibility for applying a wide range of univariate

time series models to forecast the principal component scores. The special features of this model are given below.

• 1. a penalized regression spline with partial monotonic constraint is utilized first to smooth

the log mortality rates. The following continuous smooth function ft(x) is assumed for discrete ages. mt(xi) = ft(xi) + σt(xi)ǫt,i; i = 1, . . . , p; t = 1, . . . , n (3) Here mt(xi) represents log-transformed mortality rates for each age xi in time t; σt(xi) is the noise component and ǫt,i is i.i.d. standard normal variable.

• 2. HU method allows more than one principal component unlike previously mentioned vari-

ants of LC method. Functional principal component analysis utilizes a set of continuous functions and is decomposed into functional principal components and their associated

• scores. Symbolically,

ft(x) = a(x) +

J

• j=1

bj(x)kt,j + et(x); t = 1, . . . , n Here a(x) is the mean function

• = 1

n

n

t=1 ft(x)

• ; bj(x) are set of first J functional prin-

cipal components; kt,j are set of uncorrelated principal component scores; et(x) is the residual function. It should be noted that J < n is considered for optimal number of functional principal components.

• 3. Almost every suitable univariate time series models may be utilized for principal compo-

nent scores. However, ARIMA model is suggested to have minimum AIC of the fitted model (Shang, 2012). Life tables constructed from the smoothed ft(x) consists of lower variance than that from the original mx which provides better estimates of life expectancy as well. Hyndman & Ullah 5

(2007) proposed use of weighted penalized regression splines for estimating ft(x). This weight- ing controls heterogeneity due to σt(x) and a monotonic constraint for upper ages can lead to better estimates. In this study, we applied weights equal to the approximate inverse variances wx,t = mx,tEx,t, and used weighted penalized regression splines to estimate the curve ft(x) in each year (Hyndman & Ullah, 2007). Weighted penalized regression splines are preferable due to computational time and it allows monotonicity constraints to be imposed relatively easily (Hyndman & Ullah, 2007). The details estimation procedure of interval forecast are given else- where (Hyndman & Ullah, 2007). Two more version of HU method were also proposed for special situations (Hyndman & Ullah, 2007). Those are mentioned below in brief. 2.5.1 Robust Hyndman and Ullah (2007) method Hyndman and Ullah (2007) proposed this method to forecast in presence of outliers. The Robust Hyndman and Ullah method investigates the integrated squared error for each year by calculating following measures of accuracy for the functional principal component approximation

• f the functional data.
• xp

x−1

 ft(x) − a(x) −

J

• j=1

bj(x)kt,j  

2

dx After assigning zero weight to outliers, the Robust Hyndman and Ullah method fit the mortality rates from which forecasts of age-specific life expectancies can be estimated without affect of prospective outliers. 2.5.2 Weighted Hyndman and Ullah (2007) method On of the above mentioned HU methods have one shortcoming; it doesn’t give very good fit for data which is not from very recent past. To overcome this problem, Hyndman and Ullah (2007) proposed another weighted version of the method alongside with the above mentioned

• methods. The new method can be showed symbolically as follows,

ft(x) = a∗(x) +

J

• j=1

b∗

j(x)kt,j + et(x)

(4) Here, a∗(x) is the weighted functional mean such as,

• a∗(x) =

n

• t=1

wtft(x),

n

• t=1

wt = 1, where, wt = κ(1 − κ)n−t; t = 1, . . . , n This wt is the new weights defined by Hyndman and Ullah (2007) for 0 < κ < 1, a geomet- rically decaying weight parameter. The optimal value of κ is chosen by minimizing an overall forecast error measure within the validation data set among a set of possible candidates. Details

• f the methods can be found elsewhere (Hyndman & Shang, 2009).

2.6 Bayesian framework: UN mortality and life expectacny forecast

One of the major shortcomings of all the above mentioned variants of LC model is that these methods require age-specific death rates for at least three time periods, which is not the case for many developing countries (Raftery et al, 2013). Raftery and colleagues (2013) proposed 6

an alternative approach for forecasting using Bayesian framework. They applied Bayesian Hi- erarchical Model to forecast period life expectancy directly, using a random walk model with a non-constant drift. The newly defined drift term is a nonlinear function of current life ex- pectancy and reflects the fact that life expectancy which tends to changes more slowly for the countries with the lowest and highest life expectancies, and more quickly for the countries in the middle. The UN produces estimates of age-specific mortality and period life expectancy at birth for all member countries and updates in every two years in UN World Population Prospects (United Nations 2013). The UN projects life expectancy in the next time period deterministically using the equation ℓc,t+1 = ℓc,t + g(ℓc,t) (5) Forecast of life expectancy is done by a double-logistic function of the current level of life

• expectancy. Symbolically,

g(ℓc,t|θc) = kc 1 + exp

• − A1

∆c

2 (ℓct − ∆c

1 − A2∆c 2)

• (6)

+ zc − kc 1 + exp

• − A1

∆c

4

• ℓct − 3

i=1 ∆c i − A2∆c 4

• (7)

These six parameters (∆c

1, ∆c 2, ∆c 3, ∆c 4, kc, zc) are the six parameters of the double logistic

function for country c at time t. The estimation policy changed since World Population Prospect 2012 (WPP 2012) as the UN Population Division used a probabilistic model for the first time to forecast life expectancy at birth using the methods of Raftery and colleagues (2013). Raftery and colleagues (2013) used the following hierarchical model to turned the old deterministic model into probabilistic one (with uncertainty) and hence adopted a Bayesian approach to estimate the model parameters. Hence the hierarchical model become, ℓc,t+1 = ℓc,t + g

• ℓc,t|θ(c)

+ ǫc,t+1 (8) Where, g

• ℓc,t|θ(c)

= Double-Logistic function with parameters θc θc =∆c

1, ∆c 2, ∆c 3, ∆c 4, kc, zc

∆i|σ∆i

iid

∼ Normal[0,100]

• ∆i, σ2

∆i

• ,

i = 1, . . . , 4 kc|σk

iid

∼ Normal[0,10]

• k, σ2

k

• zc|σz

iid

∼ Normal[0,1.15]

• z, σ2

z

• This hierarchical model allows pooling information about the rates of change across coun-

tries with assumption that each set of country-specific double-logistic parameters are randomly sampled from a common truncated normal distribution. Raftery and colleagues (2013) defined proper prior for all 13 parameters of the model in such a way that the prior distributions are more defuse than the posterior distributions. Thus, the above mentioned hierarchical model turned into a Bayesian Hierarchical Model. This probabilistic approach is adopted by UN for life expectancy forecast which eventually replaces the previous deterministic method; UN also adopted probabilistic methods for forecasting fertility rates or population projection as well. This method has one advantage over any other parametric or non-parametric methods: it is flexible on choosing prior to get fast, slow or medium pace for change in life expectancy level; as per usual forecasting of UN. It should be noted that this method was designed for working with data of World Population Prospects or similar formats. One of the major disadvantage is 7

that, it does not forecast considering whole life tables like the previous methods or Bayesian approach on LC method; which made it complicated to compare the outcomes with LC variants. And, as this method was innovated with aim to replace the old deterministic method of UN forecast, this method still works for five (calender) year based data. It takes single value of life expectancy for each five years and also return the forecast as a median of five (calender) years period utilizing a Bayesian Hierarchical Model.

2.7 Comparison of the existing methods

There were several studies where the point and interval forecasts were compared for developed and low mortality countries only (Hyndman & Ullah, 2007; Shang, 2012). Among the recent works, Shang and colleagues (2011) compared between major principal components methods

• nly. Another study compared all the variants of LC methods and a newly proposed weighted

version of all the LC methods (Shang, 2012). This study proposed two different types of weight- ing, both from Frequenstist and Bayesian point of view. Besides the LC and its variants using different estimation procedure; application of Bayesian framework to fit LC method is also getting popular in recent era. Several approaches in Bayesian framework exist in literature, for example, Czado et al (2005) presented a Bayesian method for estimating the Poisson log- bilinear formulation of the Brouhns et al (2002) LC model. Few studies proposed extension of LC method in Bayesian framework as well (Wi´ sniowski et al, 2015). Comparison between meth-

• ds from Frequentist and Bayesian point of view is not available till now, as several different

approaches were done in Bayesian context and all of the studies compared the results with basic LC model. The basic LC model still remains the gold-standard for judging the performance

• f newly invented forecasting model in mortality and life expectancy. As mentioned before in

Introduction, most of these comparison holds mainly for developed countries with longer time periods and better quality (both in terms of data quality and mortality pattern over the ages and years). In this study we filled the gap for comparison between existing extrapolative meth-

• ds for high mortality countries of Eastern Europe and contrasted the results with probabilistic

forecast (current UN method).

2.8 Accuracy of the fitted models and fitting period

In current study we used the whole available time period for fitting where it was possible. To estimate the forecast accuracy, we used the Mean Forecast Error (MFE) and Mean Squared Forecast Error (MSFE). This MFE is the average of errors (difference between actual and forecast) which is a measure of bias as well. Symbolically, MFE = 1 (p + 1)q

q

• j=1

p

• x=0
• yx,j −

yx,j|j−h

• ,

j = (n − 10 + h), . . . , n MSFE = 1 (p + 1)q

q

• j=1

p

• x=0
• yx,j −

yx,j|j−h 2 A low value of the MFE may conceal forecasting inaccuracy due to the offsetting effect of large positive and negative forecast errors. It’s also possible to compare the accuracy in term

• f percentage, however, we are omitting that estimation for the current study.

3 Results

The data we used for different countries and their last observed life expectancy at birth are given below (Table-1). For most of these countries the life table started from 1950 although we 8

could not utilize the whole available data for several of these countries due to lower data quality mentioned by human mortality database (Human Mortality Database, 2016). Only Bulgaria has available data from 1947. Still, to fulfill all the assumptions we started the fitting from

• 1950. So, we fitted the models for 18 series of life tables for male and female separately, for

illustration we later presented the results for females only. The fitting of Bayesian Hierarchical Model for UN forecast is not life table based like the extrapolative methods, so we discussed that in later part of results. Table 1: Fitting periods for the countries and life expectancy at birth (HMD-2016). Country Starting year End Year e0 (Female) e0 (Male) Belarus 1970 2014 78.43 67.81 Bulgaria 1950 2010 77.25 70.31 Estonia 1959 2013 81.33 72.72 Hungary 1960 2014 79.24 72.26 Latvia 1970 2013 78.73 69.26 Lithuania 1959 2013 79.37 68.52 Russia 1970 2014 76.48 65.26 Slovakia 1962 2014 80.32 73.25 Ukraine 1970 2013 76.21 66.31

3.1 Comparison of forecast: extrapolative methods

To compare the forecast of the different models we tried different future forecasting periods. The following table showed the forecast of life expectancy at birth for 20 years ahead (Table-2), which later followed by the forecast for 30 years ahead (Table-5) and 40 years ahead (Table-6). Here LC stands for basic Lee-Carter model (1992); LCPoisson stands for Lee-Carter model with Poisson regression (Brouhns et al, 2002); LM stands for Lee-Miler model (2001); Booth, Maindonald and Smith model (2002); HU stands for basic Hyndman-Ullah model (2007); HUR stands for Robust Hyndman-Ullah model (2007); HUW stands for weighted Hyndman-Ullah model (2007) and UN stands for UN life expectancy forecast using Bayesian framework (Raftery et al, 2013). Forecast and biases are presented only for the models for which the forecast is higher than last

• bserved life expectancy. These notations are used in all over the paper.

Table 2: Comparison of life expectancy forecasting for selected countries with high mortality (h = 20). Country e0 LC LCPoisson LM BMS HU HUR HUW UN Belarus 78.43

• 78.433
• 79.679

78.672 Bulgaria 77.25 81.183 78.698 79.030 78.654 77.512 79.124 80.192 78.901 Estonia 81.33 83.570 83.588 83.804 86.526 82.050 81.316 82.760 83.381 Hungary 79.24 81.173 81.055 81.578 82.901 82.222 80.847 83.860 81.165 Latvia 78.73 80.416 80.284 80.421 83.025 79.943 80.270 80.556 80.361 Lithuania 79.37 81.294 80.070 80.824 82.635 79.409 79.340 81.033 80.938 Russia 76.48

• 77.106
• 77.463

Slovakia 80.32 82.874 82.660 82.640 83.621 82.113 81.944 82.239 82.109 Ukraine 76.21

• 78.258

77.448 77.507

Results are showed for female only. e0 is the last observed life expectancy at birth of the fitting period from human mortality database.

9

Although we included the UN forecasting method to compare with our findings, still we limit the main topic of discussion only on the variants of LC methods. To analyze the performance

• f the fitted models (except the UN forecast), we further analyzed the Mean Forecast Error

(MFE) and Mean Squared Forecast Error (MFSE) of the fitted models as measure of forecast

• accuracy. The results are presented below in table 3 and 4 respectively.

Table 3: One step ahead forecast bias of the fitted models in terms of MFE. Country LC LCPoisson LM BMS HU HUR HUW Belarus

• 0.194
• 0.935

Bulgaria

• 0.010
• 0.070
• 0.076
• 0.072

0.001

• 0.394

1.181 Estonia

• 0.064
• 0.062
• 0.040
• 0.026

0.001 0.103 0.016 Hungary

• 0.023
• 0.012
• 0.007

0.001 0.001 0.346 0.590 Latvia

• 0.026
• 0.026
• 0.019
• 0.014

0.001 0.063 0.470 Lithuania

• 0.025
• 0.047
• 0.048
• 0.018

0.001 0.138 0.795 Russia

• 0.066
• Slovakia
• 0.057
• 0.054
• 0.054
• 0.016

0.001 0.064

• 0.018

Ukraine

• 0.029
• 0.014

Table 4: One step ahead forecast bias of the fitted models in terms of MSFE. Country LC LCPoisson LM BMS HU HUR HUW Belarus

• 0.151
• 0.616

Bulgaria 0.075 0.044 0.047 0.049 0.258 0.267 0.688 Estonia 0.057 0.057 0.057 0.027 0.476 0.255 0.289 Hungary 0.037 0.036 0.037 0.006 0.184 0.162 0.269 Latvia 0.019 0.019 0.020 0.012 0.629 0.437 0.776 Lithuania 0.042 0.021 0.030 0.010 0.428 0.276 0.482 Russia

• 0.091
• Slovakia

0.053 0.053 0.053 0.016 0.187 0.111 0.220 Ukraine

• 0.073

0.105 The mortality regime were notably higher for three countries; namely Belarus, Russia and

• Ukraine. All of the models gave poor fit (and forecast) for these countries. For Belarus, only

Weighted HU method worked well; all other methods extrapolate lower or equal life expectancy for Belarus. For Russia and Ukraine, Robust HU method was appropriate (for UKraine the Weighted HU method as well); all other methods provide lower or equal trend. This information was also supported by the measures of forecast biases (Table 4 & 5). We estimated one step ahead forecast bias by MFE and MSFE. For most of the cases, the models were under estimating the forecast. Although it performed well, still the highest forecast biased was observed for weighted HU method for Bulgaria. The basic LC method returned high bias for all of the

• countries. Some values of MSFE suggest over-fitting of the models for some countries. For

developed countries almost same situation was observed which settled the conclusion that no model perform uniquely well for all countries (Shang et al, 2011; Shang, 2012). The better fit

• f HU methods for higher mortality regime can be explained by the assumption of the models.

We used weighted penalized regression splines with a monotonicity constraint for all of the HU methods (Hyndman & Ullah, 2007; Hyndman, 2015). The robust HU method was defined to have better fit for country with outliers in mortality pattern. The mortality curves of Russia and Ukraine are given below for the fitting period (Figure-1). Accidental hump is clearly visible for Russia along with wide heterogeneous mortality pattern across the years for senescence

• mortality. For both of the countries, mortality rate declined rapidly for childhood.

10

Figure 1: Observed mortality rates of Russia (1970-2014) and Ukraine (1970-2013). However, better fit of weighted HU method for Belarus and Ukraine is subject to analyze. The weighted HU method was defined for countries with long time series data; which was not the case here. Due to caution notes from human mortality database we started the fitting period only from best quality mentioned in human mortality database. We also tried with all the available data from HMD for all these three countries which provided very low forecast for all of the methods (except the UN forecast). Thus, data quality remained a restriction for this comparison of the methods for high mortality countries. Estonia has comparatively better mortality scenario compare to the other countries. The life expectancy of Estonia was highest among all of these countries. This better mortality trend

• f Estonia was also reflected in the forecasting; forecasts of all the LC variants were closer for

Estonia, although HU methods seems to have slightly lower forecasts for Estonia than that

• f LC variants (LC, LCPoisson, LM, BMS). Almost same pattern was observed for Slovakia,

the second highest life expectancy was observed for Slovakia in the end of the fitting period. Besides slightly lower (but close to each other) forecasts from HU variants, the other cluster was observed for LC variants. Table 5 and 6 showed the forecast for 30 and 40 years ahead respectively for all of these methods (including UN forecast). 11

Table 5: Comparison of life expectancy forecast for selected countries with high mortality (h = 30). Country e0 LC LCPoisson LM BMS HU HUR HUW UN Belarus 78.43

• 78.433
• 80.251

79.370 Bulgaria 77.25 81.703 79.268 79.683 79.227 78.049 80.037 81.472 79.740 Estonia 81.33 84.634 84.676 84.944 88.570 83.133 81.995 83.866 84.528 Hungary 79.24 82.183 82.069 82.649 84.499 83.469 81.935 85.914 82.311 Latvia 78.73 81.167 81.043 81.191 84.776 80.659 81.070 81.468 81.384 Lithuania 79.37 81.835 80.585 81.413 83.993 79.906 79.712 81.909 81.851 Russia 76.48

• 77.477
• 78.265

Slovakia 80.32 83.922 83.680 83.682 84.997 83.047 82.813 83.235 83.273 Ukraine 76.21

• 79.247

78.053 78.371

Results are showed for female only. e0 is the last observed life expectancy at birth of the fitting period from human mortality database.

Table 6: Comparison of life expectancy forecast for selected countries with high mortality (h = 40). Country e0 LC LCPoisson LM BMS HU HUR HUW UN Belarus 78.43

• 78.433
• 80.789

79.952 Bulgaria 77.25 82.160 79.769 80.237 79.756 78.496 80.834 82.611 80.633 Estonia 81.33 85.638 85.702 86.012 90.385 84.161 82.637 84.900 85.703 Hungary 79.24 83.159 83.049 83.680 85.980 84.642 82.995 87.813 83.539 Latvia 78.73 81.881 81.763 81.919 86.369 81.341 81.830 82.333 82.432 Lithuania 79.37 82.309 81.047 81.931 85.202 80.358 80.044 82.726 82.763 Russia 76.48

• 77.818
• 79.079

Slovakia 80.32 84.911 84.645 84.661 86.235 83.934 83.638 84.167 84.436 Ukraine 76.21

• 80.190

78.637 79.243

Results are showed for female only. e0 is the last observed life expectancy at birth of the fitting period from human mortality database.

We expanded the extrapolation to 30 and 40 years ahead to check whether any convergence

• ccur in the estimates of the different methods. For the higher mortality regime, the forecast

still showed underestimation for all of the methods. Like forecast for 20 years ahead; only weighted HU method performed well for Belarus and Ukraine while robust HU method worked well for Russia and Ukraine. All of these three countries showed declining trend for 30 and 40 years ahead forecast of life expectancy for others method; which contradicts with any theory

• f demographic transition or linear trend of rise in level of life expectancy (Oeppen & Vaupel,

2002). This findings are important for future application of these extrapolation based methods for any higher mortality country. For the other countries of consideration, we did not find any closer results from different

• methods. Slower rise in longevity is indicated by all the methods as pace of rise in life expectancy

at birth is lower for 30 or 40 years ahead forecast compare to that of forecast of 20 years

• ahead. To summarize the performance of the models we may conclude that the basic LC model

underestimated the forecast compare to the later variants, whereas LM and LC model under Poisson regression performed better for other countries. The HU variant performs well almost for all of these countries although it has some possible overestimation scenario also; for example the forecast of weighted HU method for Hungary. The BMS model provided overestimated 12

forecast for all the countries except Belarus, Russia and Ukraine.

3.2 UN forecasting: Application of Bayesian Hirarchical Model

This forecast method is probabilistic rather than old deterministic methods used by UN. In the methodology it is already explained about how the life expectancy models are considered in Bayesian hierarchical structure (Reftery et al, 2013). The analysis of this method is mainly defined for data of World Population Prospects (WPP); the UN also adopted this method for life expectancy projection by switching from the old deterministic approaches. One problem of this forecast is, all the input data are on 5 (calender) years interval in WPPs (or from other sources), and the forecast we have are also for (the mid year of) 5 calender years. Thus the values we showed in Table 2 for the year 2034 is actually the forecast for time time period 2030-2035. Similarly for Table 5, the year 2044 stands for forecast for the year 2040-45 and in Table 6 the year 2054 stands for forecast for the year 2050-55. For the 9 countries considered in this study, we utilized the life expectancy at birth in 5 year basis from Human Mortality Database instead of using the UN data. We project the life expectancy at birth up to year 2100, however, for comparison with existing extrapolative methods up to 40 years; we presented the results only for 20, 30 and 40 years ahead forecast. The forecast was implemented by ‘bayesLife’ package in R. Besides the data from Human Mortality Database, we also performed the analysis using the data from WPP-2012 (imple- mented through ‘wpp2012’ package); which is later presented in Appendix section for interested

• readers. The simulation was done with 160,000 iterations (10,000 burn-in); the thinning interval

was 10 with number of chains were 3. Instead of country specific results, the posterior distri- butions and trace plots of the hyper-parameters for all the countries are given in the appendix section for interested readers. The whole simulation were conducted for all the 162 countries mentioned in WPP 2012 along with computation of hyper-parameters of all these countries. However, we kept our results limited for the desired countries of comparison only. It was already mentioned in previous section that irregular trend in mortality was visible in case of several of these countries. Thus, the trend of life expectancies in 5 year basis was not also free from it; fluctuation in trend of life expectancies remained for Belarus, Russia and Ukraine during the fitting period. Like HMD, these similar pattern was also observed in case of data from WPP-2012. Nevertheless, UN forecasts showed increasing trend in life expectancy at birth for all of these countries (both for HMD-2016 and WPP-2012). Unlike LC and HU variants, UN forecast gave better forecast for all three high mortality countries- Belarus, Russia and Ukraine. For 20 years ahead forecast, the UN forecast was lower than weighted HU method for Belarus, but still higher than last observed life expectancy at birth in 2014. It still remained lower for Belarus in 30 years ahead and 40 years ahead forecast, although the forecast became higher than weighted HU method in case of Russia and Ukraine. For Russia, the forecast produced by the UN forecast was the highest among all forecast methods. For all other countries the forecast produced by UN forecast method lied between the values

• f the other forecast methods. Some exceptions are also observed where probabilistic forecast

was lower than other methods in case of the WPP-2012 data. For example in case of Hungary and Latvia, the UN forecast were lowest for all the cases (20, 30 or 40 years ahead). Forecast for Slovakia was lowest in case of 30 and 40 years ahead (Table 9 in Appendix A2.1). One serious shortcoming of this UN forecast is, this technique is not based on life table like the previous LC or HU variants. This shortcoming limits the application of this Bayesian Hi- erarchical Model for forecasting mortality based on age-specific mortality rates. Again, as this technique was defined keeping the usual UN forecasting technique as basis, so, it can’t produce 13

precise point estimate of a single year. The mean and standard deviation of the estimated parameters (for HMD-2016 data) are attached in appendix for interested readers. To analyze the uncertainty in the context of this deterministic approach, we plotted the 95% prediction interval of forecast of life expectancy at birth for all the 9 countries below (Figure-2-5). These prediction intervals are showed here over the whole projection period (up to year 2100). It is easily observable that, higher mortality countries has more wide prediction interval than that of comparatively lower mortality countries. The diverging mortality pattern affects the forecast of life expectancy as well. Figure 2: Life expectancy at birth forecast for Belarus and Bulgaria by UN forecast. Figure 3: Life expectancy at birth forecast for Estonia and Hungary by UN forecast. 14

Figure 4: Life expectancy at birth forecast for Latvia and Lithuania by UN forecast. Figure 5: Life expectancy at birth forecast for Russia, Slovakia and Ukraine by UN forecast.

3.3 Country specific illustration: Hungary

To illustrate the performance of each of the fitted models, Hungary is considered in this section. Hungary is chosen in this study because Hungary has high mortality regime similar to the Eastern Europe, Eastern-Central Europe and Baltic states like Latvia and Lithuania. The common features of the mortality scenario of these countries can be characterized by presence

• f high level of mortality form Cardiovascular diseases, external causes of deaths (Monostori et

al, 2015). The trend of life expectancy at birth of Hungary is plotted below in figure-6 (Human Mortality Database-2016). 15

Figure 6: Life expectancy at birth of Hungary (1960-2014). In beginning of 1990s, life expectancy of the Hungarian population were among the lowest in Europe. Male life expectancy at birth in Hungary was the third lowest (after Latvia and Estonia) among the countries that currently belong to the European Union. The long trend

• f high mortality pattern among the Hungarian males did not end with the start of the socio-

economic transition, same pattern is also observed for Hungarian female (Monostori et al, 2015). As a consequence, the improvement in morality and longevity that has been observed in Hungary is not unique. Similar progress has taken place in Western European countries where the improvement got underway in the 1960s and has continued uninterrupted since then (Mackenbach, 2013; Monostori et al, 2015). As mentioned before, Hungary achieved all the gain in longevity from improve in mortality scenario from some specific causes of deaths only. These reasons are not highly attributable to infant mortality rather it came from improvement in middle-aged mortality pattern (due to several age-specific diseases and/or external causes of deaths). The cause-specific decomposition of differences in life expectancy in previous study makes it possible to identify the contribution of age groups by cause of death (Monostori et al, 2015). Decomposition of life expectancy at birth showed that the seven-year gain in male life expectancy between 1990 and 2013 are mainly attributable to decline in cardiovascular mortality which corresponds to 40% of the total gain in longevity (Monostori et al, 2015). Although steady declined in observed, still mortality due to alcohol consumption is high in Hungary (Monostori et al, 2015). The log mortality rates for Hungarian male and female from human mortality database are illustrated below (Figure 7 & 8). Irregular pattern in mortality is visible in young and senescence periods. Gender gap is also visible in different age groups during this period. 16

Figure 7: Log mortality rates of Hungarian males in 1960-2014 (Human Mortality Database- 2016). Figure 8: Log mortality rates of Hungarian females in 1960-2014 (Human Mortality Database- 2016). 3.3.1 Forecast of life expectancies We already showed the forecasts for life expectancy at birth for Hungarian female for 20, 30 and 40 years ahead. The probabilistic approach for forecasts lied among all forecast methods applied for Hungary. In this section we compared the forecasts for all the LC and HU variants

• nly for males and females of Hungary. Life expectancy at birth was 72.26 and 79.24 years

respectively for males and females in 2014. Same fitting periods are used as before for all of the models. Except the basic LC model for males, the other models perform well than that of Belarus, Russia and Ukraine. The BMS and weighted HU method provided higher forecast of life expectancy at birth for female in long run. For both male and female life tables, all the HU variants provided higher forecasts than that of LC variants. The results are given below (Table-7). 17

Table 7: Forecast of life expectancies at birth for Hungarian males and females. Forecast period LC LCPoisson LM BMS HU HUR HUW Male, 2034

• 73.927
• 77.351

77.529 77.956 Female, 2034 81.173 81.055 81.578 82.901 82.222 80.847 83.860 Male, 2044

• 74.422
• 79.328

79.695 80.416 Female, 2044 82.184 82.070 82.650 84.499 83.469 81.936 85.914 Male, 2054

• 74.821
• 81.033

81.629 82.971 Female, 2054 83.159 83.049 83.680 85.980 84.642 82.995 87.813 To see the performance of the life expectancy forecasting more closely, we extrapolated the forecast for remaining life expectancy at age 60. The forecasts for 20, 30, and 40 years ahead for both sexes are given below (Table-8). In Hungary, the remaining life expectancy at age 60 was 17.59 years and 22.30 years respectively for males and female (Human mortality database- 2016). Remaining life expectancy in the senescence period is a common indicator of longevity. Although the performance of basic LC model was not good for male life tables of Hungary, it gives reasonable forecast for age 60. However, the LC model under Poisson regression and BMS model produced underestimated forecast for remaining life expectancy of males at age 60. In

• ther cases, all the model performed well, although higher forecast is observed for weighted HU

model for females compare to the other methods. Table 8: Forecast of reamining life expectancies at age 60 years for Hungarian males and females. Forecast period LC LCPoisson LM BMS HU HUR HUW Male, 2034 22.533 17.396 19.781 17.396 20.317 20.984 20.745 Female, 2034 24.266 24.163 24.231 25.055 24.687 23.471 25.808 Male, 2044 23.993 17.800 21.026 17.800 21.590 22.396 22.538 Female, 2044 25.160 25.059 25.180 26.324 25.674 24.254 27.517 Male, 2054 25.770 18.217 22.421 18.217 22.799 23.745 24.638 Female, 2054 26.043 25.942 26.114 27.546 26.627 25.057 29.157 Moderate improvements in the mortality of people aged 60 years or over were observed in previous studies, which is not like the pattern observed in case of middle aged population (Monostori et al, 2015). More rapid increase in life expectancies were observed for middle-aged and elderly male population which indicates that there are still scope of further improvement in adult mortality. The other determinants of longevity also have vital role in case of Hungary. Previous study revealed that high disparity is available in case of life expectancy in terms of education and marital status. Life expectancy at birth for females with higher education was 5.8 years higher than that of women with primary education at 2012. For males of Hungary, the difference was 12.5 years. Married men and women have longer life expectancies and lower mortality rates from all major causes of death (Monostori et al, 2015). 3.3.2 Forecast of mortality rates The contrast of life expectancy forecast at two different ages emerges the topic of fitting of the mortality rates. As discussed before, all of these methods (except UN forecast) fit the log mortality rates in the observed time periods and hence forecast for the future mortality. In this section we will compare and contrast the forecast in log mortality rates among 7 LC variants (including HU methods). The fitted parameters with forecast of parameter kt and the observed and fitted log mortality rates for females are presented below (Figure 9 & 10). 18

Figure 9: Fitted components of basic LC method for females of Hungary. Figure 10: Observed and 20 years ahead forecast of log mortality rates for females of Hungary by basic LC method. As discussed before, the parameter kt is used for forecast; the yellow interval may be in- terpreted for interval forecasting as well. For the observed log mortality rates, the mortality scenario were improving for females of Hungary all over the life span. As mentioned before, the main source of gain in female life expectancy was result of decline in adult and senescence 19

• mortality. Previous study revealed the contribution of the female population aged 64 years and
• ver was more substantial in rise of life expectancy than that of the middle-aged who were also

affected by the economic crisis in past (Monostori et al, 2015). Previous studies explained the trend of female mortality pattern for Hungary (Monostori et al, 2015). The most remarkable change in recent mortality pattern compared to 1990s is the decline in the number of cardiovascular deaths; at least 10 percent deaths reduced for female. Prevalence of lung cancer related mortality rose sharply among females of Hungary in recent

• past. Notable change is observed for the deaths caused by diseases of the digestive system and

alcohol-related liver diseases. However, the most significant mortality decline were observed for different external causes of deaths and suicidal mortalities. For both of the cases, the number

• f deaths were half in the end of 2013 compare to that of 1990 (Monostori et al, 2015).

The other mortality models also performed well for Hungarian female life tables compare to that of male life tables. We already compared the forecasts for life expectancies at different ages, figure 11 and 12 illustrated the forecast of log mortality rates for Hungarian female by different models. The fitted components (with forecast) of all the models are given in appendix for interested readers. We illustrated the mortality forecast for 40 years ahead for all of the models to check the results converged among the models in long run; if any (Figure 11 & 12). The lighter colored lines (starting from yellow) represent the earlier period of the forecast while the darker (mainly blue) lines stands for later period of the forecast. All of the models forecast rapid and sharp decline in infant and adult mortality. Divergence in log mortality is observed for later age groups in all of the forecast models. All the LC variants showed irregular trend in decline in mortality for young ages. The case is different for HU methods due to sophisticated smoothing techniques utilized in those variants. All of the HU variants forecast regular patter over all the age groups till early senescence period. However, it may possible that these variants get

• ver-fitted for middle-aged mortality. Besides HU variants, the BMS model gave almost same

forecast for log mortality rate for adult to middle age female of Hungary. The forecast is almost similar for senescence mortality. All the models showed slower decline in later life mortality. For HU variants, these decline is more visible than LC, LM or BMS model. There are something common in all of the LC variants for age around 60s which showed future stalling in mortality level at those ages. Due to smoothing, its not so visible for HU methods, but it is common for all the old LC variants. 20

Figure 11: Comparison of log mortality rates forecast for females of Hungary: LC variants. 21

Figure 12: Comparison of log mortality rates forecast for females of Hungary: HU variants. 22

4 Conclusion

We compared several forecasting methods in line of research already done by Shang et al (2011)

• r Shang (2012) but different from that as we utilized the forecasting models for high mor-

tality countries from Eastern Europe and Russia. These countries have the characteristics of comparatively higher mortality regime; data with shorter/limited time frame for fitting period

• f the models along with quality issues. Moreover, we compared our forecasts with that of a

Bayesian Hierarchical Model used by UN. The findings from LC variants suggest that forecasts

• f countries with relatively high mortality suffer from under-fitting problems, more than what

has been found for low mortality countries (Shang et al, 2011; Shang, 2012). Some of the variants were not appropriate for too high mortality countries. Over fitting of the models were also seen for some cases. Robust and Weighted HU methods performs well for the countries with highest mortality pattern among the 9 countries, while other methods failed to do so. For comparatively better mortality regime, all the models produced fairly reasonable forecast. The probabilistic approach of forecast (UN forecast method) also produced wider prediction interval for some of the countries. So to conclude, like previous studies (Shang et al, 2011; Shang 2012); we are also unable to declare a specific model best for all of the (comparatively) high mortality

• countries. From view of limited time series data and absence of series life tables, UN forecast

(implemented through Bayesian Hierarchical Model) can be an appropriate alternative. On the

• ther hand, in presence of high mortality along with irregular mortality pattern over the fitting

period, the robust HU method can be used for best result. However, in case of long time series data with irregular trend, weighted HU method may perform well compare to other LC variants. As regarding the aim of this study (comparison of forecasting methods for limited time frame, high mortality and quality issue); the forecasts seriously suffered for irregularity in mor- tality trend and data quality. Although previous studies mentioned different mortality pattern all over the Europe, forecast for countries with similar socio-economic structure and mortality regime supposed to provide converging mortality scenario in future (Vallin & Mesl´ e, 2001). That was not the case here. Human mortality database mentioned lower data quality for some

• f the years for few of these countries due to data source; we tried to fit the models with both

including and excluding those years. The forecast was really questionable considering those years for fitting period. However, we could not fitted the models for Estonia and Lithuania by excluding the problematic years as those are almost in middle of the fitting period for Estonia and make the fitting period too short for forecast in case of Lithuania. It was mentioned in earlier studies that long time series is preferable for fitting of these models, that condition could not be hold for many of these countries due to data problem (Booth et al, 2002). Availability of reliable data made us unable to analyze the performances of these models to other developing countries as well; that might be another future extension of this work. We may suggest some further extension of this line of research on comparison of forecasting

• methods. The basic LC models can be estimated using different algorithms like weighted least

square (WLS), singular value decomposition technique (SVD), maximum likelihood estimation (MLE). We estimated the model only from SVD in this study; the others estimation procedures may provide different results than that we had. For the HU variants; we only used the weighted penalized regression splines to smooth ft(x); use of others splines methods or changing the basic smoothing techniques may provide more better insight to understand and forecast the mortality

• regime. Another shortcoming of current study is, we did not compare the interval forecast for

these methods in current study, which could be another topic of further research. An important point regarding model fitting was that we considered all the countries in setting provided by human mortality database, considering certain censoring in ages like 0 to 100+ or 0 to 90+ may provide different results than that we had. 23

Bayesian framework is getting popular in recent era for forecast of mortality and life ex- pectancy (Raftery et al, 2013). We used the data of Human Mortality Database and World Population Prospect-2012; application of World Population Prospect 2015 (WPP-2015) may produce different results than ours. However, a serious drawback of this Bayesian hierarchical model is its computational complexity; it took almost 3 days to converge the mcmc in current

• study. Also, it is not based on whole life table, rather it works only on time series data of

life expectancy at birth. Several approaches exist on extending basic LC method in Bayesian

• framework. It is possible to have better forecast from understanding the underlying probabil-

ity distribution properly than the extrapolative approaches. However, more in-depth study is required to compare the forecast among Frequentist and Bayesian approaches. Nevertheless, invention of new forecasting method is advisable from findings of current study as none of the models were undoubtedly suitable for high mortality countries.

Acknowledgement

The analysis performed in this study can be implemented by ‘Demography’ package of R for all

• f the LC and HU variants. For Bayesian forecast we used the ‘bayesLife’ package in R.

References

[1] Booth, H.; Maindonald, J.; Smith, L. Applying Lee-Carter under conditions of variable mortality decline. Population Studies, 56(3), 325-336, 2002. [2] Brouhns, N.; Denuit, M.; Vermunt, J. A Poisson log-bilinear regression approach to the construction of projected lifetables. Insurance: Mathematics and Economics, 31, 373-393, 2002. [3] Girosi, F.; King, G. Demographic Forecasting. Cambridge University Press, Cambridge, 2006. [4] Human Mortality Database. University of California, Berkeley (USA), and Max Planck Institute for Demographic Research (Germany), 2016. [5] Hyndman, R.J. Demography: Forecasting mortality, fertility, migration and population

• data. R package version 1.18, 2015.

[6] Hyndman, R.J.; Ullah, M.S. Robust forecasting of mortality and fertility rates: A functional data approach. Computational Statistics and Data Analysis, 51(10), 4942-4956, 2007. [7] Hyndman, R.J.; Shang, H.L. Forecasting functional time series (with discussion). Journal

• f the Korean Statistical Society, 38(3), 199-221, 2009.

[8] Lee, R.D.; Carter, L.R. Modeling and forecasting U.S. mortality. Journal of the American Statistical Association, 87(419), 659-671, 1992. [9] Lee, R.D.; Miller, T. Evaluating the performance of the Lee-Carter method for forecasting

• mortality. Demography, 38(4), 537-549, 2001.

[10] Keyfitz, N. Experiments in the projection of mortality. Canadian Studies in Population, 18(2), 1-17, 1991. [11] Mackenbach, J. P. Convergence and divergence of life expectancy in Europe: a centennial

• view. European Journal of Epidemiology 28, 229-240, 2013.

24

[12] Monostori, J.; ¨ Ori, P.; Sp´ eder, Z. Demographic Portrait of Hungary 2015. Hungarian De- mographic Research Institute, Budapest, 2015. [13] Oeppen, J.; Vaupel, J.W. Broken limits to life expectancy. Science, 296: 1029-1031, 2002. [14] Pollard, J.H. (1987). Projection of age-specific mortality rates. Population Bulletin of the United Nations, 21-22: 55-69, 1987. [15] Raftery, A.E.; Chunn, J.L.; Garland, P.; ˇ Sevˇ c´ ıkov´ a, H. Bayesian Probabilistic Projections

• f Life Expectancy for All Countries. Demography, 50, 777-801, 2013.

[16] Renshaw, A.E.; Haberman, S. Lee-Carter mortality forecasting with age-specific enhance-

• ment. Insurance: Mathematics and Economics, 33: 255-272, 2003.

[17] ˇ Sevˇ c´ ıkov´ a, H., Raftery, A.E. bayesLife: Bayesian Projection of Life Expectancy. R package version 3.0-0, 2015. [18] Shang, H.L. Point and interval forecasts of age-specific life expectancies: A model averaging

• approach. Demographic Reseach, 27(21), 593-644, 2012.

[19] Shang, H.L.; Booth, H.; Hyndman, R.J. Point and interval forecasts of mortality rates and life expectancy: A comparison of ten principal component methods. Demographic Research 25, 173-214. [20] Stoeldraijer, L.; Duin, C.V.; Wissen, L.; Janssen, F. Impact of different mortality forecast- ing methods and explicit assumptions on projected future life expectancy: The case of the

• Netherlands. Demographic Research, 29(13): 323-354, 2013.

[21] United Nations. World Population Prospects: The 2008 Revision. New York: United Na- tions, Department of Economic and Social Affairs (DESA), Population Division, Population Estimates and Projections Section, 2009. [22] United Nations. World Population Prospects: The 2012 Revision. New York: United Na- tions, Department of Economic and Social Affairs (DESA), Population Division, Population Estimates and Projections Section, 2013. [23] Vallin, J.; Mesl´ e, F. Trends in mortality in Europe since 1950: age-, sex- and cause-specific

• mortality. In Trends in mortality and differential mortality, ed. J. Vallin, F. Mesl´

e and T.

• Valkonen. Strasbourg: Council of Europe Publishing, 31-184, 2001.

[24] Wi´ sniowski, A., Smith, P., Bijak, J., Raymer, J., Forster, J. Bayesian population forecast- ing: extending the Lee-Carter method. Demography, 29(13): 323-354, 2013. 25

Appendix

• A1. Mortality forecast for Hungarian female: components of various fitted

models

Figure 13: Fitted components of LC model under Poisson regression for females of Hungary. Figure 14: Fitted components of LM model for females of Hungary. 26

Figure 15: Fitted components of BMS model for females of Hungary. 27

Figure 16: Fitted components of HU model for females of Hungary. 28

Figure 17: Fitted components of Robust HU model for females of Hungary. 29

Figure 18: Fitted components of Weighted HU model for females of Hungary. 30

• A2. Findings from Bayesian Hierarchical Model (UN forecast)

A2.1 Comparison between forecasts for different data source (HMD-2016 & WPP- 2012) Table 9: Forecast of life expectancies at birth using UN forecast method for HMD-2016 and WPP-2012†. Data & Forecast year Belarus Bulgaria Estonia Hungary Latvia Lithuania Russia Slovakia Ukraine HMD, 2034 78.672 78.901 83.381 81.165 80.361 80.938 77.463 82.109 77.507 WPP, 2034 77.225 78.914 81.395 80.600 78.960 79.807 76.223 81.071 75.621 HMD, 2044 79.370 79.740 84.528 82.311 81.384 81.851 78.265 83.273 78.371 WPP, 2044 78.022 79.858 82.371 81.592 79.724 80.645 77.204 82.057 76.219 HMD, 2054 79.952 80.633 85.703 83.539 82.432 82.763 79.079 84.436 79.243 WPP, 2054 78.812 80.729 83.315 82.707 80.564 81.428 77.991 83.053 76.906

†For WPP-2012 data, the life expectancies were available from 1873 to 2015. The results are showed for females only.

A2.2 Fitted parameters Table 10: Fitted parameters from Bayesian Hirerchical Model. Parameter Mean SD ∆1 12.891 3.765 ∆2 40.655 3.985 ∆3 4.003 10.116 ∆4 17.911 5.503 k 3.969 0.252 z 0.634 0.018 λ1 0.014 0.005 λ2 0.011 0.003 λ3 0.008 0.004 λ4 0.007 0.003 λk 0.684 0.132 λz 32.860 9.344 ω 1.558 0.036

Data: Human Mortality Database-2016.

A2.3 Posterior distributions and traceplots The posterior distributions of each of the parameters and trace plots for the parameters are given below in figure-19 & 20 respectively. It should be noted that The R package ‘bayesLife’ denote the symbols ∆i by ‘Trianglei’ and λi by ‘lambdai’. 31

32 Figure 19: Posterior ditributions of the parameters of Bayesian Hierarchical Model.

33 Figure 20: Trace plots of hyper parameters of mcmc.