Multiple imputation of systematically missing data on gait speed in the Swedish National Study on Aging and Care

Study population and study variables

This simulation study is grounded in data from SNAC, which is an ongoing longitudinal cohort study based on samples of the Swedish elderly population launched in 2001. A detailed description of the study structure and methods can be found elsewhere [11]. In brief, it consists of four study sites, namely Kungsholmen, Skåne, Nordanstig, and Blekinge. Data collection includes information on health determinants, disease outcomes, functional capacity, and social conditions [12]. A key predictor of all-cause mortality – gait speed (≤0.8, 0.8-1.2, >1.2 meters per second) – is systematically missing at one study site, Blekinge. To evaluate potential imputation methods to impute the systematically missing variable we defined the outcome as all-cause mortality (yes/no) within 5 years from the examination date. Further, four key health indicators were chosen based on previous analysis [7] that included, severe disability (measured as the number of personal activities of daily living (ADL) a person was unable to perform independently, categorised into 0 and ≥1), mild disability (measured as the number of instrumental activities of daily living (IADL), categorised into 0 and ≥1), cognitive status measured with the Mini-Mental State Examination (MMSE) ranging from 30 (best possible score) to 0, categorised into ≤20, 20-25 and >25), and the number of chronic diseases (count of chronic diseases performed by a clinical examination (ICD-10 diagnostic criteria); categorised into ≤2, 2-4, >4 chronic diseases). Two demographic factors were also included: sex (female/male), and age (59-70, 70-80, 80-90, 90+ years). The cut-offs were chosen based on a combination of avoiding numerical problems during the simulations (i.e., having a sufficient number of subjects in each category of the variables) and on clinical cut-offs used in Santoni et al. [7].

The following sections outline the structure of the simulation study, including a) the data generating mechanism, b) imputation methods, c) analytical methods, and d) estimands and performance measures [13].

Data generating mechanism

In this paragraph we describe how we simulated one IPDMA. We created four synthetic data [14] sets based on the original data from the four SNAC sites (Kungsholmen, Skåne, Nordanstig, Blekinge), keeping their main statistical properties. The synthetic data sets included the five above mentioned predictors of 5-years mortality (gait speed, ADL, IADL, MMSE, number of comorbidities) and two demographic factors (sex and age).

Our simulation strategy consisted of first reproducing marginal and conditional relationships of all the predictors separately within each study; and second, generating individual binary outcomes according to a common set of regression coefficients across all four studies.

All predictors of mortality were randomly generated from a multivariate normal distribution given a set of observed means and variance/covariances. The variables were discretised using the inverse cumulative distribution function method based on observed frequencies [15, 16]. To simulate gait speed in the study with missing information (Blekinge), we used the inverse cumulative distribution function method based on an arithmetic average of observed frequencies available in the other studies (Kungsholmen, Skåne, Nordanstig).

We denote with i the index for the studies included in the prospective MA data. In our simulated scenario, the index i ranges from 1 to 4 representing the Kungsholmen, Skåne, Nordanstig, and Blekinge studies, respectively. Data on seven possibly correlated predictors X (gait speed, ADL, IADL, MMSE, comorbidities, age, and sex) of 5-year mortality risk were randomly generated from single multivariate normal distribution [15]

$$X_i~N({\mu }_i,{\Sigma }_i)$$

where Σ_i is the symmetric observed variance/covariance matrix and μ_i is the observed vector of means for the i-th study. The variance covariance matrices are based on the real data with the same number of discretised variables (eAppendix A in the Supplement). Next, each variable in X_i was discretised using the inverse of the normal cumulative distribution given the empirical mean μ_i, its standard deviation (square root of the diagonal elements of Σ_i), and the observed probabilities shown in Table 1 of the original data [15, 16].

The individual binary outcome, 5-year mortality status, was randomly generated according to a Bernoulli distribution. The outcome probability varied conditionally on all the predictors modelled with indicator variables and a common set of regression coefficients. Independent observations within each study were described by the statistical model Y_i|X_i ~ Bernoulli(π_i) and ${\pi }_i=e^{X_i\beta }\text{/}(1+e^{X_i\beta })$ with X_iβ, the linear predictor of the logit (log odds) of the 5-years mortality probability:

$$\begin{array}{c}{X}_i\beta =-0.786-0.599\text{I}(\text{gaitspeed}=2) \\ -1.237\text{ I}(\text{gaitspeed}=3)\\ -0.791\text{ I}(\text{MMSE}=2)-1.13\text{ I}(\text{MMSE}=3)\\ +0.297(\text{ADL})+0.520(\text{IADL})\\ -0.0197\text{ I}(\text{Cormobidities}=2)\\ +0.297\text{ I}(\text{Comorbidities}=3)\\ -0.725(\text{Sex})+0.949\text{I}(\text{Age}=2)\\ +1.382\text{ I}(\text{Age}=3)\\ +2.17\text{ I}(\text{Age}=4)\end{array}$$

