ARMA MODELS FOR MORTALITY FORECAST

Abstract. In the last several decades, many countries have been paying a lot of attention to mortality forecasting because of high longevity risk. The purpose of this paper is to analyze mortality characteristics of Baltic countries and make predictions using ARMA models. Research showed that mortality rate distribution is almost the same in Lithuania, Latvia and Estonia and all of them represent longevity trends. It means that men and women, children and adults have the same mortality structure in all Baltic countries and live longer than before.


Introduction
Time series analysis is one of the most popular methods to look through the data of different types with alternate characteristics.We will use an autoregressive moving average (ARMA) model as it helps to represent the data structure of the mortality rate in Baltic countries.The main objectives of this research are to analyze mortality rate data and its structure, compare mortality distribution by sex and age.We will answer the question whether mortality rate distribution is the same in all Baltic countries and forecast future mortality rates in case ARMA models work correctly.This information and the results obtained could be used by the social insurance systems of Baltic countries, insurance companies and pension funds as longevity trends can cause a lot of problems to their businesses.Many attempts were made to understand and explain these changes in different countries.Some of them are discussed by Tan K.S., Blake D. and MacMinn R. (see [11]).The application of AR-ARCH models for mortality analysis were discussed by Giacometti R., Bertocchi M., Rahev S.T. and Fabozzi F.J. (see [4]).An in-depth longevity analysis can represent changes in recent mortality rates and help make accurate forecasts.
The paper is organized in the following way.In section 2, we describe the main tools which we use for the mortality data analysis.Section 3 is dealt with the mortality forecast for populations of Lithuania, Latvia and Estonia.In sections 4, we compare the empirical data with the forecast results.Finally, in section 5, we present calculations of the remaining life-time for the population members under condition that mortality rates vary according to ARMA processes.

Autoregressive models
Let X = (X t ,t ∈ T ), T ⊂ R, be a stochastic process.Then X t with EX 2 t < ∞, t ∈ T , is said to be (weakly) stationary if, for all t 1 ,t 2 ∈ T and for such h > 0 that t 1 + h,t 2 + h ∈ T , two following equalities are satisfied: To analyze stationary time series, different models are used.One of the most popular models is a discretetime stationary process called autoregressive moving average process (ARMA).The model ARMA(p, q) is a generalization of AR(p) and MA(q) models.The process X t , t ∈ Z, is an autoregressive process AR(p) of order p if X t is stationary and where µ and α 1 , . . ., α p are parameters.Here and everywhere W t ∼ W N(0, σ 2 ) is a stochastic process called white noise that has further characteristics: It means that model AR(p) tries to approximate the process value at the moment t as linear dependence on the previous values.The process X t , t ∈ Z, is called the moving average process MA(q) of order q if where ν and β 1 , . . ., β q are parameters and W t ∼ W N(0, σ 2 ).Process X t , t ∈ Z, is called the ARMA(p, q) process with mean a if X t is stationary and ARMA models are widely used for almost all types of data as they not only explain regression in data but also add a stochastic factor to a model.

Tests for stationarity
As we are going to use ARMA models, we should check whether our data is suitable for this kind of analysis, i.e. whether the time series is stationary.To make such analysis, we use two different tests.
The first one is an augmented Dickey-Fuller test.The ADF test examines unit root existence for the process with test hypothesis The null hypothesis H 0 states that there is a unit root, so time series is not stationary.To find the answer, we use test statistic γ/SE(γ) (for details, see [2]) and then compare the result with the relevant critical value for the Dickey-Fuller test.If the test statistic is less than the critical value then the null hypothesis is rejected.Sometimes the ADF test states that analyzed data has a unit root.However, it is not the final answer to our question because the ADF test suggests only an approximate p-value and may be inaccurate.The shortest way to check results is to test data with different tests.We will do this with Kwiatkowski-Phillips-Schmidt-Shin (KPSS) test.
The KPSS test is very practical as modeling the ARMA process with a function auto.arima in R software this test is done automatically.Using the KPSS test it is important not to confuse given results.The null hypothesis H 0 in the KPSS test is opposite to that in the ADF test.Hypothesis H 0 of the KPSS test states that the process is stationary.For more information about the theoretical base and practical usage of these tests, see [2] and [6].So, we will use both ADF and KPSS tests in order to check whether our analyzed data is stationary.The confidence level selected for our tests is 95%.We will require that at least one of the tests results states that the analysed data is stationary.

