Email updates

Keep up to date with the latest news and content from BMC Medical Research Methodology and BioMed Central.

Open Access Research article

Attributable risk from distributed lag models

Antonio Gasparrini1* and Michela Leone2

Author Affiliations

1 Department Medical Statistics, London School of Hygiene and Tropical Medicine, Keppel Street, WC1E 7HT, London, UK

2 Department of Epidemiology, Lazio Regional Health Service, Via di Santa Costanza 53, 00198 Rome, Italy

For all author emails, please log on.

BMC Medical Research Methodology 2014, 14:55  doi:10.1186/1471-2288-14-55


The electronic version of this article is the complete one and can be found online at: http://www.biomedcentral.com/1471-2288/14/55


Received:28 January 2014
Accepted:8 April 2014
Published:23 April 2014

© 2014 Gasparrini and Leone; licensee BioMed Central Ltd.

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly credited. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Abstract

Background

Measures of attributable risk are an integral part of epidemiological analyses, particularly when aimed at the planning and evaluation of public health interventions. However, the current definition of such measures does not consider any temporal relationships between exposure and risk. In this contribution, we propose extended definitions of attributable risk within the framework of distributed lag non-linear models, an approach recently proposed for modelling delayed associations in either linear or non-linear exposure-response associations.

Methods

We classify versions of attributable number and fraction expressed using either a forward or backward perspective. The former specifies the future burden due to a given exposure event, while the latter summarizes the current burden due to the set of exposure events experienced in the past. In addition, we illustrate how the components related to sub-ranges of the exposure can be separated.

Results

We apply these methods for estimating the mortality risk attributable to outdoor temperature in two cities, London and Rome, using time series data for the periods 1993–2006 and 1992–2010, respectively. The analysis provides estimates of the overall mortality burden attributable to temperature, and then computes the components attributable to cold and heat and then mild and extreme temperatures.

Conclusions

These extended definitions of attributable risk account for the additional temporal dimension which characterizes exposure-response associations, providing more appropriate attributable measures in the presence of dependencies characterized by potentially complex temporal patterns.

Keywords:
Attributable risk; Attributable fraction; Distributed lag models

Background

Epidemiological studies usually rely on effect summaries based on ratio measures, such as relative risk, odds ratio or rate ratio, with the choice depending on the specific study design [1]. Although these measures are ideal for summarizing the association of interest, they offer limited information on the actual impact of the exposure. This information is critical for the planning and evaluation of public health interventions, and it is better provided by relative excess measures such as the attributable fraction (AF), or absolute excess measures such as the attributable number (AN). Steenland and Armstrong offer a thorough overview on the topic [2]. Here we generally refer to these summaries as attributable risk measures.

Problems in the definition of these measures may arise in the presence of delayed associations, occurring when an exposure generates a risk lasting well beyond the exposure period. Researchers in different fields, during the last thirty years, have proposed approaches for modelling this type of association [3-6]. In time series analysis, a popular approach is based on distributed lag models (DLMs) [7,8], generalized to distributed lag non-linear models (DLNMs) when including non-linear exposure-response associations [9,10]. Recently, the DLNM framework has been extended beyond time series data for modelling such dependencies, defined exposure-lag-response associations, in different study designs [11]. However, attributable risk measures have not been developed for the DLM and DLNM class, with the result that their current definitions do not take into account the additional temporal structure in the exposure-response association. Previous work investigated the issue in time series analysis, although without producing a general approach [12,13].

In this contribution, we extend this research and attempt to formalize the definition of attributable risk measures within the DLNM modelling framework. In particular, we illustrate how complex non-linear and temporal patterns can be accounted for in the computation of the attributable risk. Also, we show how attributable components related to different exposure ranges can be separated. We propose an example on the estimation of the health burden attributable to temperature with time series data, a problem which has received quite a lot of interest recently due to climate change predictions. However, the approach can be easily applied to other exposure-lag-response associations. The method is implemented in simple functions developed within the R software package and provided as additional files.

Methods

Attributable risk measures

A general definition of the attributable fraction AFx and number ANx for a given exposure x can be provided by:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M1','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M1">View MathML</a>

(1a)

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M2','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M2">View MathML</a>

(1b)