The values of the regression coefficients above were obtained by computing the inverse-variance weighted average of the regression coefficients estimated in the three studies with complete data. The common adjusted effects of gait speed on 5-year mortality risk comparing 0.8-1.2 vs ≤0.8 m/s was $O{R}_1=e^{{\beta }_1}=e^{-0.599}=0.55$ and comparing >1.2 vs ≤0.8 m/s was $O{R}_2=e^{{\beta }_2}=e^{-1.237}=0.29$. These parameter values served as a benchmark to evaluate the performance of the different imputation strategies. Once the individual mortality outcomes were generated based on the above realistic values of the regression coefficients and data, we set gait speed to systematically missing in the study SNAC-Blekinge. Each of the three imputation strategies described in the following paragraph were then used for the same IPDMA to impute gait speed.

Imputation methods

This paragraph describes the three evaluated imputation methods. Based on an IPDMA of four studies, we evaluated two standard imputation methods of systematically missing discrete data, fully conditional specification (FCS) and multivariate normal (MVN). In addition, we evaluate a method based on conditional quantiles (CQI). We imputed the systematically missing variable 100 times for each imputation method given that the Monte Carlo Error (MCE) for 100 imputations and 1 000 repetitions is expected to be very small [17].

Fully conditional specification (FCS)

FCS was first described in detail by van Buuren, Boshuizen and Knook [18]. It identifies a suitable conditional imputation model for each incomplete variable and iteratively imputes until convergence [10, 17–19]. We used FCS for one systematic missing discrete variable. As an imputation model for the discrete missing variable gait speed, we used multinomial logistic regression on joined (appended) datasets [20]. The model included all of the previously mentioned predictors: ADL, IADL, MMSE, number of comorbidities, sex, age, mortality outcome, and study-level indicator variables to identify the original structure of the data.

Multivariate normal (MVN)

MVN imputation assumes that variables being imputed follow a multivariate normal distribution [21]. The method uses an iterative Markov Chain Monte Carlo method to impute missing values [22]. The performance of MVN has also been investigated and evaluated for binary variables [23–25]. The MVN imputation model for gait speed included the same predictors as mentioned for the FCS.

Conditional quantile imputation

This paragraph provides a brief description of a conditional quantile imputation (CQI) approach [26]. A more detailed description of the imputation method can be found in Bottai et al. [27]. Based in Bottai et al. [27] the rationale of the imputation method consists of three main steps: I) quantification of the association between the missing variable, gait speed in our example, and any other observed variable in the three studies (Kungsholmen, Skåne, Nordanstig) where the missing variable is available; II) prediction of the probabilities of any level of the discrete missing variable conditionally on the observed variables in the study with missing data (Blekinge) based on the estimated average relationships obtained in the previous step; and III) imputation of individual missing values of the discrete variable by inverting the cumulative distribution function of a random uniform with quantiles equal to the cumulative predicted conditional probabilities. More technical details of the steps involved in CQI are noted in eAppendix B in the Supplement of this paper.

Analytical methods

The results of the imputation methods were analysed with the following analytical model. MI methods that heavily rely on random-effect models in the case of limited number of studies are difficult to estimate [8]. Moreover, the harmonisation in measurement and data collection between the four study sites of SNAC allow us to assume a common effect as opposed to a heterogenous (random) effect. Thus, the multivariable logistic model described above to predict 5-years mortality risk was estimated for any simulated IPDMA using a two-stage common effect meta-analysis. Estimates of the regression coefficients were combined across imputations using Rubin’s rules [28].

Simulation estimands and performance measure

We simulated the mechanism described above 1 000 times to obtain a sampling distribution of the adjusted effect of gait speed, the predictor that is systematically missing, on mortality risk. Gait speed was modelled with two indicator variables. The performance of the three imputation methods was assessed for the two corresponding regression coefficients (the conditional log odds ratios) of gait speed. The key numerical quantity used to assess the performance was the average relative bias comparing the estimated regression coefficients $({\widehat{\beta }}_1,{\widehat{\beta }}_2)$ with the parameter values in the outcome model previously specified (β₁ = −0.599 and β₂ = −1.237). In addition, we estimated the following performance measures including their Monte Carlo Error (MCE) described in Morris et al. [13] for all three methods: i) bias in point estimate, ii) model-based standard error (SE) (the mean of the SE from the 1 000 repetitions), iii) empirical SE (the standard deviation of the 1 000 estimates from the 1 000 repetitions), and iv) nominal coverage level (proportion of CIs covering the reference value).

Application

We applied the three imputation methods to the original SNAC data. In addition to the adjusted odds ratio and 95% confidence interval, we calculated the predictive capacity of the model based on the area under the curve (AUC) [29]. The imputation and analytical strategy followed the same procedure as in the simulations.

