A Study on Link Functions for Modelling and Forecasting Old-Age Survival Probabilities of Australia and New Zealand

Abstract

Forecasting survival probabilities and life expectancies is an important exercise for actuaries, demographers, and social planners. In this paper, we examine extensively a number of link functions on survival probabilities and model the evolution of period survival curves of lives aged 60 over time for the elderly populations in Australasia. The link functions under examination include the newly proposed gevit and gevmin, which are compared against the traditional ones like probit, complementary log-log, and logit. We project the model parameters and so the survival probabilities into the future, from which life expectancies at old ages can be forecasted. We find that some of these models on survival probabilities, particularly those based on the new links, can provide superior fitting results and forecasting performances when compared to the more conventional approach of modelling mortality rates. Furthermore, we demonstrate how these survival probability models can be extended to incorporate extra explanatory variables such as macroeconomic factors in order to further improve the forecasting performance.

1. Introduction

The average human lifespan has been increasing consistently throughout the developed nations in the last hundred years or so. Oeppen and Vaupel (2002) reported that the highest female period life expectancy at birth around the world each year has been growing by about 0.24 years annually for more than a century. In Australia, period life expectancy at age 60 has grown from 16.5 years in 1950 to 25.0 years in 2017; in New Zealand, it has increased from 16.9 years in 1950 to 24.0 years in 2013. Life expectancy is one of the major indicators of a country’s wellbeing. Forecasting life expectancies accurately is of critical importance for government social planners as well as demographic researchers and industry practitioners (e.g., Lee et al. 1995; Nikolaevich 2019).

The continuous mortality decline, together with the lack of clear improvement in the maximum lifespan, has caused the common phenomenon of a “rectangularisation” of the survival curve. This kind of observations has been well documented in earlier studies (e.g., Fries 1980; Cheung et al. 2005). As shown in Figure 1, it can be characterised by an upward and rightward shift of the period survival curve across successive time periods. The increasing area under the survival curve over time refers to the extent of rectangularisation. This concept is one of the major analytical frameworks in demographic research and can actually be further adapted to describe past changes in survival levels and predict future survival rates. The recent developments in mortality projection methods can provide a useful reference for this direction of exploiting the trends in survival probabilities.

While there is a vast literature on modelling and projecting mortality rates (e.g., Lee and Carter 1992; Booth and Tickle 2008), relatively less attention has been paid on forecasting survival probabilities directly. Amongst the few previous works, De Jong and Marshall (2007) applied the probit link function to the survival probabilities of Australian females and assumed that it is driven by a single time trend. Hatzopoulos and Haberman (2015) used the complementary log-log link function and age-cohort effects within the GLM (generalised linear model) framework for the cohort survival probabilities of a few European countries. Wong and Tsui (2015) proposed a new survival function for US women and men and modelled the changes of its parameters over time by autoregressive processes. Tan et al. (2016) constructed a hybrid survival curve and applied the logit link function to the annualised survival probabilities with two or three sets of time-varying parameters using Swedish and Bulgarian data.

There are potential advantages of modelling and projecting the survival probabilities rather than the mortality rates. First, when producing life expectancy forecasts, it would be more convenient to work with the survival probabilities directly without having to compound the mortality rates to form the survival rates. In a similar vein, when pricing pensions or annuities, the future probabilities of survival are a major input for the calculation process, and so it would be more natural to build a projection model that focuses on the survival probabilities. Moreover, as shown in the following sections, the survival probability patterns and trends are generally quite stable, which make their forecasting more straightforward than otherwise. We will show in our empirical analysis that some of these survival rate projection models can actually outperform the more conventional approach of a mortality rate projection model.

As mentioned above, forecasting survival probabilities are much less explored in the literature when compared to forecasting mortality rates. In this paper, we attempt to reduce this knowledge gap by examining extensively a number of link functions on survival probabilities and modelling the evolution of the parameters and so the survival rates over time. The link functions considered include the probit, complementary log-log, logit, gevit functions, and a new link function based on the theory of minima, and both the age and period effects are incorporated. Furthermore, we illustrate how the survival rate projection models can be extended to include additional explanatory variables such as the GDP (gross domestic product) per capita. Some previous examples of incorporating macroeconomic factors into mortality rate projection models are Hanewald (2011); Niu and Melenberg (2014); French and O’Hare (2014); and Seklecka et al. (2019). To the best of our knowledge, our paper provides the first attempt to incorporate macroeconomic factors into the survival rate projection models.

The remainder of the papers is as follows. Section 2 gives an introduction of the various survival rate projection models being considered. Section 3 compares the fitting performances of the models for Australian and New Zealand data. Section 4 sets forth an out-of-sample analysis for measuring the forecasting performances of the models. Section 5 studies the potential link between the survival probabilities and the economic growth. Section 6 concludes.

2. Survival Rate Projection Models