with n as the total number of cases. The parameter βx used in Eq. (1a) represents the risk associated with the exposure, and it usually corresponds to the logarithm of a ratio measure such as relative risk, relative rate or odds ratio. It is generally obtained from regression models while adjusting for potential confounders. The general definition of βx used here refers to the association with a specific exposure intensity x compared to a reference value x0. For linear exposure-response relationships, the association can also be reported as β·x, where in this case β refer to a unit increase in x. For binary variables reporting presence/absence of the exposure, Eq. (1a) simplifies to AF=(RR−1)/RR, with RR as relative risk, as reported by Steenland and Armstrong [2]. We keep the more general definition of βx, which is easily applicable to non-linear exposure-response relationships, throughout the manuscript.

The theoretical nature of these effect measures is based on a counterfactual, where the observed condition is compared with a reference state which never occurred. This state postulates that the same population is followed in an identical situation where only the exposure level changes to the reference value x0. Typically, such a reference is represented by the absence of association, meaning x0=0 and <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M3','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M3">View MathML</a>. However, different counterfactual conditions can be used, for example a lower exposure which can be determined by an intervention. In this case the quantity βx can be simply re-parameterized as <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M4">View MathML</a>, and Eq. (1) still applies.

Eq. (1a) can be extended to define the risk attributable to multiple exposures x1,…,xp:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M5','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M5">View MathML</a>

(2)

with <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M6','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M6">View MathML</a> obtained by substituting Eq. (2) in Eq. (1b) [2]. For the specific form of Eq. (2), it should be noted that <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M7','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M7">View MathML</a>, i.e. the sum of the attributable risk measured for individual exposures is usually higher than their concurrent attributable risk.

A review of the DLNM modelling framework

The basic idea underpinning the development of DLNMs is that the risk at time t can be described as the weighted sum of effects cumulated from a series of exposures <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M8','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M8">View MathML</a> experienced in the past over the lag period =0,…,L, with 0 and L corresponding to minimum and maximum lags, respectively. The risk can be described by the function f(x), determining the exposure-response, and the function w(), specifying the lag-response, related to the weights given to exposures at different lags . These functions are combined in a bi-dimensional exposure-lag-response functionf · w(x,). Algebraically, the risk is defined by a function s(x,t;η), written in terms of parameters η as:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M9','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M9">View MathML</a>

(3)

The function s(x,t) is computed as the approximate integral of the exposure-lag-response function over the lag dimension, representing the cumulated risk over the lag period. The parameterization in the final step of Eq. (3) is obtained through a cross-basis, involving a tensor product between the basis chosen for f(x) and w(), generating the transformed variables wx,t linearly combined with the parameters η. Simpler DLMs are defined by Eq. (3) by assuming f(x) as linear. Algebraic details and additional information are provided elsewhere [11]. The cross-basis is specified with a reference value x0 used later as a centering point for the function f(x), which is used to define the counterfactual condition.

The complex parameterization of exposure-lag-response associations provided by Eq. (3) can be more easily interpreted by computing effect summaries from the original parameters η. Specifically, the bi-dimensional exposure-lag-response risk surface modelled through f · w(x,) can be expressed by a grid of effect summaries βx,, each interpreted as the association with an exposure x at lag versus the reference value x0. For a given time t, the cross-basis parameterization in (3) can be re-expressed as:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M10','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M10">View MathML</a>

(4)

This overall cumulative association is composed of the sum of contributions βx, from exposures <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M11','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M11">View MathML</a> experienced within the lag period. Algebraic definitions have been previously provided [11].

Forward and backward perspectives

The term βx, for each intensity x can be interpreted using two complementary perspectives, illustrated graphically in Figure 1. From a forward standpoint, looking from current exposure to future risks, the terms βx, are the contributions from the exposure xt occurring at time t to the risk at times t+0,…,t+L, identified by green circles. From a backward standpoint, looking from current risk to past exposures, the terms βx, are the contributions to the risk at time t from exposures <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M12','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M12">View MathML</a> experienced at t0,…,tL, identified by yellow squares. The underlying curve in Figure 1 depicting such associations is called the lag-response curve related to a given exposure intensity x. The sum of these contributions over the whole lag period can be interpreted as the overall cumulative risk.

thumbnailFigure 1. Conceptual model for the interpretation of exposure-lag-response associations: forward (left panel) and backward (right panel) perspectives.