Availability of data and materials

The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request.

Simulation results

Table 2 describes the combined adjusted estimates of the levels of gait speed on 5-years mortality and performance measures for the three imputation methods based on 1 000 simulations. We used the first level of gait speed as a reference (≤0.8 m/s). The average combined estimate for the mortality odds ratio and the relative bias (%) for the second level of gait speed (0.8-1.2 vs ≤0.8 m/s) were highest for the FCS method ${\widehat{OR}}_1=0.579$; relative bias = 8.2%), and lowest for the CQI method (${\widehat{OR}}_1=0.55$; relative bias = 0.7%). For the third level of gait speed (>1.2 vs ≤0.8 m/s), estimates were highest for the MVN method (${\widehat{OR}}_2=0.33$; relative bias = 9.9%). Compared to FCS and MVN, CQI seems less efficient due what can be seen in a higher empirical SE (0.096 compared to 0.074 and 0.079 for MVN and FCS, respectively. The fraction of simulated studies in which the parameter values were included in the confidence intervals were 95.80%, 96.70%, and 94.90% for FCS, MVN, and CQI, respectively.

For the third level of gait speed, the FCS method had a lower relative bias and more precise point estimate compared to its performance for the second level of gait speed (${\widehat{OR}}_2=0.30$; relative bias = 2.5%). Again, the lowest bias was shown for the CQI method with an average combined estimate for the mortality odds ratio of ${\widehat{OR}}_2=0.29$ and a relative bias of 0.9%. The average estimated standard error of the combined estimates was similar across all three methods. A slightly higher average estimated spread could be found for the CQI method for the second and third level of gait speed (SE = 0.096). Similar to the second level of gait speed, CQI indicates a less efficient performance when comparing empirical SE. The nominal coverage seems sufficient for FCS and CQI, however is only at 77.10% for MVN. This is not surprising however, given the large bias in point estimate (0.124).

Figure 1 shows the approximately symmetric and bell-shaped simulated distribution of the combined adjusted effect estimates of the three levels of gait speed on 5-years mortality. For comparing the second vs first level of gait speed (0.8-1.2 vs ≤ 0.8 m/s), all three methods share a similar distribution with a substantial overlap. CQI method indicates a better precision and is centred around the estimated common effect of gait speed on mortality (OR = 0.548). Both the FCS and MVN share a particularly similar distribution with less spread compared to the CQI method. The effect estimates for the MVN method show a larger divergence in distribution from the estimated common effect of gait speed on mortality for both levels of gait speed (0.8-1.2 vs ≤ 0.8 m/s and >1.2 vs ≤ 0.8 m/s).

Application to SNAC data

We applied the investigated MI methods to the four original data sets of the SNAC studies where gait speed was systematically missing at the study site in Blekinge. Table 3 shows that the three MI approaches show comparable effect sizes. As expected from the simulation results, there are no large differences between the three methods on the adjusted effect estimates of gait speed on 5-years mortality. For the second level of gait speed (0.8-1.2 vs ≤ 0.8 m/s) the odds ratios for the FCS, MVN, and CQI method are 0.571, 0.568, and 0.563, respectively. In addition, for the third level of gait speed (>1.2 vs ≤ 0.8 m/s) the odds ratios for the FCS, MVN, and CQI method are 0.309, 0.318, and 0.301, respectively. All three methods share a similar predictive capacity around an AUC of 0.82 with CQI being slightly higher compared to FCS and MVN. There was no substantial difference between the MI methods and the complete case analysis including only the three complete data sets.

Strengths and limitations

We investigated one specific scenario of multiple imputation of a systematically missing discrete variable in an IPDMA with only four studies. The simulations were based on unique data from SNAC and have a practical relevance for researchers involved in working with SNAC and similar data. Further, this is one of the first specific simulations that investigate multiple imputation of systematic missing variables in the context of IPDMA with only a small (<5) number of studies. Last, the imputation method based on condition quantiles operated sufficiently well for small IPDMA and should be further explored.

We acknowledge several limitations. First, the simulations in this paper may have restricted generalisability. Our simulations were specifically tailored to the SNAC studies and might perform differently in other scenarios. We chose to relate the simulations to a specific example to have a direct impact on people working with SNAC data. General conclusions are reasonable to be extended to other IPDMA with only a small number of large observational studies. Second, we only explored one specific scenario as a probabilistic sensitivity analysis of multiple imputation for systematically missing data in IPDMA. Future research in this area should be directed towards a more general approach, including scenarios on varying clusters and sample sizes, varying levels for missing predictors, and combining systematically with sporadically missing data. Last, in our simulation settings we assumed a common effect of gait speed on mortality across the four studies. Although the homogeneity of effects can be easily relaxed, it would be very difficult to derive good estimates of variability across studies based on a limited number of studies.