Suppose $q_{x,t}$ is the mortality rate at age $x$ in year $t$, and ${}_{n}p{}_{x,t}={\displaystyle {\prod }_{i=1}^n(1-q_{x+i-1,t})}$ is the corresponding period survival probability in year $t$ for a surviving period of at least $n$ years. The first link function we consider is the probit function used by De Jong and Marshall (2007). We set the model structure as ${\mathsf{\Phi }}^{-1}({}_{n}p{}_{x_0,t})=h(x,t)$, in which ${\mathsf{\Phi }}^{-1}$ is the inverse standard normal cdf (cumulative distribution function) and $h(x,t)$ is a regression structure allowing for the age and period effects with $x=x_0+n$. If one treats $1-{}_{n}p{}_{x_0,t}$ as the cdf of the future lifetime (within $n$ years) of a life aged $x_0$ in year $t$, ${}_{n}p{}_{x_0,t}$ as the survival function of that future lifetime, and $h(x,t_2)-h(x,t_1)=\lambda$ as a certain constant for $t_2>t_1$, it can be deduced that ${}_{n}p{}_{x_0,t_2}=\mathsf{\Phi }({\mathsf{\Phi }}^{-1}({}_{n}p{}_{x_0,t_2}))$$=\mathsf{\Phi }(h(x,t_2))=\mathsf{\Phi }(h(x,t_1)+\lambda )$$=\mathsf{\Phi }({\mathsf{\Phi }}^{-1}({}_{n}p{}_{x_0,t_1})+\lambda )$ and so $1-{}_{n}p{}_{x_0,t_2}$$=\mathsf{\Phi }({\mathsf{\Phi }}^{-1}(1-{}_{n}p{}_{x_0,t_1})-\lambda )$. It then means that under this probit model structure, the future death rates can be seen as a Wang transform (Wang 2000) of the past death rates, with the parameter $\lambda$ capturing the mortality decline. Note that the probit link function ensures that the survival rates are between 0 and 1, regardless of the values of $h(x,t)$, and that it is a symmetric link as $\mathsf{\Phi }(z)$ approaches 0 and 1 at the same pace.

The next one is the complementary log-log link function used by Hatzopoulos and Haberman (2015). We follow their model structure as $\ln (-\ln (1-(1-{}_{n}p{}_{x_0,t})))$$=\ln (-\ln ({}_{n}p{}_{x_0,t}))=h(x,t)$, that is, the link function is being applied on $1-{}_{n}p{}_{x_0,t}$, but not ${}_{n}p{}_{x_0,t}$. Unlike the probit function, the complementary log-log function is an asymmetric link, which would be useful for the rectangularisation patterns as in Figure 1. This asymmetry may suit survival modelling better when compared to a symmetric one. The survival rates under this link function are between 0 and 1.

We apply the logit link function from Tan et al. (2016) as $\ln ({}_{n}p{}_{x_0,t}/(1-{}_{n}p{}_{x_0,t}))=h(x,t)$. It is a symmetric link like the probit function, and it constrains the survival rates between 0 and 1. Its inverse function $\exp (z)/(1+\exp (z))$ is actually the cdf of the logistic distribution with the location parameter equal to 0 and the scale parameter equal to 1. Both the probit and logit functions are used extensively in binary response models.

Recently, Medford and Vaupel (2019) proposed the so-called gevit link function for modelling mortality rates. We apply this link function differently here as $[{(-\ln ({}_{n}p{}_{x_0,t}))}^{-\xi }-1]/\xi =h(x,t)$; that is, it is being applied on the survival probability but not on the death rate. This link is asymmetric and it constrains the survival rates like those above. Its inverse function $\exp (-{(1+\xi z)}^{-1/\xi })$ is indeed the cdf of the GEV (generalised extreme value) distribution for maxima with the location parameter equal to 0 and the scale parameter equal to 1. The extra shape parameter $\xi$ provides more flexibility to manage the extent of asymmetry.

Inspired by the gevit link function, which is based on maxima, we exploit the GEV distribution for “minima” instead (e.g., Liu and Li 2019) as an alternative. The cdf of the minima GEV distribution is specified as $1-\exp (-{(1-\xi z)}^{-1/\xi })$. Accordingly, we consider a new model structure $[1-{(-\ln (1-{}_{n}p{}_{x_0,t}))}^{-\xi }]/\xi =h(x,t)$. It is straightforward to deduce that this new link is also asymmetric and the resulting survival rates are within the valid range. We refer to it as the “gevmin” model structure in the following analysis.

If this new link is applied on $1-{}_{n}p{}_{x_0,t}$ instead of ${}_{n}p{}_{x_0,t}$, the model structure becomes $[1-{(-\ln (1-(1-{}_{n}p{}_{x,t})))}^{-\xi }]/\xi$$=[1-{(-\ln ({}_{n}p{}_{x_0,t}))}^{-\xi }]/\xi =h(x,t)$. Then if $\xi =0$, it reduces to $\ln (-\ln ({}_{n}p{}_{x,t}))=h(x,t)$, which means that the complementary log-log model structure above can effectively be seen as a specific example of this alternative model structure from minima.

