Estimating interactions and subgroup‐specific treatment effects in meta‐analysis without aggregation bias: A within‐trial framework

Abstract Estimation of within‐trial interactions in meta‐analysis is crucial for reliable assessment of how treatment effects vary across participant subgroups. However, current methods have various limitations. Patients, clinicians and policy‐makers need reliable estimates of treatment effects within specific covariate subgroups, on relative and absolute scales, in order to target treatments appropriately—which estimation of an interaction effect does not in itself provide. Also, the focus has been on covariates with only two subgroups, and may exclude relevant data if only a single subgroup is reported. Therefore, in this article we further develop the “within‐trial” framework by providing practical methods to (1) estimate within‐trial interactions across two or more subgroups; (2) estimate subgroup‐specific (“floating”) treatment effects that are compatible with the within‐trial interactions and make maximum use of available data; and (3) clearly present this data using novel implementation of forest plots. We described the steps involved and apply the methods to two examples taken from previously published meta‐analyses, and demonstrate a straightforward implementation in Stata based upon existing code for multivariate meta‐analysis. We discuss how the within‐trial framework and plots can be utilised with aggregate (or “published”) source data, as well as with individual participant data, to effectively demonstrate how treatment effects differ across participant subgroups.


| INTRODUCTION
Estimation of within-trial interactions in meta-analysis is considered crucial for reliable assessment of how treatment effects vary across participant subgroups. 1 In recent years, a strong recommendation has emerged for a focus on covariate interactions derived within trials. 2,3 However, the within-trial approach as previously described 2 has various limitations. Where covariate interactions are identified, clinicians and policy makers ultimately need reliable estimates of effect specific to each patient subgroup-ideally expressed on an absolute scale-to allow treatments to be targeted at those who might benefit most and to inform patient choice. Furthermore, the existing method is unable to incorporate information from trials in which interactions cannot be estimated, such as those with all participants belonging to a single subgroup. Older approaches, whereby meta-analyses are carried out separately for each subgroup and then compared, allow the straightforward estimation of such quantities and their presentation on forest plots. However, such approaches conflate within-and across-trial information and are therefore prone to aggregation (or "ecological") bias. 1,4,5 In addition, while the described procedure for two subgroups was straightforward, 2 it was not extended to the case of more than two subgroups, or for subgroups without a clear ordering or choice of reference group. Such cases require the estimation of more than one parameter, necessitating multivariate meta-analysis. Importantly, this imposed unwanted restrictions on the range of patient covariates, potentially of clinical importance, that could be investigated in an appropriate way.
In this article, we extend our within-trial interaction methodology 2 to address these challenges. We intend these methods to be applied with observed effect sizes that are either extracted or calculated from aggregated source data, 6,7 although they may also be derived directly from individual participant data (IPD) as the first stage of a two-stage analysis. 8,9 We provide a fully flexible framework to: (1) estimate within-trial interactions across two or more subgroups, ordered or unordered, for categorical covariates; (2) estimate subgroup-specific treatment effects (hereafter described as "floating" estimates) that make maximum use of available data and are compatible with the within-trial interactions in the sense which we describe below; and (3) clearly present this data using novel implementations of forest plots.
The article is structured as follows: in Section 2 we describe motivating examples for this work, from two previously published meta-analyses. We then proceed in Section 3 to describe the within-trial framework and present the steps involved in estimating within-trial interactions and compatible floating subgroup-specific treatment effects. In Section 4 we apply these methods to data from our motivating examples. A discussion of the capabilities, further considerations and future developments of the framework is given in Section 5, and we finish with a brief conclusion in Section 6.

| MOTIVATING EXAMPLES
Our motivating examples both come from published meta-analyses in which detailed aggregate data were available, thus enabling us to apply the within-trial framework.

| Effects of interleukin-6 antagonists on 28-day all-cause mortality in hospitalised COVID-19 patients
Our first example is a collaborative, prospective metaanalysis of the effects of interleukin-6 antagonists, including tocilizumab, on outcomes for patients hospitalised with COVID-19. 10 One aim was to investigate whether there was a difference in the effect of tocilizumab on 28-day all-cause mortality in subgroups of patients who were and were not assigned to concomitant corticosteroids at baseline. A total of 15 trials were included: 1 that did not use corticosteroids, 2 that gave corticosteroids to all patients and 12 that administered corticosteroids to varying proportions of patients.

