How to estimate mortality trends from grouped vital statistics.

BACKGROUND: Mortality data at the population level are often aggregated in age classes, for example 5-year age groups with an open-ended interval for the elderly aged 85+. Capturing detailed age-specific mortality patterns and mortality time trends from such coarsely grouped data can be problematic at older ages, especially where open-ended intervals are used.
METHODS: We illustrate the penalized composite link model (PCLM) for ungrouping to model cancer mortality surfaces. Smooth age-specific distributions from data grouped in age classes of adjacent calendar years were estimated by constructing a two-dimensional regression, based on B-splines, and maximizing a penalized likelihood. We show the applicability of the proposed model, analysing age-at-death distributions from cancers of all sites in Denmark from 1980 to 2014. Data were retrieved from the Danish Cancer Society and the Human Mortality Database.
RESULTS: The main trends captured by PCLM are: (i) a decrease in cancer mortality rates after the 1990s for ages 50-75; (ii) a decrease in cancer mortality in later cohorts for young ages, especially, and very advanced ages. Comparing the raw data by single year of age, with the PCLM-ungrouped distributions, we clearly illustrate that the model fits the data with a high level of accuracy.
CONCLUSIONS: The PCLM produces detailed smooth mortality surfaces from death counts observed in coarse age groups with modest assumptions, that is Poisson distributed counts and smoothness of the estimated distribution. Hence, the method has great potential for use within epidemiological research when information is to be gained from aggregated data, because it avoids strict assumptions about the actual distributional shape.

Introduction

Vital statistics data often lack fine resolution when publicly released. Death counts are usually provided in aggregated form. Databases of the World Health Organization (WHO) and NORDCAN, a Nordic tool for cancer information, planning, quality control and research, are two examples in which deaths are grouped in age classes, for example age groups of 5 years with a coarse open-ended interval that gathers together the elderly aged 85 and over. The reasons for such aggregation might be privacy protection and the wish for a compact and easy illustration. With increasing longevity in Western countries, wide grouping schemes at older ages are particularly problematic for analysis, given that an increasing proportion of the population dies in the given open-ended age group, for example people aged 85 and over. Records from the Human Mortality Database (HMD) show, for example, that deaths in Denmark of people aged 85 and over rose from 29% of total deaths in the year 2000 up to almost 35% in the year 2014. The analysis of mortality and leading causes of death at older ages is therefore problematic, not only because of the known variation in mortality reporting by physicians, but also because of coarse age groupings. Additionally, grouping schemes might differ across time and countries, which makes comparability challenging.

The most widely used methods to estimate univariate distributions, for example age-specific deaths for a single calendar year, from grouped data are flexible parametric models (with, however, a high number of parameters that lead to overfitting) and non-parametric methods such as kernel density estimators and spline interpolations.

When both detailed age-specific mortality patterns and the underlying mortality time trends are of interest, it is necessary to model mortality surfaces, for example a bivariate distribution. We propose using the penalized composite link model (PCLM) to estimate detailed age-at-death distributions of adjacent calendar years from coarsely grouped death counts. The PCLM is a flexible tool that requires modest assumptions. The coarsely grouped counts observed by the researcher can be regarded as indirect observations of a latent sequence of expected counts. These expected counts represent the distribution on a fine resolution that we seek to estimate. This distribution is assumed to be smooth, that is expected counts on the detailed grid are similar to each other and can be estimated from the composite observed data by maximizing a penalized Poisson likelihood. We use a two-dimensional regression analysis, following the P-spline method by Currie et al. 2004, for ungrouping the age-specific distributions from the coarsely grouped data and smoothing across adjacent calendar years. Thus, the suggested methodology combines two approaches: the PCLM for ungrouping in one dimension, and two-dimensional smoothing with P-splines.,

The PCLM is a powerful tool to model aggregated epidemiological and demographic data. In a one-dimensional setting it has been found to model age-specific grouped data, outperforming kernel density estimator and spline interpolation methods in the presence of open-ended intervals especially, and data suffering from digit preferences. In two- or three-dimensional settings it has been applied to aggregated spatial counts and to fertility rates grouped by age, time and birth order. The methodology was also explored in a Bayesian framework. The use of this model is as yet low. Here, we focus on smooth estimations of mortality surfaces from data grouped in age classes, with particular attention to open-ended intervals at the right-hand tail of the age-at-death distributions.