Figure 2 compares the symmetry of the probit and logit functions with the asymmetry of the complementary log-log, gevit, and gevmin functions as described above. For the symmetric ones, the response approaches 0 and 1 at the same pace. For the asymmetric complementary log-log and gevit functions, the response approaches 1 slower than reaching 0. By contrast, under the new gevmin model structure, the response of the (inverse) link approaches 1 faster than reaching 0, which is an opposite situation. As reflected in Figure 1, as the rectangularisation continues to occur, the survival rates of more and more ages $x$ rise above 0.5 and get closer to 1 over time $t$, while the rates drop more and more sharply at the progressively narrowing highest end. This phenomenon may make one or more of the asymmetric candidates more suitable for modelling how the survival rates evolve over time. More details of the maxima and minima GEV distributions are given in the Appendix A.

In addition to ${}_{n}p{}_{x_0,t}$, we also follow Tan et al. (2016) and consider another response ${({}_{n}p{}_{x_0,t})}^{\frac{1}{n}}$, that is, the “annualised” survival probability for comparison. Figure 3 shows that the resulting annualised survive curve also displays an upward and rightward shift over successive time periods. Regarding the regression structure $h(x,t)$, we employ both the Lee and Carter (1992) structure and the Cairns et al. (2009) structure to allow for the age and period effects. The first one is taken as $h(x,t)=a_x+b_x{k}_t$, and the second one as $h(x,t)=k_{t,1}+k_{t,2}(x-\overline{x})$$+k_{t,3}({(x-\overline{x})}^2-{\sigma }^2)$, where $a_x$ is the age effect, $k_t$ is the “survival index” with age-specific sensitivity $b_x$, $k_{t,1}$ to $k_{t,3}$ are three time-varying parameters, $\overline{x}$ is the average of the age range, and ${\sigma }^2$ is the average of $(x-\overline{x})$ over the age range considered. Altogether, there are 5 (links) × 2 (responses) × 2 (structures) = 20 combinations under our consideration. This coverage is much more comprehensive than those of the few earlier papers on projecting the survival rates directly. Subsequently, we will also explore adding macroeconomic factors into the regression structure to see whether it can improve the performances. The Appendix A provides the parameter estimation methods for the models tested.

Note that some of these stochastic projection models can be seen as modifications of those standard static mortality or survival functions. Suppose ${\mu }_x$ and ${}_{t}p{}_x$ represent the force of mortality and the t-year survival rate of a life aged x for a static one-year life table. For instance, the Gompertz law states that ${\mu }_x=\alpha {\beta }^x$, which can be transformed into the survival function ${}_{t}p{}_x=\exp (-\alpha {\beta }^x({\beta }^t-1)/\ln \beta )$ and a regression structure $\ln {\mu }_x=\ln \alpha +x\ln \beta$. The old-age component of the Heligman–Pollard curve can be taken as $q_x/p_x=\alpha {\beta }^x$, which can be expressed as the survival function $p_x={(1+\alpha {\beta }^x)}^{-1}$ and a regression structure $\mathrm{logit}{q}_x=\ln \alpha +x\ln \beta$. The famous Lee–Carter model (Lee and Carter 1992) and the CBD model (Cairns et al. 2006) can be seen in some way as adding time-varying components into the two regression structures above (log and logit) in order to turn them from being static to stochastic. While the stochastic Lee–Carter and CBD models deal with mortality rates, as compared to our approach of modelling survival probabilities, we will include them in the following analysis for comparison.

3. Fitting Performances

In this section, we apply the 20 models (combinations) to female and male populations in Australia and New Zealand. The mortality data of ages 60 to 99 and years 1970 to 2017 of Australia (1970 to 2013 for New Zealand) are extracted from the Human Mortality Database (HMD 2020). For demonstration purposes, we consider the survival probabilities for a life aged 60 and set $x_0=60$ (for surviving periods $n$ = 1, 2, 3, …, 40), as the mortality rates below age 60 have already reached very low levels and have relatively little impact on the overall survival rates. In fact, the survival rate of a newborn for a surviving period up to 60 years is now very close to one, and the corresponding movements over time look too insignificant for a meaningful projection. By contrast, there is still much room for the survival rates at ages 60 and above to improve, and as shown below, the resulting patterns and trends are rather stable and can readily be projected into the future. Moreover, longevity of retirees is a serious concern for governments, insurers selling annuities, and pension funds because of the increasing financial burden. It would be of high practical interest to focus on the mortality of retirement ages.