ARMA model selection. AIC criterion
The main purpose of the analysis is to select the ARMA(p, q) model that approximates empirical data best.To achieve our objectives, we will use the Akaike information criterion (AIC): where σ2 W is the estimation of the white noise variance σ 2 W , calculated by the maximum likelihood method.It ensures that the optimized ARMA(p, q) model is selected according to the smallest AIC value.Using the function auto.arima in R software, we do not need to check the criterion as it selects the best model for the analyzed data automatically.As a result, we have ARMA(p, q) models where complexity and accuracy match best.

Empirical data analysis
Historical mortality data of three Baltic countries -Lithuania, Latvia, and Estonia -is represented in this research.All the required mortality tables have been taken from the Human Mortality Database, University of California [12].Baltic countries data is collected from year 1959 to 2013.1959-1979 years of data will be ignored because of high probability of inaccuracy; 1980-2008 years of data will be used for analysis and predictions.Consequently, 2009-2013 years empirical data will be compared with the predictions made.Research will be done for people aged 0-90.
It should be noted that modeling will be done for a traditional variable m x,t , called mortality rate: , where µ x,t is the empirical force of mortality.
Most graphs show the log mortality rate ln(m x,t ), as this form of data has a better level of representation differences.Figure 1 represents historical log mortality rates for the total population.From the first point of view, Baltic countries seem to have similar variance level and trends.In further analysis, all time series will be broken down by sex and year, i.e. we consider m (M) x,t for men and m From Figure 2, it can be seen that Latvian and Estonian data do not show a lot of fluctuation till before the year 1995, but then similar downward trends are presented that prove longevity existence.However, Lithuanian mortality data fluctuates depending on time without any exact trend.
Besides, if we look at Lithuanian men's and women's data separately (see Figure 3), we will see absolutely different graphs according to sex.So, we will try to find these changes and forecast them using ARMA models.

Forecast of m x,t for Lithuanian data
We will perform a minute analysis of 0-, 18-, 40-and 65-year-old people mortality rate series.All other results will be summarized and presented in graphs and tables.
Let's analyze newborn girls' data series m (W ) 0,t ,t = 1980, ..., 2008.In order to check whether the series is stationary, we will use KPSS test results.The KPSS test p-value is p > 0.1, so our analyzed data is considered to be stationary.As the stationarity requirement is fulfilled, we can apply an ARMA model.The suggested model is AR (1).Therefore, one step prediction is the following m(W) 0,t = 0.0095 + 0.8886m We also looked through the 18-year-old girls' series m 18,t ,t = 1980, ..., 2008.Again, the KPSS test suggests p-value p > 0.1, so we consider the data analyzed to be stationary.As a result, the model AR(1) is selected.Hence, m(W) 18,t = 0.0006 − 0.3268m The same analysis was done for 40-year-old women, and based on the ADF test result (p = 0.0196), the AR(1) model was selected.Therefore, m(W) 40,t = 0.0021 + 0.3761m The series of 65-year-old women's data is more complex as no stationary models were found.So, we are going to analyze time series of differences: ∆m The stationarity of 65-year-old men's data is stated only by the ADF test with p = 0.0279.We decided that this is enough to find the data approximation.In the case, the model selected is AR (1).So, according to this selection m(M) 65,t = 0.0353 + 0.639m The same analysis was done for people of different ages.Using suggested ARMA models for all the data, forecasts were made for the years 2009-2013.The results of ln( mx,t ) broken down by sex can be seen in Figure 4.