Attributable risk from DLNMs

The effect summaries provided above can be used for defining attributable risk measures within the DLNM framework. The idea is to treat the associations with exposures at different lags as independent contributions to the risk. A neat definition can be developed using a backward perspective, assuming the risk at time t as attributable to a series of exposure events in the past. The backward attributable fraction b-AFx,t and number b-ANx,t at time t are obtained by substituting Eq. (4) in Eq. (2):

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M13','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M13">View MathML</a>

(5a)

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M14','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M14">View MathML</a>

(5b)

with nt as the number of cases at time t. This structure is consistent with the configuration of the regression model usually applied to fit the data, where the risk at time t is associated with lagged exposures at times t. The definition of backward attributable risk requires an extended version of the counterfactual condition accounting for the additional lag dimension: b-ANx,t and b-AFx,t are interpreted as the number of cases and the related fraction at time t attributable to past exposures to x in the period t0,…,tL, compared to a constant exposure x0 throughout the same period.

An alternative version can be obtained using a forward perspective. Among other possible definitions, forward attributable number f-ANx,t and fraction f-AFx,t can be defined as:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M15','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M15">View MathML</a>

(6a)

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M16','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M16">View MathML</a>

(6b)

This alternative version has some advantages if compared to the backward definition. First, the counterfactual condition is simpler: f-AFx,t and f-ANx,t are interpreted as the fraction and number of future cases in the period t+0,…,t+L attributable to the single exposure x occurring at time t, compared to x0. Moreover, the overall cumulative risk <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M17','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M17">View MathML</a> for a given exposure xt in (6a) is available also when the bi-dimensional exposure-lag-response is reduced to uni-dimensional exposure-response relationship, a step often needed in multi-site studies [14]. In contrast, all the lag-specific contributions are needed to compute <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M18','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M18">View MathML</a> in (5a) for the backward counterpart.

However, the forward version also has an important limitation, related to the fact that the contributions are associated to risks measured at different times. The attributable number f-ANx,t in (6b) is computed by averaging the total counts experienced in the next 0,…,L times, thus only approximating the lag structure of risks. This approximation is likely to produce some bias, which is expected as an underestimation of the attributable number if compared to the backward version.

Separating attributable components

The definitions provided in Eq. (5)–(6) can be extended to separate the attributable components related to specific exposures or exposure ranges. This will be used later in the example to single out the contributions from cold and heat in temperature-health associations. Let’s define a range r=[l,h] between low and high exposure limits l and h. The definition of forward attributable number <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M19','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M19">View MathML</a> and fraction <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M20','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M20">View MathML</a> limited to exposures within the range r is clear-cut, as they are either equal to the quantities reported in Eq. (6) if xr or zero otherwise. Adopting a backward perspective, a similar definition of <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M21','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M21">View MathML</a> can be obtained by modifying Eq. (5a) as:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M22','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M22">View MathML</a>

(7)

simply selecting the risk contributions from past exposures included in the range r. The related attributable number <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M23','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M23">View MathML</a> is computed by substituting Eq. (7) into Eq. (5b). Attributable components referring to different ranges can be summed up, as all are defined using the same counterfactual condition of a constant exposure x=x0 for the whole lag period =0,…,L.

The forward version has the additional advantage that for two non-overlapping ranges r1 and r2 the sum of the components is equal to the overall attributable risk, namely <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M24','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M24">View MathML</a>. In contrast, adopting a backward perspective <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M25','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M25">View MathML</a>, as the risks are simultaneously computed for the same time t in the like of Eq. (2).

Total attributable risk

The attributable risk measures provided above can be computed for each of the i=1,…,m observations in a data set. An estimate of the total attributable number ANtot and fraction AFtot is provided by:

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M26','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M26">View MathML</a>

(8a)

<a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M27','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M27">View MathML</a>

(8b)

The equations above can be applied either to forward or backward attributable risk and to separate components, simply substituting the related attributable numbers in Eq. (8a).

Computing uncertainty intervals