Figure 4 plots the survival index $k_t$ of the first regression structure and also the time-varying parameters $k_{t,1}$ to $k_{t,3}$ of the second regression structure based on the five different link functions for the survival probabilities of Australian females. It is interesting to see that all the temporal parameters of $k_t$ and $k_{t,1}$ demonstrate a strong linearly increasing trend. It reflects clearly the continually improving survival rates at old ages as a whole. (For the complementary log-log model structure, Figure 2 displays a negative relationship between the response and the argument in contrast to the others. The resulting major trends of $k_t$ and $k_{t,1}$ in Figure 4 are then inverted, but the implication on mortality decline is the same.) The time-varying parameters $k_{t,2}$ and $k_{t,3}$ of the second regression structure refer to the slope and curvature in year $t$. While the directions of their trends are different because of the differences in how the link functions operate on the survival rates, a high level of linearity can largely be observed. We can then use the (univariate or multivariate) random walk with drift to project all these linear trends. Compared to the time-varying trends usually seen when modelling the mortality rates (e.g., Cairns et al. 2009), modelling the survival rates here generates more linear trends, which make the use of the random walk with drift more justifiable than otherwise. The observations for Australian male and New Zealand populations (not shown here) exhibit similar patterns.

Table 1. MAPE values (%) on fitted survival probabilities at age 60 from fitting 22 models to Australian and New Zealand data (1970–2017 and 1970–2013 respectively).

Model	Australian Females		Australian Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
$\mathrm{probit}-k_t$	1.24	1.47	2.03	2.53
$\log \text{-}\log -k_t$	1.59	1.59	2.75	2.75
$\mathrm{logit}-k_t$	1.39	1.59	1.94	2.72
$\mathrm{gevit}-k_t$	1.11	1.09	1.91	2.00
$\mathrm{gevmin}-k_t$	1.12	1.14	1.82	1.87
$\mathrm{probit}-k_{t,1}$$-k_{t,3}$	4.38	1.14	5.76	1.12
$\log \text{-}\log -k_{t,1}$$-k_{t,3}$	7.20	2.17	11.74	1.80
$\mathrm{logit}-k_{t,1}$$-k_{t,3}$	8.49	1.96	11.76	1.60
$\mathrm{gevit}-k_{t,1}$$-k_{t,3}$	2.71	0.85	3.48	1.01
$\mathrm{gevmin}-k_{t,1}$$-k_{t,3}$	2.95	0.73	3.71	0.93
Lee-Carter	1.39		2.42
CBD	1.62		1.57
Model	New Zealand Females		New Zealand Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
$\mathrm{probit}-k_t$	2.23	2.61	2.76	3.11
$\log \text{-}\log -k_t$	2.82	2.82	3.27	3.27
$\mathrm{logit}-k_t$	2.36	2.81	2.77	3.24
$\mathrm{gevit}-k_t$	2.06	2.02	2.70	2.73
$\mathrm{gevmin}-k_t$	2.09	2.11	2.66	2.70
$\mathrm{probit}-k_{t,1}$$-k_{t,3}$	4.80	1.45	5.95	1.86
$\log \text{-}\log -k_{t,1}$$-k_{t,3}$	8.21	2.24	12.65	2.35
$\mathrm{logit}-k_{t,1}$$-k_{t,3}$	9.22	2.06	12.19	2.20
$\mathrm{gevit}-k_{t,1}$$-k_{t,3}$	3.09	1.26	3.53	1.71
$\mathrm{gevmin}-k_{t,1}$$-k_{t,3}$	3.38	1.20	3.71	1.68
Lee-Carter	2.44		3.02
CBD	2.00		2.15

reports the MAPE (mean absolute percentage error) values on ${}_{n}p{}_{60,t}$ of fitting the 20 models, with the original Lee–Carter model (Lee and Carter 1992) and the CBD model with curvature (Cairns et al. 2009) for the death rates also included for comparison. The major observations are stated below:

(1)

For the first regression structure, the MAPE values are often smaller when the response is ${}_{n}p{}_{x_0,t}$. By contrast, for the second regression structure, the MAPEs are much smaller when the response is ${({}_{n}p{}_{x_0,t})}^{\frac{1}{n}}$.

(2)

When the response is ${}_{n}p{}_{x_0,t}$, the MAPE values from the first regression structure are clearly smaller. However, when the response is ${({}_{n}p{}_{x_0,t})}^{\frac{1}{n}}$, the situation is mostly reversed.

(3)

For the complementary log-log model structure with the first regression form, the MAPE remains the same regardless of the response (i.e., ${}_{n}p{}_{x_0,t}$ or ${({}_{n}p{}_{x_0,t})}^{\frac{1}{n}}$). The underlying reason is that $\ln (-\ln ({({}_{n}p{}_{x_0,t})}^{\frac{1}{n}}))=a_x+b_x{k}_t$ is indeed equivalent to $\ln (-\ln ({}_{n}p{}_{x_0,t}))$$=a_x+b_x{k}_t+\ln n$$=a_x^{\ast }+b_x{k}_t$, where $a_x^{\ast }=a_x+\ln n$. Hence, they produce the same fitted values of ${}_{n}p{}_{60,t}$ and so the same MAPEs.

(4)