We will first introduce the PCLM for ungrouping coarse age-at-death distributions of adjacent calendar years and, secondly, analyse age-at-death distributions from cancers of all sites in Denmark from 1980 to 2014, using the proposed methodology. We conclude with a discussion and give a demo R-code of the model in the Appendix and in the Supplementary material, with the aim of bridging the statistical methodology with a view to its application in the public health field.

Methods

The statistical method for ungrouping coarse age-at-death distributions of adjacent calendar years

Suppose that one could observe death counts by single year of age for several adjacent calendar years. Let us denote these death counts on a fine grid by the matrix of dimension, where corresponds to the ages in a single year step and to adjacent calendar years. Death counts are assumed to follow a Poisson process, with expected value In practice, we are only able to observe these death counts grouped in coarser age classes. We denote the observed death counts aggregated in coarse age groups for adjacent calendar years by the matrix of dimension, where corresponds to the age classes and to adjacent calendar years. Death counts follow a Poisson distribution, with expected value, where can be seen as expected values composed of the latent unobserved expectations. Hence, is the expected mortality surface on a fine grid, which we aim to estimate from the composite means . This is done by maximizing a penalized Poisson likelihood via maximum likelihood estimation.

For the purpose of the regression, we arrange the data by column vector, that is , and . The PCLM for ungrouping coarse age-at-death distributions of adjacent calendar years is then given by: where is the sequence of detailed expected death counts that we aim to estimate, is a composition matrix, a B-spline basis with being the corresponding coefficient. The matrix is a block matrix that ‘composes’ from and describes how the latent distribution was mixed before generating the data. The matrix is derived from the Kronecker product of the two marginal composition matrices, where is an identity matrix of dimension and of dimension . Elements of are zero, except for those that indicate the elements of the expected counts γ that are summed to get the aggregated expected counts µ for each calendar year . In order to model rates instead of counts, each column of is multiplied by the corresponding exposure to risk, that is person-years at risk of dying, by single year of age (if the exposures are also grouped, one can produce detailed estimates again with the PCLM).

In the PCLM, is overparametrized with respect to . This is solved by two regularizations following the P-splines methodology by Eilers and Marx (1996). First, the set of is restricted by using a B-spline basis of lower dimension than; second, the linear space of the coefficients is restricted by introducing a penalty matrix .

As the first step, we use as a two-dimensional regression matrix a B-spline basis , which reduces the number of parameters to be estimated, since the number of elements in is very large. Following, where is the marginal B-spline basis in the calendar year direction of dimension and is the marginal B-spline basis in the age direction of dimension, are equal to the number of knots of the corresponding marginal B-spline basis, plus the degree of the spline. Throughout the paper cubic splines are used.

As the second step, we assume that the latent distribution is smooth, that is adjacent elements of are similar, and implement this assumption by a roughness penalty on the parameters, with θ being the vector of coefficients associated with the regression matrix of length . Then, the penalty matrix follows as , with and the smoothing parameters in the age and calendar year directions, respectively, that control the amount of smoothness, and and the matrices that compute the d differences of the regression coefficients . In the paper, the roughness of the coefficients is measured by second-order differences following Currie et al. (2004).

Estimation procedure

The PCLM for ungrouping coarse age-at-death distributions of adjacent calendar years can be estimated by maximum likelihood estimation, by maximizing the following penalized Poisson log likelihood:

Maximizing (2.2) leads to a system of equation that can be solved by an appropriately modified version of the iteratively reweighted least squares (IRWLS) algorithm. The system of equations becomes: where which are diagonal matrices of weights, with tilde indicating the approximation to the solution of a specific iteration, and the penalty matrix. For a given value of the smoothing parameters, the system can be solved. To obtain the optimal values of smoothing, we computed the Akaike Information Criterion (AIC) for a two-dimensional grid of -values (on a log scale) and chose the optimal smoothing corresponding to its minimal value. Standard errors for the estimated can be obtained from a sandwich-estimator for , given by or via a Bayesian approach, such as . Standard errors are obtained by taking the square root of the diagonal elements as , and confidence intervals for the estimated are consequently given by . Standard errors do not reflect the uncertainty introduced by the choice of the smoothing parameters and that are treated as fixed. The estimating procedure was implemented in R version 3.2.2, and the code is provided in the Appendix to this paper.

Data

Age-specific cancer deaths in adjacent calendar years

Age-specific cancer mortality of all sites, including non-melanoma skin cancer, men and women together, in Denmark from 1980 to 2014, is analysed with the proposed method. Raw death counts were retrieved from the Danish Cancer Society from the database for NORDCAN. NORDCAN is a database that contains incidence, mortality, prevalence and survival statistics for 50 major cancers in Nordic countries. Danish data for cancer mortality in the NORDCAN database stem from the Danish Registry of Causes of Death. An internet application provides access to summary data with graphic and tabulation facilities. From the Danish Cancer Society, overall cancer death counts were obtained on a finer resolution, that is from age 0 up to the last age of recorded events in the time window 1980–2014 for each consecutive calendar year and for both sexes combined. Deaths were classified according to the International Classification of Diseases, Tenth Revision (ICD-10, codes CXX.X + D09.0–1+ D30.1–9 + D35.2–4 + D41.1–9 + D32-33 + D42-43 + D44.3–5 + D46-47). The quality of this classification relies mainly upon the reporting of physicians; the analyses of overall cancer mortality reduces the impact of possible misclassifications or discontinuity in registrations of specific cancer sites., The total number of cancer deaths among the Danish population in the 35 years considered was 529 511. The year 1995 had the highest number of cancer deaths, with 15 914 cases, whereas the age of 77 was the most affected age across the entire time window, with 17 675 deaths. The age-at-death distributions for adjacent years from the raw data are illustrated in Figure 1.

Figure 1.

Age-at-death distributions from cancers of all sites, including non-melanoma skin cancer, in Denmark from 1980 to 2014, men and women together from raw data. Source: Danish Cancer Society.

Person-years at risk of dying, by single year of age for each calendar year, were taken from the Human Mortality Database, which is a free database containing data on death from all causes, exposure to risk, and population size, by age and calendar year for 38 countries. Danish data in the Human Mortality Database stem from Statistics Denmark. We used exposures to risk provided for Denmark for each calendar year at 1-year age intervals from 0 to 109 years, and then a final group of 110+.

Application and results

The high-quality data described in the section above, allowed us to investigate the performance of the proposed method, by comparing the raw detailed death rates with the PCLM ungrouped estimates. To construct the mortality surface on a fine grid (Figure 2, right-hand panel), death counts for each age and calendar year, obtained from the Danish Cancer Society, were divided by the corresponding exposures to risk retrieved from the Human Mortality Database.

Figure 2.

Age-specific death rates from cancers of all sites, including non-melanoma skin cancer, in Denmark from 1980 to 2014. Left-hand panel: smooth cancer mortality surface estimated by the PCLM for ungrouping age-at-death distributions of adjacent calendar years from coarsely aggregated data. Right-hand panel: detailed cancer mortality surface from raw data. Sources: Danish Cancer Society and Human Mortality Database.