Forecast of m x,t for Latvian data
Again, we will analyze 0-, 18-, 40-and 65-year-old women's and men's data time series.We will start with women's data.Let us take data on newborn girls m (W ) 0,t ,t = 1980, . . ., 2008.This time series is considered to be stationary by KPSS test (p > 0.1).The selected model is AR (1).Hence, we can suppose that m(W) 0,t = 0.0107 + 0.8206m As series m (W ) 18,t ,t = 1980, . . ., 2008 is considered to be non-stationary by both tests, we will study time series of the data differences ∆m The stationarity of 40-year-old men's data is proven by both ADF (p = 0.027) and KPSS (p > 0.1) tests, and the AR(2) model is selected, we have that m(M) 40,t = 0.0076 + 1.0242m The AR(1) model is suggested for 65-year-old men's data and stationarity is checked by the KPSS test (p = 0.0914).Due to the obtained result we have m(M) 65,t = 0.0404 + 0.7064m After the same analysis is done for all the ages and appropriate ARMA models are chosen, 5-year forecasts for women's and men's data are presented in Figure 5.

Forecast of m x,t for Estonian data
We will repeat the same analysis for Estonian mortality data and present some of the results in more detail.Let us take time series of newborn girls.As m (W ) 0,t ,t = 1980, . . ., 2008 is not stationary itself, we will look through its differences.Time series ∆m (W ) 0,t ,t = 1981, . . ., 2008 is considered to be stationary by both ADF (p < 0.01) and KPSS (p > 0.01) tests.So, after processing of the data we select the MA(1).So we have that ∆ m(W) 0,t = 0.0004.
Then we will analyze 18-year-old women's series of differences.The p-value of the KPSS test for data 18,t is p > 0.1, and the ADF test suggests p = 0.0194.So, time series of differences is stationary, the best model for data is MA(1), and, therefore, ∆ m(W) 18,t = 0.00094.The AR(1) model is selected for 40-year-old women's data.The stationarity is checked by the ADF test (p = 0.0173).Therefore, the following formula should be used for one step prediction m(W) 40,t = 0.002 + 0.4067m The 65-year-old women's time series is considered to be stationary by both tests.The p-value of the ADF test is 0.015, and p-value of the KPSS test is p > 0.1.As a result, the ARMA ( The stationarity of 65-year-old men's data is proven by the KPSS test (p > 0.0894).As the stationarity condition is fulfilled, we can proceed the approximation by some ARMA model.As a result, the AR(1) model is selected, and we can suppose that m(M) 65,t = 0.0386 + 0.323m The forecasts for both women's and men's data are presented in Figure 6.With reference to all the results and forecasts made for Baltic countries, the following findings were made: (i) the infant mortality rate is relatively high; (ii) the mortality rate distribution for women and men in all Baltic countries is almost the same; (iii) a considerable difference is seen in women's and men's mortality rates of average years: women's data surface is convex upward, men's -downward.

Comparison of forecasts
As it was stated before, a comparison of children's and adults' data would be made separately.So, we use empirical data for the years 2009-2013 and compare it with forecasts made by ARMA models.We calculate the mean square error by the formula: In our case, For better understanding of the results, let us look at figures 9 and 10, where a comparison of the empirical mortality rate and forecasts is presented for the years 2009 and 2013.According to the information that can be seen from the graphs, we can make some insights: Forecasts for the year 2009 have a high level of accuracy, while forecasts for the year 2013 are not so precise, in comparison with empirical mortality data.This is a typical situation.As information is aging and forecasts depend on previous values, the results grow worse as the number of forecast steps increases; (ii) During analysis, it was seen that women's and men's results have the same structure of √ MSE.As time increases, we see more errors.Still ,the situation of children's data is even more complex.
√ MSE calculated from the children mortality data is much higher than for adults.Moreover, the accuracy of forecast has a strong fluctuation depending on time.
(iii) ARMA models predict adult mortality rate more accurately, in comparison with results for children.Now we are going to discuss forecasts for the year 2016.From Figure 11, we can see that mortality rates in all Baltic countries are identical for all cuts.Again, boys' and men's mortality rates are much higher than girls' and women's respectively.The most significant difference is noticed in adults data forecasts.
Unfortunately, after all the analysis, it was found that no long-term forecasts could be done.They go to the average too fast because of small ARMA model order.So, forecasts for the years 2016 and 2026 have almost no differences.In order to represent longer-term forecasts with a higher level of accuracy, we must use ARMA models of greater order or choose other time series models.