Analytical formulae for confidence intervals of attributable risk measures are not easily produced [15], and this also applies to the extended versions developed here. Although approximated estimators have been proposed [15,16], in this context the most straightforward approach is to rely on interval estimation obtained empirically through Monte Carlo simulations [17,18]. Basically, we take random samples η(j) of the original parameters η of the cross-basis in Eq. (3) from the assumed multivariate normal distribution with point estimate <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M28','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M28">View MathML</a> and (co)variance matrix <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M29','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M29">View MathML</a> derived from the regression model. These samples η(j) are used to compute <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M30','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M30">View MathML</a> for =0,…,L and each intensity x, empirically reconstructing the distributions of the attributable measures defined in Eq. (5)–(8). The related 2.5th and 97.5th percentiles of such distributions are interpreted as 95% empirical confidence intervals (eCI).

Results

The methods illustrated in the previous section are applied to estimate the all-cause mortality risk attributable to temperature, using daily time series from two cities, London and Rome, in the periods 1993-2006 and 1992-2010 respectively. R scripts and data implementing the method and partly replicating the results are provided as Additional files 1, 2, 3, 4, 5 and 6.

Additional file 1. R script implementing the function to compute attributable risk measures.

Format: ZIP Size: 2KB Download fileOpen Data

Additional file 2. Documentation of the attrdl function.

Format: PDF Size: 69KB Download file

This file can be viewed with: Adobe Acrobat ReaderOpen Data

Additional file 3. R script for fitting the regression models in the example.

Format: ZIP Size: 1KB Download fileOpen Data

Additional file 4. R script for computing attributable risk in the example.

Format: ZIP Size: 1KB Download fileOpen Data

Additional file 5. R script for producing the graphs in the example.

Format: ZIP Size: 2KB Download fileOpen Data

Additional file 6. File with data used in the example.

Format: ZIP Size: 271KB Download fileOpen Data

Modelling strategy

We fitted a standard time series Poisson model allowing for overdispersion, controlling for seasonal and long term trends and day of the week, using a 10 df/year spline and indicator variables, respectively. Model selection is still an issue of current research within the DLNM framework, although simulation studies indicate a good performance of methods based on the Akaike Information criterion (AIC) [11]. Considering the illustrative purpose of the example, we selected a-priori the cross-basis function in Eq. (3) for representing the association between mean daily temperature and mortality, basing our choice on previous analyses. Specifically, the cross-basis is composed of a quadratic B-spline with two equally-spaced knots as the exposure-response function f(x), and a natural cubic B-spline with three equally-spaced knots in the log-scale as the lag-response function w() over lags 0–25.

In the specific case of temperature where a null exposure condition cannot be defined, a reasonable choice is to center the cross-basis in Eq. (3) to the temperature of minimum risk, as suggested in previous publications [13]. This optimal temperature corresponds to 20°C and 21°C for London and Rome respectively, and it represents the reference point x0 for the computation of the attributable risk measures. These are obtained for the whole temperature range, and then for cold and heat contributions by separating the associations with temperatures lower or higher than x0. In addition, the attributable components are separated further in mild and extreme cold and heat by selecting as cut-off values the 1st and 99th percentiles of city-specific distributions, corresponding to 0.4°C and 23.7°C in London and 2.6°C and 28.6°C in Rome.

We derived empirical confidence intervals for backward total attributable numbers and fractions, computed overall and for separated components, by simulating 5,000 samples from the assumed distribution of <a onClick="popup('http://www.biomedcentral.com/1471-2288/14/55/mathml/M31','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2288/14/55/mathml/M31">View MathML</a>.

Risk attributable to temperature

The estimated associations between temperature and all-cause mortality in the two cities are illustrated in Figure 2. The left panels show the bi-dimensional exposure-lag-response surfaces, while the right panels display the overall cumulative exposure-response curves, interpreted as the risk cumulated over the entire lag period of 0–25 days. The associations are represented in the RR scale, with the centering point and the cut-off values for defining extreme cold and heat displayed as dotted and dashed vertical lines respectively.

thumbnailFigure 2. Association between temperature and all-cause mortality.Left panels: 3-D graphs of the exposure-lag-response risk surface. Right panels: overall cumulative exposure-response associations with temperature distributions. London 1993–2006 and Rome 1992–2010.