Overall, the combination of the response ${({}_{n}p{}_{x_0,t})}^{\frac{1}{n}}$ and the second regression structure (i.e., $k_{t,1}+k_{t,2}(x-\overline{x})$$+k_{t,3}({(x-\overline{x})}^2-{\sigma }^2)$) gives the better MAPEs. In particular, the gevmin model structure, based on the newly proposed gevmin link, leads to the smallest MAPEs consistently for all the populations (0.73, 0.93, 1.20, 1.68) considered.

(5)

The best gevmin model structures noted above outperform the more traditional approaches of the Lee–Carter and CBD models in modelling the mortality rates (with MAPEs of 1.39 to 3.02).

Compared with the traditional link functions, the additional shape parameter $\xi$ in the gevmin and gevit link functions makes them a lot more flexible in capturing different degrees of asymmetry. For each population in Table 1, the two smallest MAPEs are all generated from either the gevmin or gevit model structure. The empirical advantage of using these two newer links is apparent here. Note also that the response ${({}_{n}p{}_{x_0,t})}^{\frac{1}{n}}$ (see Figure 3 again) has a simpler shape over $n$ than the response ${}_{n}p{}_{x_0,t}$ (see Figure 1). Hence, a simple regression structure in terms of merely $x$ (i.e., $k_{t,1}+k_{t,2}(x-\overline{x})$$+k_{t,3}({(x-\overline{x})}^2-{\sigma }^2)$) would suffice for the former, while a more dedicated allowance for the age effect (i.e., $a_x+b_x{k}_t$) would be needed for the latter.

As noted earlier, if the new gevmin link is applied on $1-{}_{n}p{}_{x_0,t}$ rather than ${}_{n}p{}_{x_0,t}$, the model structure turns into $[1-{(-\ln ({}_{n}p{}_{x_0,t}))}^{-\xi }]/\xi =h(x,t)$. It can further be rearranged as $[{(-\ln ({}_{n}p{}_{x_0,t}))}^{-\xi }-1]/\xi =-h(x,t)=h^{\ast }(x,t)$, which then becomes the gevit model structure effectively. Consequently, they generate the same fitted values of ${}_{n}p{}_{60,t}$ and the same MAPEs, though the resulting signs and trends of their parameters are opposite to each other. However, the increasing survival index based on the gevmin link on ${}_{n}p{}_{x_0,t}$ has a more natural interpretation in terms of improving survival over time.

4. Forecasting Performances

In this section, we perform an out-of-sample analysis to assess the forecast accuracy of the survival rate projection models. We apply the models to four fitting periods of 1970 to 1989 (20 years), 1970 to 1994 (25 years), 1970 to 1999 (30 years), and 1970 to 2004 (35 years) and then forecast the survival rates for the remaining periods until the very last year of available data. Here we focus on the annualised response and the second regression structure since they give the better fitting performances in Table 1. In addition to continuing to use the Lee–Carter and CBD models for comparison, we also apply the multivariate random walk with drift to the survival probabilities as a benchmark, naive model. The MAPEs of not only the projected ${}_{n}p{}_{60,t}$ but also the projected life expectancies at age 60 are examined to compare the performances. They are calculated as $\frac{1}{n_d}{\displaystyle {\sum }_{x,t}|{}_n\widehat{p}{}_{60,t}}-{}_n\dot{p}{}_{60,t}|/{}_n\dot{p}{}_{60,t}$ and $\frac{1}{n_t}{\displaystyle {\sum }_t|{\widehat{e}}_{60,t}}-{\dot{e}}_{60,t}|/{\dot{e}}_{60,t}$ respectively, where ${}_n\widehat{p}{}_{60,t}$ and ${}_n\dot{p}{}_{60,t}$ (${\widehat{e}}_{60,t}$ and ${\dot{e}}_{60,t}$) are the projected and observed survival probabilities (life expectancies) in year $t$, $n_d$ is the number of data points, and $n_t$ is the number of years in the testing period. For each case (column) in

Table 2. MAPE values (%) on projected survival probabilities at age 60 using Australian and New Zealand data for four different fitting periods.

	Australian Females				Australian Males
Model	1970–1989	1970–1994	1970–1999	1970–2004	1970–1989	1970–1994	1970–1999	1970–2004
probit	5.25	1.87	3.03	3.03	14.65	9.92	2.70	1.65
log-log	6.08	2.05	2.65	2.02	15.45	10.79	3.49	1.54
logit	5.95	1.97	2.64	2.77	15.30	10.63	3.30	1.97
gevit	4.43	2.61	4.21	3.94	13.59	8.91	2.42	2.08
gevmin	4.42	2.44	3.96	3.72	13.98	9.26	2.40	1.77
LC	3.96	3.31	3.24	2.41	10.65	7.51	3.97	3.20
CBD	4.62	1.87	3.28	3.51	15.45	11.01	3.75	2.39
MRW	8.61	4.44	2.51	2.05	18.05	14.03	7.70	5.50
	New Zealand Females				New Zealand Males