| Effects of postoperative radiotherapy on survival of patients with non-small cell lung cancer
Our second example is an IPD meta-analysis which showed that post-operative radiotherapy (PORT) was associated with poorer survival of patients with nonsmall cell lung cancer. In the original analyses of these data, 11,12 using the older approach of pooling separately within subgroups, there was strong evidence that the negative effect of PORT increased with the number of affected lymph nodes. However, in the most recent update, the data were re-analysed using a within-trial approach, which showed much weaker evidence for such a trend. 13 This meta-analysis included 11 trials that provided data on patients' lymph node status, categorised as either N 0 , N 1 or N 2 /N 3 . Only three of the trials recruited patients across all three lymph node categories, and four recruited patients only within one lymph node category. Although this was an IPD meta-analysis, hazard ratios and confidence intervals were reported for each nodal status subgroup for each trial, allowing subgroup-specific effect estimates and standard errors to be extracted.

| Model setup
Suppose we have n studies, within which patient observations are split into k disjoint subgroups. Let b β ji represent the treatment effect estimated within subgroup j of trial i, and b β i be the vector of estimated subgroup-specific treatment effects from trial i, with covariance matrix S i . We suppose, at least for now, that the full set of subgroupspecific treatment effects are observed for every trial, so there is no missing data. Then the standard multivariate meta-analysis model is: where β is the vector of subgroup-specific treatment effects and Σ β is the between-trial heterogeneity covariance matrix associated with β, to be estimated from the data. Structures for Σ β are discussed in Section 3.2 below. Under this model alone, the standard estimator and variance gives results consistent with the older method of pooling within each subgroup separately. Now, the interaction between treatment and the covariate from trial i may be represented by the k À 1 treatment effect contrasts b γ i with respect to a reference subgroup, with covariance matrix V i . These effect vectors, which we refer to as "within-trial interaction" vectors, are derived from b β i and S i via a simple linear combination, represented by the "contrast matrix" M. We may write another multivariate meta-analysis model, similar to Equation (1), for the contrasts b γ i : Typically, covariances between the b β ji at the trial level will be unreported but may be assumed to be negligible, so that S i is diagonal and we may write, where J kÀ1 ð Þ is the k À 1 by k À 1 matrix full of ones: However, unless explicitly stated, nothing in what follows is affected by the precise form of V i ; the equations above demonstrate that non-diagonal matrices S i are easily incorporated. Under our proposed within-trial framework, we wish for our subgroup-specific estimates b β to be compatible with the pooled within-trial interactions b γ, in the sense that the contrasts of b β are equal to b γ. Let b θ represent the subgroup-specific treatment effect estimate in the reference subgroup, which without loss of generality we choose as j = 1. Then b β is parameterised by b θ plus the k À 1 estimated contrasts b γ with respect to that reference, under the relationships shown in Equation (4). In this scenario, we describe b β as the "floating" subgroupspecific effect estimates.
If, as stated above, the reference subgroup is identified by ; that is, a column of length k À 1 full of 'À1' adjacent to the identity matrix of order k À 1. Similarly, we have the design matrix Z ¼ 0 kÀ1

| Random effects considerations
We propose three basic forms for the heterogeneity covariance matrices Σ β and Σ γ , described below in order of increasing complexity:

| Common-effect
A fully common-effect model-that is, with the interaction term(s) and the subgroup-specific treatment effects common across studies-is obtained simply by setting both Σ γ and Σ β to zero.

| Exchangeable random-effects
An exchangeable structure implies that heterogeneity variances and covariances do not depend on which subgroups are being compared. Given the constraint , it can be shown that: That is, we have a single heterogeneity parameter τ 2 γ shared across all estimated interaction contrasts; since contrasts are correlated due to their dependence on the reference subgroup, the correlation between each pair must be one-half. 14,15 Similarly, we have a shared subgroup heterogeneity variance, τ 2 β . Note that Equation (5) includes the special case where τ 2 γ is set to zero, implying Σ β ¼ τ 2 β J, so that the contrasts γ are estimated under a common-effect model with a shared heterogeneity variance τ 2 β both across and within subgroups.

| Unstructured random-effects
Finally, we could allow Σ γ and Σ β to be unstructured, subject to the constraint ). This formulation allows a different heterogeneity variance to be estimated within each subgroup.