Table 1 reports the estimates of the total backward and forward attributable fraction b-AFtot and f-AFtot, with 95%eCI. The backward version indicates that overall 13.59% and 12.58% deaths are attributable to temperature in London and Rome within the study periods, respectively. As expected, the corresponding estimates computed forward are affected by a certain degree of negative bias associated to the averaging of future mortality within the lag period, as described above. Nonetheless, the difference is not substantial in this case.

Table 1. Total mortality fraction (%) attributable to temperature, computed backward (b-AFtot) and forward (f-AFtot), reported as overall, hot and cold components with 95% empirical confidence intervals (eCI)

Cold, heat and extreme components

The total backward attributable risk is then separated into components due to cold and hot temperatures, defined as those below and above the optimal temperature, respectively. The estimates, computed using Eq. (7), are reported in Table 1. The comparison of the two contributions clearly indicates that cold is responsible for most of the mortality attributable to temperature, with b-AFtot equal to 12.95% and 10.84%, compared to 0.66% and 1.74% for heat, in the two cities. Estimates of forward attributable risk are very similar, and as expected their sum is equal to the overall burden, differently than for the backward version.

The analysis is extended by further separating the attributable components into contributions from mild and extreme temperatures, with results summarized in Table 2. For cold temperatures, the contribution from extreme days accounts only for a minimal part of the mortality burden, while still 12.38%–12.48% and 10.27%–10.37% of deaths in London and Rome, respectively, are attributable to mild cold days. This result is expected when inspecting the right panels of Figure 2, with the definition of mild cold including the majority of days in the series for both cities.

Table 2. Total mortality fraction (%) attributable to temperature, computed backward (b-AFtot) and forward (f-AFtot), reported as components from mild and extreme hot and cold contributions with 95% empirical confidence intervals (eCI)

In contrast, the comparison between the two cities is rather different for the components attributable to mild and extreme hot temperatures. In spite of the stronger risk in London, the attributable fraction is similar for extreme heat and even higher in Rome for mild heat (1.32%–1.45% versus 0.25%–0.33%). This apparent contradiction is explained by the different temperature distribution, and in particular the percentile corresponding to the optimal temperature, corresponding to 93.6th and 72.5th in London and Rome. This result suggests the hypothesis that the population in Rome is more adapted to the range of temperatures corresponding to extreme hot if compared to London, where the population experienced only a few days of unusually high temperatures.

The harvesting paradox

Accounting for the additional lag dimension in exposure-lag-response associations involves further complexities in the interpretation of attributable risk measures. We now focus our attention to the association with hot temperature in Rome. The left panel of Figure 3 shows the estimated lag-response curves at various temperatures, computed versus the reference optimal temperatures of 21°C. The graph indicates a strong risk in the first lags, then followed by a protective association at longer lags. This pattern is consistent with the harvesting hypothesis, which assumes that the initial risk is partly discounted by the depletion of the pool of susceptible after an extreme heat event [19]. Also, the extent of harvesting seems more pronounced for extreme temperatures.

thumbnailFigure 3. Lag structure and harvesting paradox. Left panel: lag-response associations between various temperatures and all-cause mortality. Rome 1992–2010. Right panel: daily number of deaths attributable to heat, computed forward (green circles) and backward (yellow squares), and temperature trend. Rome July-Sept 1992.

This phenomenon has interesting implications. An example is offered by the right panel of Figure 3, illustrating the estimated daily deaths b-ANx,t and f-ANx,t attributable to heat, computed backward and forward for the first summer in the time series for Rome, with related temperature trend. As expected, a substantial number of deaths are attributable to temperatures above the optimal value, represented by the horizontal dotted line, in the period mid-July to mid-August. The trend of forward attributable deaths f-ANx,t closely follows the daily temperatures, consistently with the definition of number of deaths attributable to the temperature in day t cumulated in the next L days. In contrast, the backward attributable number b-ANx,t decreases to zero and even becomes negative in late summer days, although the overall cumulative exposure-response in Figure 2 (bottom-right panel) does not show a RR below 1 for any temperature.