Model	1970–1989	1970–1994	1970–1999	1970–2004	1970–1989	1970–1994	1970–1999	1970–2004
probit	4.78	5.70	3.34	2.77	13.23	8.48	9.53	8.15
log-log	5.13	5.04	3.53	3.23	13.65	8.89	9.97	8.47
logit	5.04	5.15	3.47	3.12	13.57	8.81	9.88	8.41
gevit	4.79	6.49	3.54	2.88	12.57	7.95	9.03	7.94
gevmin	4.80	6.47	3.53	2.85	12.81	8.17	9.20	7.97
LC	6.60	7.28	5.09	3.56	15.95	9.77	8.70	6.26
CBD	5.66	8.10	3.77	2.90	14.09	9.19	9.99	7.77
MRW	6.61	3.28	3.60	3.37	15.60	10.81	11.25	9.74

, the three lowest MAPE values are highlighted (in grey). We notice the following patterns within the results in Table 2:

(1)

Out of the 16 cases (4 fitting periods × 2 countries × 2 sexes), the gevit and gevmin model structures produce the three lowest MAPEs in 12 cases. Their performances are the most consistent ones amongst all the candidates.

(2)

For the gevit and gevmin model structures, the average MAPE is 6.10. For the probit, complementary log-log, and logit model structures, the average MAPE is 6.29. For the LC and CBD models, the average MAPE is 6.53. For the naive random walk model, the average MAPE is 7.95.

(3)

The MAPE values tend to be lower for females (4.04 on average) than for males (8.98 on average).

(4)

The MAPE values tend to be lower for Australia (5.65 on average) than for New Zealand (7.36 on average).

(5)

The naive random walk model leads to the highest MAPE values in 10 cases.

Overall, the gevit and gevmin model structures deliver the best and the most consistent forecasting performances among the eight alternatives. Moreover, the survival rate projection models as a whole also tend to outperform the usual Lee–Carter and CBD models in this out-of-sample analysis of the survival rates. These results highlight the potential usefulness of projecting the survival rates directly when the focus is on the survival probabilities. The generally poor performances by the naive random walk model also emphasise the importance of having a proper model structure for modelling the survival rates.

Regarding the MAPEs of the projected life expectancies in

Table 3. MAPE values (%) on projected life expectancies at age 60 using Australian and New Zealand data for four different fitting periods.

	Australian Females				Australian Males
Model	1970–1989	1970–1994	1970–1999	1970–2004	1970–1989	1970–1994	1970–1999	1970–2004
probit	1.67	0.69	1.12	1.20	5.34	3.68	0.84	0.55
log-log	2.00	0.71	0.88	1.03	5.90	4.23	1.34	0.39
logit	1.97	0.71	0.90	1.05	5.83	4.17	1.28	0.40
gevit	1.24	0.95	1.61	1.52	4.34	2.70	0.66	1.03
gevmin	1.34	0.82	1.42	1.39	4.85	3.22	0.66	0.77
LC	2.04	0.74	0.79	0.82	6.11	4.13	1.54	0.70
CBD	2.23	0.78	0.78	1.00	6.39	4.52	1.69	0.51
MRW	1.33	0.79	1.43	1.49	5.05	3.40	0.74	0.68
	New Zealand Females				New Zealand Males
Model	1970–1989	1970–1994	1970–1999	1970–2004	1970–1989	1970–1994	1970–1999	1970–2004
probit	2.83	0.75	1.78	1.17	7.16	4.29	4.09	2.20
log-log	2.93	0.81	1.83	1.14	7.43	4.57	4.31	2.33
logit	2.92	0.80	1.83	1.15	7.39	4.53	4.29	2.32
gevit	2.63	0.78	1.64	1.15	6.60	3.78	3.70	1.91
gevmin	2.70	0.75	1.70	1.18	6.86	4.03	3.90	2.07
LC	3.50	0.91	2.38	1.12	9.32	5.28	4.52	2.24
CBD	3.12	0.80	1.78	1.16	7.34	4.78	4.36	2.54
MRW	3.13	0.79	1.74	1.14	7.18	4.25	3.97	2.00

, the gevit and gevmin model structures continue to perform well relative to the others. Their average MAPE is 2.31, compared to 2.56 for the probit, complementary log-log, and logit model structures and 2.81 for the LC and CBD models. It highlights again the flexibility offered by the extra shape parameter and the potential benefits. Moreover, the survival rate projection models still perform better than the Lee–Carter and CBD models in general. Projecting the survival rates directly would be a useful strategy or alternative if the main task is to forecast future life expectancies. It is also interesting to note that the average MAPE for the naive random walk model is only 2.44. We conjecture that their larger errors as shown in Table 2 may have coincidently offset one another to some extent across different ages as life expectancy is effectively an aggregate measure of survival probabilities. Anyway, the advantages of imposing an appropriate model structure are obvious in modelling the survival rates.