| Estimation
In order to derive floating subgroup-specific treatment effects, we first estimate the vector of pooled contrasts γ using a within-trial approach. By subtracting this from the set of trial-level subgroup effects b β i , we proceed to estimate the treatment effect θ in the reference subgroup. Importantly, this ensures that all available data is used, and that the trial-specific weighting is identical across subgroups (see Appendix 1.2). Finally, we associate the floating subgroup-specific treatment effect in the reference subgroup β 1 with b θ, and derive the remaining floating effects as shown in Equation (4). Hence, the floating subgroup-specific effects are constrained to differ by the values of the pooled contrasts γ. We present this approach below, broken down into three steps which are straightforward to implement. In general, the facility to estimate multivariate metaanalysis models is required, such as "mvmeta" in Stata 14,16 or R 17 ; sample Stata code to implement the framework under a common-effect model is presented in Appendix 2, and a full Stata package "metafloat" is available via GitHub (https://github.com/UCL/metafloat). However, with a binary subgroup (k = 2) and under a common-effect model, the approach simplifies considerably (see Appendix 1.5) such that the floating subgroup estimates and variances may be evaluated using spreadsheet software, or even by hand.

| Step 1. Estimate the within-trial interaction(s)
Our first step is to apply the multivariate model for the interactions (Equation 2). In the general case, the standard estimator and variance under this model is: This is a straightforward generalisation of the existing concept of pooling within-trial covariate interactions. 2 The heterogeneity covariance matrix Σ γ may be estimated using restricted maximum likelihood (REML) either with an exchangeable structure or unstructured, as described in Section 3.2. A global test of interaction may be performed via a simultaneous Wald test of all elements of b γ being equal to zero. With three or more subgroups in a natural ordering, it may be useful to perform an additional test for linear trend across the elements of b γ. This may be achieved with a simple modification to Equations (2) and (6); for further details see Appendix 1.4.

| Step 2. Estimate the floating subgroup-specific treatment effects
We now subtract the fitted treatment effect contrasts b γ from the observed subgroup-specific treatment effects b β 2i … b β ki in the non-reference subgroup(s) in each trial i.
Together with the effects b β 1i in the reference subgroup, we then have a set of values which we may interpret as estimates of the true underlying treatment effect θ in the reference subgroup. The standard estimator under the model described by Equations (1) and (2) As stated in Section 3.1, the identification of the reference subgroup by j = 1 is arbitrary; see Appendix 1.3 for suggested modifications to the presented equations for a choice of reference other than the first subgroup. Again, the matrix Σ β may be estimated using REML; but now constraints are required on the values of certain matrix elements, dependent upon values taken from Σ γ (see Section 3.2; further details are given in Appendix 1.7). Having estimated b θ, we form the vector of floating subgroup estimates b β by adding to b θ the appropriate elements of b γ; see Equation (4). Note that at no point are different subgroup estimates compared between studies, so avoiding the risk of introducing aggregation bias. 2 3.3.3 | Step 3. Correct the variances of the floating subgroup-specific treatment effects Although we have now obtained the point estimates for the floating effects, we have a little more work to do to obtain their variance. Observe that not only is b β a random function of both b θ and b γ, but that b θ is itself a random function of b γ. Therefore, although Step 1 gives us Var b γ ð Þ, the apparent variance given by Step 2 only gives us Var b θ γ ð Þ evaluated at the fixed value γ ¼ b γ. To resolve this, note from Equation (7) that b θ may be separated into two parts, only one of which is dependent on b γ. Therefore, with reference to Equation (4), we can express the vector of floating effects b β as a linear combination of the b β i and b γ. As these are independent (see Appendix 1.1), the variance of b β works out to take the following form: where the k by k À 1 matrix A is in terms of fixed, known quantities (see Appendix 1.1 for further details).

| Accommodating studies where not all subgroups are observed
If any estimates b β ji are unobserved for trial i, then b β i and b γ i contain estimates for just the observed subgroups. Similarly, the right-hand side of Equations (1) and (2) contain the corresponding subvectors of β and γ, and submatrices of S i þ Σ β and of V i þ Σ γ . This does not cause any technical issues for estimation. Rather, the unobserved estimates may be considered to be very imprecisely estimated, for example by assigning to them a value of zero for the effect size and a variance much exceeding the largest observed variance (e.g., in our examples below, a value of 10,000 was used). A similar approach has been suggested in related contexts 14 ; we echo their recommendation to check that alternative values of the assigned variance give near-identical results.

| APPLICATION TO EXAMPLE DATA
We now demonstrate this methodology in practice using our two examples described in Section 2. The first example will primarily demonstrate the simplicity of calculations involved in the case of a binary covariate, whilst the second example demonstrates how to proceed in the case of covariate with three levels.

| Effects of interleukin-6 antagonists on 28-day all-cause mortality in hospitalised COVID-19 patients
For this example, we use indices 1 and 2 to represent the "no corticosteroids" and "corticosteroids" groups respectively. The protocol [PROSPERO Identifier: CRD42021230155] specified a priori the use of a common-effect model, which is used here in the first instance. As suggested in Section 2.2, to assist in calculation we assign a value of 0 to the unobserved subgroups for the effect size and 10,000 for the variance. Box 1 contains a detailed walk-through of the within-trial framework under a common-effect model. Stata code to replicate this example is provided in Appendix 2.
In Figure 1, we demonstrate presentation of the effect estimates for each subgroup within each trial on the left-hand side, accompanied by the corresponding pooled floating subgroup-specific treatment effects as derived above. The within-trial interactions for each trial are shown on the righthand side of the plot, with the pooled interaction odds ratio (OR) of 0.69 (95% CI 0.52-0.91, p = 0.01), suggesting a greater benefit from Tocilizumab for patients that did receive corticosteroids compared to those who did not. Our floating subgroup-specific treatment effect estimates for β 1 and β 2 translate to OR 1.10 (95% CI 0.87-1.39) and OR 0.76 (95% CI 0.67-0.86) respectively; note that their ratio exactly equals the pooled interaction estimate of 0.69. In this case, the results are similar to additional published results, which used the older method of pooling within each subgroup separately (Figure 2 of WHO REACT Group 10 ), so here study conclusions remain unchanged. Supplementary random-effects modelling showed no evidence of statistical heterogeneity in either interactions or subgroups (τ 2 β and τ 2 γ both <0.0001), and so results are not presented here. 4.2 | Effects of postoperative radiotherapy on survival of patients with non-small cell lung cancer As described in Section 2.2, data from this IPD meta-analysis had already been analysed using a within-trial approach, assuming a linear trend across the three subgroups to estimate a single interaction. 13 We now extend this example by estimating separate within-trial interactions that compare effects both between participants with N 0 (no affected lymph nodes) and N 1 (one lymph node affected) disease to the reference group N 2/3 (2 to 3 lymph nodes affected), and estimate floating subgroup-specific treatment effects compatible with these interactions. We use the most up-to-date aggregate data, presented in the 2016 Cochrane review, 13 rather than re-analysing the IPD. We choose the N 2/3 subgroup as the reference, as it is the best-represented category among trials including more than one subgroup (see Appendix 1.3 for details of implementation when the reference is not the first subgroup). The floating subgroup-specific treatment effects are presented in Figure 2, below the results for subgroups within each trial (as in Figure 1), and with within-trial F I G U R E 2 Post-operative radiotherapy (PORT) for non-small lung cancer meta-analysis: Effect of PORT treatment on overall survival by nodal status. Green interactions represent the interaction between N 0 and N 2/3 , light blue interactions represent the interaction between N 1 and N 2/3 . Refer to Figure 1 for detailed descriptions on the various markers, except instead of odds ratios, these markers refer in the same way to hazard ratios. Note also that there are two pooled interaction hazard ratios represented by dark blue circles. interactions with respect to the reference N 2/3 subgroup on the right-hand side. Whereas the analysis using the older method of pooling within each subgroup separately suggested that PORT adversely affects the survival of participants with N 0 disease, our floating estimate for N 0 is near to null, similar to the N 2/3 subgroup (see Supplementary Figure 1). Hence, where an across-and-within-trial approach strongly suggests a trend (χ 2 = 4.51 on 1 d.f., p = 0.034), there is no such suggestion under our within-trial framework (χ 2 = 0.31 on 1 d.f, p = 0.58). This difference is likely due to the Lille 18 and Belgium 19 trials, both of which included only N 0 patients and reported large detrimental treatment effects. Under the older method, these trials heavily and directly influence the "naïve" pooled subgroup estimates; whereas under the within-trial framework their influence is more indirect, as they do not contribute to the pooled interactions.
Although common-effect modelling was specified in the PORT protocol (available on request), we fit random-effects models here as a demonstrative example (see Figure 3). Fitting an exchangeable structure results in a shared subgroup heterogeneity variance of τ 2 ¼ 0:062, such that the floating estimates all have moderately wider (compared to common-effect) confidence intervals, whilst their point estimates remain roughly the same. However, fitting an unstructured model shows that there is substantially more heterogeneity within the N 0 subgroup (τ 2 ¼ 0:795) than in the other two subgroups (τ 2 ¼ 0:017 and τ 2 ¼ 0:033 for N 1 and N 2/3 respectively). Hence, the N 0 subgroup is now estimated noticeably differently, with an extremely wide confidence interval and with a point estimate shifted to the left (as a result of extreme estimates from small trials which are now given relatively greater weight). As a result, the floating estimates for N 1 and N 2/3 are also shifted moderately to the left, although in this particular example this does not affect the overall interpretation.

| DISCUSSION
In this article, we extend the previously recommended 2 within-trial approach to a simple yet flexible framework for estimating treatment-covariate interactions in meta-analysis. Both ordered and unordered categorical covariates may be analysed, and heterogeneity covariance is handled straightforwardly. Crucially, we also provide practical methods of estimating floating treatment effects for each level of a covariate.
Previous criticism of testing for interactions using a within-trial approach has included the potential for loss of power relative to older methods, for example from exclusion of "single-subgroup" studies. By contrast, the within-trial framework discussed here allows the entirety of information from all trials to be incorporated, whilst also taking account of the magnitude, direction and precision of the within-trial pooled interaction. The floating subgroup-specific treatment effects from this framework are compatible with within-trial interactions-that is, each contrast is equal in magnitude to F I G U R E 3 Post-operative radiotherapy (PORT) for non-small lung cancer meta-analysis: Effect of PORT treatment on overall survival by nodal status. Results from a common-effect model, and from random-effects models with exchangeable and unstructured heterogeneity covariances Floating subgroup-specific treatment effects (left panel) and pooled within-trial interactions (right panel). The top dark blue open circle represents a meta-analysis of the interaction hazard ratios between N 0 and N 2/3 , with the bottom dark blue open circle representing a meta-analysis of the interaction hazard ratios between N 1 and N 2/3 . Note only summary information is presented on this plot and not trial-level information as in Figures 1 and 2. [Colour figure can be viewed at wileyonlinelibrary.com] that of the corresponding interaction-and each is estimated using the same set of trial-specific weights which reflect the totality of information in each trial. By comparison, pooling separately within subgroups assigns a different set of weights within each subgroup, making interpretation of subgroup contrasts difficult. Note that, when the subgroup of interest is fully observed and "balanced" in all trials (in the sense that S i ¼ s i I 8i), then subgroup-specific treatment effects from our within-trial framework and from the "pooling within subgroups" approach can be shown to coincide under a common-effect model (see Appendix 1.6). Hence, we propose our framework as giving a more appropriate representation of comparative subgroup-specific treatment effects. As such, where clinical decision-making is best served by deriving additional statistics, such as translation (where appropriate) from relative to absolute scale, 20,21 this may be done most effectively via floating subgroup estimates.

| Further considerations
We assume that treatment effects for the different subgroups within each trial can be extracted or calculated from aggregated source data, 6,7 or are derived directly from IPD. The results of a previous methodological review 2 suggested that the availability of suitable aggregate data from peer-reviewed articles was uncommon. However, this can be obtained through working collaboratively with trialists, 22,23 which allows prospective agreement on the important patient-level covariates to analyse, as demonstrated in recent meta-analyses in prostate cancer [24][25][26]27,28 As discussed previously in both our examples, whilst studies including only a single participant subgroup cannot contribute to the within-trial interaction, we can use their information in the within-trial framework when estimating floating subgroup-specific treatment effects compatible with this interaction. Doing so requires that we make the assumption of transitivity across subgroups: that is, that any non-observed subgroup-specific treatment effect could in principle have been observed, and that its true value (or its true distribution under a random-effects model) would be identical to those of the remaining studies. If such studies (either singly, or taken together) are assigned relatively large weights, then this assumption may have a substantial impact upon the floating subgroup estimates. Conversely, if a "single-subgroup" estimate is extreme relative to the remaining data, then it may be questionable whether the pooled interaction (and therefore the transitivity assumption) is applicable to that trial. In the context of our within-trial framework, we therefore strongly recommend that reviewers critically evaluate the design and setting of "single-subgroup" trials to assess whether the transitivity assumption holds, and whether construction of floating subgroup estimates including such trials would involve excessive extrapolation outside the range of observed data. Whilst this is not an issue for the examples presented in this paper, the illustrative example in Fisher et al. 2 and the data presented in eFigure 1 of Sterne et al. 27 are examples of situations where such evaluation might be recommended. One possible statistical approach here might be to include within the model described by Equation (4) an additional parameter, identifying "singlesubgroup" studies-or, more generally, the set of subgroups provided by each trial-analogously to identification of "designs" and inconsistency parameters in network meta-analysis. 29 As a sensitivity analysis, this approach might be used to test the effect of removing such studies from the estimation procedure.
The choice between common-effect and randomeffects modelling is discussed elsewhere. 30,31 In the multivariate case, as here, the most general random-effects model uses an unstructured covariance matrix, estimating a separate heterogeneity for each subgroup and thereby allows conclusions to be drawn regarding the source of heterogeneity or whether results from some subgroups are more reliably estimated than others. However, the number of parameters to be estimated increases sharply with the number of subgroups, and models might not converge with smaller datasets. Conversely, an exchangeable structure estimates a single heterogeneity parameter across all subgroups, with a common pairwise covariance. This simpler structure borrows strength across subgroups for the estimation of τ 2 β . The special case of τ 2 γ ¼ 0 places random-effects on the subgroup estimates but not on the interaction estimates, and is simpler still. Note that we do not recommend random-effects formulations which violate the compatibility relationships given in Equation (4).
A further consideration when utilising the withintrial framework is presentation of the results on a forest plot. Where possible, we advocate presentation of the within-trial interactions alongside the within-trial subgroup-specific treatment effects, as in Figures 1 and 2. However, for covariates with three levels or more and/or with large numbers of included studies (as in Figure 2), such plots become challenging to present clearly. Furthermore, although the choice of reference subgroup does not affect the pooled interactions or floating subgroup estimates themselves, it can have considerable impact on presentation. Specifically, if a trial does not contain the reference subgroup, then interactions cannot be plotted on the second panel of the plot for that trial. Typically, the reference subgroup should contain a relatively large amount of information; but it is possible that a suitable choice does not always exist. In such cases, we would recommend checking the "transitivity assumption" as described above; otherwise, displaying only the summary results (as in Figure 3 and Supplementary Figure 1), or only the left-hand panel, may be a suitable alternative.

| Future developments
The methods presented here may be generalised to encompass two-stage IPD meta-analysis 8 which, among other things, would allow interactions with continuous covariates to be analysed in a within-trial framework. However, our methods are not directly compatible with a one-stage IPD meta-analysis, in which a single generalised linear model is fitted to all trial data simultaneously. It has been noted 32 that simply including a global treatment-covariate interaction term into a one-stage model would incorporate across-trial information and is at risk of aggregation bias, 2,3 and so it has been recommended that within-and across-trial information should be separated out to provide the equivalent one-stage within-trial interaction. 32 Further possibilities to estimate patient-level interactions exist for a one-stage IPD metaanalysis, with shrinkage methods suggested for use when a large number of treatment-covariate interactions are to be modelled simultaneously. 33 Further work is needed to extend the within-trial framework into a one-stage IPD meta-analysis setting in order to estimate subgroupspecific treatment effects that are compatible with patient-level interactions and free from aggregation bias, and also in settings where "personalised" treatment effects utilising multiple covariates are desired.

| CONCLUSIONS
Our within-trial framework allows straightforward estimation of a range of within-trial treatment-covariate interactions, compatible subgroup-specific treatment effects and corresponding forest plots, to clearly and effectively demonstrate how treatment effects differ across patient subgroups.
AUTHOR CONTRIBUTIONS Peter J. Godolphin: methodology; analysis and interpretation of data; writing -original draft; writing -review & editing. Ian R. White: conceptualization; methodology; software; writing -review & editing. Jayne F. Tierney: conceptualization; acquisition of data; writing -review & editing; project administration. David J. Fisher: methodology; software; analysis and interpretation of data; writing -review & editing; supervision.