This paradox is explained by the counterfactual condition associated with the backward perspective. Specifically, each b-ANx,t compares the association with the observed temperatures in the past L days to a constant exposure x0. In the presence of harvesting, the observed population becomes ‘healthier’ than the counterfactual population after a series of heat days, due to the depletion of the susceptible pool. This explains the negative attributable numbers for specific combinations of lagged exposures. This fact emphasises that harvesting should not be interpreted as a true protective association at longer lags, but rather as an artefact due to a change in the underlying population following a stress, which affects the counterfactual condition. This issue is relevant when using backward attributable risk measures b-ANx,t and b-AFx,t to assess the contribution of specific days. However, similarly to the net overall cumulative risk, the total attributable number b-ANtot and fraction b-ANtot, produced by Eq. (8) and reported in Tables 1 and 2, account for the discount by summing the contributions over the whole series.

Discussion and conclusions

In this contribution we illustrate an extended definition of attributable risk measures based on the DLNM framework. Consistently with this class of models, such a definition accounts for the complex pattern of potentially non-linear and delayed associations described through exposure-lag-response associations.

Two alternative definitions of attributable risk are proposed, assuming backward or forward perspectives. The former provides more consistent estimators which naturally arise from the structure of the regression model, where distributed lag terms at times t contributes to the risk at time t. The forward attributable measures, in contrast, are affected from a negative bias related to the averaging of future counts, which nonetheless is likely to be relatively low. On the other hand, the forward version is well suited for separating the risk in components attributable to different ranges, as their sum matches the overall risk. Furthermore the forward perspective, looking from current exposure to future risk, seems more appropriate for quantifying the health burden due to specific exposure occurrences, as it is based on a more coherent counterfactual condition. Corrections have been proposed in previous works on risk attributable to multiple exposures [20-22], and can be applied to the backward version.

Strictly speaking, the definition given in Eq. 1a is interpreted as the attributable fraction among the sub-population of exposed subjects. In the setting of time series analysis for environmental stressors, the whole population is usually considered as exposed, and this definition can be more generally interpreted as the population attributable fraction. If only a subset is instead exposed, Eq. (5)–(8) can be easily extended using the equations proposed by Steenland and Armstrong [2] for population attributable risk.

Previous papers suggested approaches for producing attributable risk from distributed lag models when applied to heat-mortality associations. Baccini and colleagues applied DLMs and computed attributable risk measures, specifically addressing the issue of harvesting [12]. Honda and colleagues illustrated an analysis on the mortality burden due to heat using DLNMs [13]. However both approaches are limited, as the former assumes a linear threshold form of the exposure-response, while the latter averages the non-linear risk across the whole temperature range. In this paper we offer a formal and more consistent definition of such attributable risk measures.

An advantage of the proposed method is the provision of estimates for separate components of the attributable risk, associated with different exposure ranges. In the specific case of temperature-health associations, this allows the separation of attributable risks from cold and heat, and further from mild and extreme temperatures. The estimates reported in the example highlights how the simple analysis of exposure-response curves can be misleading in the attribution of risk, and that most of the mortality in the two cities is in fact attributable to mild cold temperatures, in spite of the relatively low RR.

The availability of attributable risk measures, complementary to estimates of exposure-response associations, is essential for the identification and planning of public health interventions. Their extension to exposure-lag-response associations allows the computation of such measures from dependencies showing potentially complex non-linear and temporal patterns.

Abbreviations

AF: Attributable fraction; AN: Attributable number of cases; RR: Relative risk; DLM: Distributed lag models; DLNM: Distributed lag non-linear models; eCI: empirical confidence interval.

Competing interests

The authors declare that they have no competing interests.

Authors’ contributions

AG conceived the idea of attributable risk measures for DLNMs and worked out their algebraic definitions. AG and ML developed the final version of the measures, planned the example and carried out the analysis, drafted the final version of the manuscript. AG provided the software implementation through the R scripts. Both authors read and approved the final manuscript.

Acknowledgements

AG is supported through a Methodology Research fellowship awarded by Medical Research Council-UK (grant ID G1002296).