To study the performance of the proposed model, we grouped the raw death counts from cancer into age classes of different lengths, that is 5-year age classes with an open-ended interval starting at age 85, following the age grouping of most of the freely available databases., After grouping the raw cancer deaths, we applied the PCLM to model the mortality surface at a fine resolution, that is from age zero up to age 110+, for the time window 1980–2014, in 1-year age groups. Age 110+ was set as the maximum age, assuming that no cancer deaths are expected after that age. In the present application, the age groups are equal to , whereas are the detailed ages and the calendar years analysed. of dimension denotes the raw cancer deaths with expected value and of dimension denotes the aggregated cancer death counts with expected value . We aimed at estimating from , with . The two-dimensional regression basis has a dimension of , as it results from the Kronecker product of the two marginal B-spline bases: of dimension , where 35 is the number of calendar years analysed, and of dimension , where 111 is the number of detailed ages. The number of columns of the marginal B-spline bases is the sum of the number of knots chosen and the degree of the spline. Following the rule of thumb described in Currie et al. (2004), one knot every eight elements of equally spaced and cubic splines are used, giving seven columns, (floor(35/8)+3) for and 16 columns, , (floor(111/8)+3) for . The resulting model has parameters ) in order to estimate the mortality surface, that is . The composition matrix of dimension is the Kronecker product of the two marginal matrices and of dimensions and respectively. is an identity matrix with as many rows and columns as the calendar years analysed, whereas has as many rows as the age classes and as many columns as the ages on the fine resolution. When modelling mortality rates, the composition matrix is further multiplied by the respective exposures as offset so that the estimated is a mortality surface. From the Human Mortality Database, detailed exposure figures, that is for single year of age from zero up to 110+ for all calendar years analysed, were retrieved. In order to find the solution of (2.3), the optimal values of the smoothing parameters are selected via AIC: in the age direction and in the calendar year direction. The so obtained smooth cancer mortality surface, of dimension , is illustrated in Figure 2, left-hand panel.

Comparing the raw data by single year of age with the PCLM-ungrouped distributions, we find that the model fit is good. This is captured visually in Figure 2 and it is further quantified by the root mean squared error (RMSE) reported in Table 1.

Table 1.

Root mean squared error (RMSE) calculated in the year and age dimensions. RMSE measures the differences between estimates and raw data. Values close to 0 indicate a good fit. The RMSE of the PCLM is compared with the RMSE derived by only smoothing the raw detailed data. Results do not differ much, indicating that the ungrouping procedure fits the data well

	RMSE PCLM	RMSE smoothing
Year	0.0151	0.0088
Age	0.0034	0.0013

When modelling a smooth mortality surface, the estimated trends in the age and year directions are of particular interest. To better show such patterns modelled by the PCLM, we report estimates versus raw data in a cross-sectional cut in the age and year direction in Figures 3 and 4, respectively. The level of uncertainty of the estimates is captured by the width of the confidence intervals (Figure 3; and Figure A1 in the Appendix): uncertainty increases at high ages, particularly after 100 years of age, where death counts and exposures are of very small numbers.

Figure 3.

Death rates from all cancers, including melanoma skin cancer, in Denmark from 1980 to 2014 for ages 40, 50, 60, 70, 80 and 90 years. Smooth lines are PCLM smooth estimates with confidence intervals (dashed lines), and dots are raw death rates. Sources: Danish Cancer Society and Human Mortality Database.

Figure 4.

Age-specific death rates from all cancers, including melanoma skin cancer, in Denmark from 1980 to 2014. Different calendar years are displayed with different shades of grey: lighter grey corresponds to earlier years and progressively darker grey to more recent years. Left-hand panel: PCLM smooth estimates. Right-hand panel: raw death rates. Sources: Danish Cancer Society and Human Mortality Database.

Figure A1.

Age-specific death rates from all cancers including melanoma skin cancer in Denmark in 1980, 1990, 2000 and 2010. Smooth lines are PCLM smooth estimates with confidence intervals (dashed lines) from death counts grouped in age classes of 5 years, with open-ended age group 85+. Dots are raw death rates. Sources: Danish Cancer Society and Human Mortality Database.

The main trends captured by PCLM are: (i) a decrease in cancer mortality rates after the 1990s for ages 50–75 (Figure 3; and for more details, see Figure A2 in the Appendix); (ii) a decrease in cancer mortality for later cohorts, particularly for young ages and very advanced ages (Figure 4); (iii) a peak of cancer mortality around age 90, followed by a levelling off (Figure 4; and Figure A1). These results are in line with the existing literature on cancer mortality patterns in Denmark, and similar to those observed in international studies. Modelling cancer mortality surfaces with the PCLM allows us to gain information about the mortality patterns, of the oldest old especially, from data that are otherwise aggregated at the higher ages. The levelling off of cancer mortality around age 90 could not be captured with data coarsely grouped in the traditional 85+ open-ended age class. Smoothing in the calendar year direction allows us at the same time to capture the trend across time and to correct for fluctuations or noise in data collection.