In addition to the MAPE values, the sMAPE (symmetric mean absolute percentage error) values are also calculated and given in the Appendix A. The major implications are largely the same as those discussed above.

5. The Effect of Economic Growth

There have been previous studies in the area of demography and macroeconomics showing the impact of economic fluctuations on mortality levels (e.g., Ruhm 2000; Brenner 2005). A few earlier papers in mortality projection that incorporated exogenous economic factors include Hanewald (2011), Niu and Melenberg (2014), and Boonen and Li (2017). They showed that the Gross Domestic Product (GDP) of a nation may serve as an explanatory factor of a country’s mortality rates. Furthermore, they demonstrated that it may be integrated into extrapolative mortality models such as the Lee–Carter model or its extensions to improve the fitting quality. In this section, we experiment with embedding such macroeconomic factor into the survival rate projection models. As a preliminary analysis, we find some positive correlation ($\rho$ = 0.17 for Australian females, 0.10 for Australian males, 0.08 for New Zealand females, 0.15 for New Zealand males) between the annual change in survival index $k_t$ and the annual growth in real GDP per capita for the period from 1970 onwards (for Australia and New Zealand respectively). The survival index $k_t$ can be seen as an underlying driver of the survival probabilities and so an indicator of how life expectancy changes, and the GDP trend is a very important indicator of economic growth. This observation provides some incentive for putting the survival trend and the economic trend together into the modelling framework. By doing so, the extended model may produce more interpretable implications on the projected survival rates based on the projected economic growth. In effect, such a model can be regarded as partly extrapolative (using past survival trends) and partly explanatory (using exogenous GDP figures). Note that the GDP can be regarded as an overall indicator of the quality of life. As the focus of our study is the use of new link functions in modelling survival probabilities and the aggregate data used are of country level, we deem that the GPD itself would be sufficient for our analysis1.

Accordingly, we modify the first regression structure as $h(x,t)=a_x+b_x{k}_t+c_x{g}_t$, in which $g_t$ is the (log) real GDP per capita in year $t$ and $c_x$ is the (age-specific) associated coefficient. The values of $g_t$ are adjusted such that ${\displaystyle {\sum }_t{g}_t}=0$, and so $a_x$ still refers to the overall (average over time) age effect. The real GDP per capita data are obtained from World Development Indicators (WDI 2020), the primary database maintained by World Bank on comparative development indicators across countries worldwide. Figure 5 (upper panel) illustrates the increasing trend of the GDP in recent decades, as well as some potential co-movements with the life expectancy trend. Following the Box and Jenkins (1976) approach, we find that the ARIMA(1,1,0) process (autoregressive integrated moving-average) provides an adequate description of the dynamics of the GDP process. Figure 5 (lower panel) also shows that the fitted GDP values from ARIMA(1,1,0) are very close to the observed values.

Figure 6 plots the survival indices under the modified (with GDP) first regression structure with the gevit and gevmin links. While the new factor $c_x{g}_t$ is supposed to explain part of the rising survival trend, the computed survival indices still exhibit a clear linearly increasing trend, though their drifts (around 0.37) are smaller than previously (about 0.43 without the GDP factor). It is also interesting to see in Figure 6 that the age-sensitivity $c_x$ is high at ages 60 to 80 but it drops drastically after around age 80, reflecting that the very-old-age mortality is still much less responsive to economic growth (likely due to current technological and medical limitations).

Table 4. MAPE values (%) on fitted survival probabilities at age 60 from fitting 2 models with real GDP per capita as explanatory variable to Australian and New Zealand data (1970–2017 and 1970–2013 respectively).

Model	Australian Females		Australian Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
with GDP
$\mathrm{gevit}-k_t$	0.86	0.79	1.69	1.61
$\mathrm{gevmin}-k_t$	0.87	0.73	1.66	1.41
without GDP
$\mathrm{gevit}-k_t$	1.11	1.09	1.91	2.00
$\mathrm{gevmin}-k_t$	1.12	1.14	1.82	1.87
Model	New Zealand Females		New Zealand Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
with GDP
$\mathrm{gevit}-k_t$	1.76	1.37	2.54	2.55
$\mathrm{gevmin}-k_t$	1.70	1.40	2.52	2.50
without GDP
$\mathrm{gevit}-k_t$	2.06	2.02	2.70	2.73
$\mathrm{gevmin}-k_t$	2.09	2.11	2.66	2.70

compares the fitting performances with and without the GDP factor. Interestingly, all the MAPE values become smaller after incorporating the GDP covariate. Furthermore,

Table 5. MAPE values (%) on projected survival probabilities at age 60 from using 2 models with real GDP per capita as explanatory variable on Australian and New Zealand data (left value—fitting period 1970–1989; right value—fitting period 1970–1999).