References

  1. Rothman KJ, Greenland S, Lash TL: Modern Epidemiology . Philadelphia: Lipcott Williams & Wilkins; 2008. OpenURL

  2. Steenland K, Armstrong B: An overview of methods for calculating the burden of disease due to specific risk factors.

    Epidemiology 2006, 17(5):512-519. PubMed Abstract | Publisher Full Text OpenURL

  3. Thomas DC: Statistical methods for analyzing effects of temporal patterns of exposure on cancer risks.

    Scand J Work Environ Health 1983, 9(4):353-366. PubMed Abstract | Publisher Full Text OpenURL

  4. Breslow NL, Day NE:

    Statistical Methods in Cancer Research. Vol. II: The Desing and Analysis of Cohort Studies .

    Lyon: International Agency for Reasearch on Cancer (IARC); 1987:232–271. Chap. 6: Modelling the relationship between risk, dose and time

    OpenURL

  5. Sylvestre MP, Abrahamowicz M: Flexible modeling of the cumulative effects of time-dependent exposures on the hazard.

    Stat Med 2009, 28(27):3437-3453. PubMed Abstract | Publisher Full Text OpenURL

  6. Richardson DB: Latency models for analyses of protracted exposures.

    Epidemiology 2009, 20(3):395-399. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  7. Almon S: The distributed lag between capital appropriations and expenditures.

    Econometrica 1965, 33:178-196. Publisher Full Text OpenURL

  8. Schwartz J: The distributed lag between air pollution and daily deaths.

    Epidemiology 2000, 11(3):320-326. PubMed Abstract | Publisher Full Text OpenURL

  9. Armstrong B: Models for the relationship between ambient temperature and daily mortality.

    Epidemiology 2006, 17(6):624-631. PubMed Abstract | Publisher Full Text OpenURL

  10. Gasparrini A, Armstrong B, Kenward MG: Distributed lag non-linear models.

    Stat Med 2010, 29(21):2224-2234. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  11. Gasparrini A: Modeling exposure-lag-response associations with distributed lag non-linear models.

    Stat Med 2014, 33(5):881-899. PubMed Abstract | Publisher Full Text OpenURL

  12. Baccini M, Kosatsky T, Biggeri A: Impact of summer heat on urban population mortality in Europe during the 1990s: an evaluation of years of life lost adjusted for harvesting.

    PloS One 2013, 8(7):69638. Publisher Full Text OpenURL

  13. Honda Y, Kondo M, McGregor G, Kim H, Guo YL, Hijioka Y, Yoshikawa M, Oka K, Takano S, Hales S, Kovats RS: Heat-related mortality risk model for climate change impact projection.

    Environ Health Prev Med 2013.

    (Epub ahead of print. doi:10.1007/s12199-013-0354-6)

    OpenURL

  14. Gasparrini A, Armstrong B: Reducing and meta-analyzing estimates from distributed lag non-linear models.

    BMC Med Res Methodol 2013, 13(1):1. PubMed Abstract | BioMed Central Full Text | PubMed Central Full Text OpenURL

  15. Graubard BI, Fears TR: Standard errors for attributable risk for simple and complex sample designs.

    Biometrics 2005, 61(3):847-855. PubMed Abstract | Publisher Full Text OpenURL

  16. Cox C, Li X: Model-based estimation of the attributable risk: a loglinear approach.

    Comput Stat Data Anal 2012, 56(12):4180-4189. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  17. Greenland S: Interval estimation by simulation as an alternative to and extension of confidence intervals.

    Int J Epidemiol 2004, 33(6):1389-1397. PubMed Abstract | Publisher Full Text OpenURL

  18. Wood SN: Generalized Additive Models: An Introduction with R . Boca Raton: Chapman & Hall/CRC; 2006. OpenURL

  19. Schwartz J: Is there harvesting in the association of airborne particles with daily deaths and hospital admissions?

    Epidemiology 2001, 12(1):55-61. PubMed Abstract | Publisher Full Text OpenURL

  20. Eide GE, Gefeller O: Sequential and average attributable fractions as aids in the selection of preventive strategies.

    J Clin Epidemiol 1995, 48(5):645-655. PubMed Abstract | Publisher Full Text OpenURL

  21. Land M, Vogel C, Gefeller O: Partitioning methods for multifactorial risk attribution.

    Stat Methods Med Res 2001, 10(3):217-230. PubMed Abstract | Publisher Full Text OpenURL

  22. Hamel JF, Fouquet N, Ha C, Goldberg M, Roquelaure Y: Software for unbiased estimation of attributable risk.

    Epidemiology 2012, 23(4):646-647. PubMed Abstract | Publisher Full Text OpenURL

Pre-publication history

The pre-publication history for this paper can be accessed here:

http://www.biomedcentral.com/1471-2288/14/55/prepub