Figure A2.

Death rates from all cancers, including melanoma skin cancer, in Denmark from 1980 to 2014 for ages 50–75 years. The mortality trajectories over calendar years are displayed with different shades of grey, according to the ages: lighter grey corresponds to earlier ages, starting at age 50, and progressively darker grey corresponds to older ages, up to age 75. Left-hand panel: PCLM smooth estimates. Right-hand panel: raw death rates. Sources: Danish Cancer Society and Human Mortality Database.

Discussion

We have presented the PCLM to model smooth mortality surfaces from age grouped data. In our application, we were able to reconstruct from aggregated data a cancer mortality surface on a fine age grid equivalent to the raw ungrouped data. At the same time we model a smooth time trend across calendar years.

In the PCLM setting, the user has the freedom of selecting the number of knots in the marginal B-spline bases, the degree of the spline and the order of the penalty plus the smoothing parameters. To choose the optimal number of knots in the B-spline basis, no information criterion is adopted: the idea is to choose a rich enough basis and let the smoothing parameters control the optimal level of smoothness. Following Currie et al. (2004), we suggest using one knot about every eight elements of the latent and detailed distribution, cubic splines and quadratic penalties. This choice is not crucial, since a different degree of the spline or order of the penalty leads to a similar fit. The optimal smoothing parameters are obtained using a grid search of possible values of and selecting those that minimize AIC. In our application, the range of smoothing parameter values is set as. We found that allowing too small values of leads to singularity of the system of equations: This occurs because the parameter is free to vary where it is also not informative.

The main assumptions of the PCLM for ungrouping in two dimensions are smoothness of the estimated surface and Poisson distributed counts. When consecutive observations are expected to be characterized by drastic differences in magnitude, the smoothness assumption is violated. In age-at-death analysis appreciable differences are found, for example between mortality levels at age 0 years and following ages. Explicitly including a point mass for infant mortality can solve this limitation. Another assumption that does not always hold for epidemiological data is the Poisson distribution of counts, that is . In some circumstances data are overdispersed, that is , and therefore the variability of the estimates might be underestimated.

In our application, we chose the same age grouping for all calendar years analysed. However, in practice, the age grouping scheme can differ across calendar years. This is particularly true for historical data or data of the past few decades, aggregated, for example, in 10-year age classes with open-ended intervals starting at age 75 or 80. By changing the composition matrix accordingly, the PCLM can handle different age groupings. In population-based data, deaths are usually reported by each single calendar year. However, in some cases, time might also be grouped, for example in classes of 5 calendar years. Modifying the marginal composition matrix accordingly, the PCLM can be extended to model smooth mortality surfaces from both grouped ages and calendar years.

The proposed methodology can be applied to deaths from various causes—since no assumption about the shape of the sequence γ is made. Because of its flexibility, it can also serve other epidemiological applications with bivariate data that can be assumed to be Poisson distributed, for exampple grouped data on individuals’ body mass index (BMI) and height.

Conclusion

The PCLM has proven useful in this study in gaining detailed information when vital statistics are aggregated in coarse age groups. The PCLM applied to age-specific mortality data allows us to capture the underlying trends in two dimensions, that is age and calendar time, when limited information is to hand.

Vital statistics, such as disease incidence and cause-of-death data, are often provided in coarse age groups, as is the case for the elderly aged 85+. Such aggregated data hide the age-specific trajectories at higher ages and consequently hinder accurate data analysis in an ageing population. Moreover, ungrouping is needed so that data that follow different groupings can be comparable.

The two-dimensional penalized composite link model for ungrouping has a double goal: to estimate detailed age-specific trajectories (of mortality or incidence of various diseases) from grouped data and to model a smooth trend over calendar time to control for small fluctuations or noise and misreporting in the data collection.

Application of this approach in epidemiology and public health research can improve the detail of vital statistics used for further analysis.

Conflict of interest: None declared.

Click here for additional data file.