Model	Australian Females		Australian Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
with GDP
$\mathrm{gevit}-k_t$	4.31/4.67	4.12/6.21	11.67/3.35	11.65/4.79
$\mathrm{gevmin}-k_t$	3.99/2.88	4.06/6.74	12.27/3.80	11.79/4.20
without GDP
$\mathrm{gevit}-k_t$	4.89/5.62	5.21/8.02	12.50/2.81	9.41/7.21
$\mathrm{gevmin}-k_t$	5.24/5.45	5.83/8.08	12.49/2.14	9.53/5.97
Model	New Zealand Females		New Zealand Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
with GDP
$\mathrm{gevit}-k_t$	6.79/4.18	5.94/4.34	12.93/8.67	12.47/7.00
$\mathrm{gevmin}-k_t$	6.53/3.71	6.37/4.94	12.60/8.67	12.64/6.95
without GDP
$\mathrm{gevit}-k_t$	7.38/5.32	7.51/5.36	13.19/8.60	12.34/7.52
$\mathrm{gevmin}-k_t$	7.25/5.23	6.68/5.73	13.00/8.51	12.49/7.53

and

Table 6. MAPE values (%) on projected life expectancies at age 60 from using 2 models with real GDP per capita as explanatory variable on Australian and New Zealand data (left value—fitting period 1970–1989; right value—fitting period 1970–1999).

Model	Australian Females		Australian Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
with GDP
$\mathrm{gevit}-k_t$	2.16/1.78	1.61/3.34	5.81/0.91	3.07/3.90
$\mathrm{gevmin}-k_t$	2.07/1.29	1.55/3.18	5.56/0.58	4.11/2.83
without GDP
$\mathrm{gevit}-k_t$	2.55/2.04	2.27/3.48	7.16/1.49	4.94/4.96
$\mathrm{gevmin}-k_t$	2.76/1.76	2.24/2.83	7.35/0.94	5.74/2.79
Model	New Zealand Females		New Zealand Males
	Non-Annualised	Annualised	Non-Annualised	Annualised
with GDP
$\mathrm{gevit}-k_t$	3.71/1.87	2.61/1.14	7.91/4.23	7.03/2.20
$\mathrm{gevmin}-k_t$	3.70/1.76	2.67/1.03	7.38/4.11	7.34/2.72
without GDP
$\mathrm{gevit}-k_t$	4.01/2.32	3.57/1.97	8.28/4.67	7.65/3.87
$\mathrm{gevmin}-k_t$	4.00/2.27	3.37/2.08	8.16/4.63	7.86/4.16

present the forecasting performances before and after embedding the GDP covariate under the same out-of-sample test setting as in Section 4. For both the projected survival probabilities and the projected life expectancies at age 60, adding the GDP factor leads to a clear improvement in forecast accuracy in most of the cases. It is obvious that besides the selection of appropriate link function, age and period effects, and annualisation of the survival rates, integrating a GDP explanatory factor within the model has much potential in further enhancing the performance of a survival rate projection model.

A final note is that while our regression structure here is designed to examine how the economic trend would affect and/or move synchronously along with the survival trend over the long term, we acknowledge that the underlying process can be more complicated in terms of the direction of causation and the lagging effect. These matters are beyond the scope of this paper and we leave them for future research.

6. Concluding Remarks

In this paper, we have conducted a thorough examination on several old and new link functions for their applications in modelling how old-age survival probabilities evolve across time in Australasia. The link functions under investigation include the newly proposed gevit and gevmin links, which are compared against the traditional ones like probit, complementary log-log, and logit links. We use the random walk with drift to project the temporal parameters due to their highly linear trends in the past decades. Future survival probabilities and life expectancies are then projected using the estimated parameters. We notice that many of these survival rate projection models can produce better fitting and forecasting performances than the more conventional approach of modelling mortality rates, in terms of achieving lower MAPE and sMAPE values. Hence, projecting the survival rates directly would serve as a useful strategy or alternative if the objective is to forecast future life expectancies. In particular, the gevmin and gevit links are found to be able to offer extra flexibility in depicting the (annualised or non-annualised) survival curve patterns and improve the model performance. This extra flexibility comes from the additional shape parameter $\xi$ that copes with any extent of asymmetry. Lastly, we illustrate how these survival rate projection models can be modified to incorporate additional explanatory variables, in an attempt to further enhance the model results. We notice that adding a covariate on the real GDP per capita can improve the general fitting performances as well as many of the forecasting performances. Note that our approach can be applied similarly to other ages above 60 for $x_0$. The focus of this work is on forecasting survival probabilities at ages 60 and above, as the current death rates below age 60 are very low and their impact on the overall survival rates has become much less important. If one is interested in modelling the entire age range, the Lee–Carter model and its various extensions in modelling mortality rates can be used instead.

There are two potential areas for future research. First, since a link function is often constructed from the inverse of a cdf, the use of other cdfs can be explored. There are a range of skewed distributions which may be useful for tackling different levels of asymmetry. Second, while we have used two particular regression forms in our analysis, there are other possible structures that also allow for the age and period effects. Some examples include adding extra time-varying parameters and co-modelling the survival rates of related subpopulations.