Average time of life remaining
In this chapter, we will use forecasts for the year 2016 to calculate the average time of life remaining for individuals of different sex and age, i.e. e (W ) x,2016 , e As k p x,t = p x,t p x+1,t+1 p x+2,t+2 . . .p x+k−1,t+k−1 , for all possible x, k,t and p x,t = 1 − q x,t ≈ 1 − m x,t we have: êx,2016 All the results given are shown in Table 2 and Figure 12.It is one more approval that mortality distribution in Baltic countries is very similar and women's and men's mortality rates have differences.In addition, these results for the Lithuanian population were compared with the data presented in the research ( [5]).It was noticed that forecasts for Lithuanian data are very similar.Even more, Lithuanian men live by almost 1 year longer in the year 2016 than 4 years ago.It is one more proof that the average number of years remaining is still growing for Lithuanian people.

Concluding remarks
Mortality rates of all Baltic countries were compared and forecasts were presented.After data analysis, it can be proven that the mortality rate distribution of three Baltic countries is almost identical.Latvian population lives a little bit shorter, Estonian population lives a little bit longer, in comparison with the Lithuanian population.The men's death rate is higher than women's death rate; children's mortality rate is lower than infant, but it rises as children grow up.In general, children's data has a strong fluctuation depending on age, while adults' mortality rate faces a visible upward trend.
During the research, some problems were met.Time series for every year's category was not so long: it had only 29 observations and is not so much in order to find a stationary ARMA model that fits the data best.As data had a high level of fluctuation, especially in recent years, the accuracy of the forecasts dropped.As a result, we had quite big mean square errors, especially in children's data analysis.For the same reason, long time predictions do not tell us any additional information, which could be reliable for the strict statements.
In order to achieve a high level of forecast accuracy, different ARMA model modifications should be used, or one may approximate mortality rate series by another models.Nevertheless, it was enough to state the longevity trend in all Baltic countries and see the forecasts still representing it.

Figure 9 .Figure 10 .
Figure 9.Comparison of the empirical log mortality rate and forecast for the year 2009
−1 ,t = 1981, . . ., 2008.This time series is considered to be stationary by both tests (the p-value of the ADF test is p = 0.027 and the KPSS test p-value satisfies p > 0.1).The suggested model is the ARMA(1,1) model.Consequently, one step prediction is Now let us turn to time series analysis for men's data.The AR(1) model is selected for newborn boys, as both ADF (p = 0.0254) and KPSS (p > 0.1) tests state stationarity of the data.According to such selection (1)ng analysis for 18-year-old boys' data we get that the KPSS test is suggested with p > 0.1, and the model selected is MA(1).So, m(M) 18,t = 0.0019.The stationarity of 40-year-old men's data is proven by both ADF (p = 0.0351) and KPSS (p > 0.1) tests and the model suggested is AR(2).Hence, (1,1) = 1981, ..., 2008.This series is stationary as the KPSS test suggested p > 0.1 and p-value of the ADF test is p = 0.0187.The model ARMA(2, 1) is chosen for the data.The analysis of 40-year-old women's data shows that this time series is stationary.It is proven by both tests: the ADF test p-value is 0.0367 and KPSS test p-value is p > 0.1.The selected model is AR(1).2008 is not stationary, so let us look through the time series differences which are considered to be stationary by the ADF test (p = 0.0194).Then, the model chosen is ARMA(1,1).According to this we have ,t = 1981, . . ., 2008.This time series is stationary because the ADF test suggests pvalue 0.0321, and the KPSS test p-value is p > 0.1.So, after processing ARMA(1,1) model we get that (W ) 65,t ,t = 1980, . . ., (M) 18,t ,t = 1980, . . ., 2008 is not stationary, so we will analyze data differences ∆m

Table 2 .
Estimation of the average number of years remaining e x,2016 for Baltic countries Figure 12.Estimation of the average number of years remaining êx,2016