A Method to adjust for measurement error in multiple exposure variables measured with correlated errors in the absence of an internal validation study

Data and study design

In this work, we use a subset data from a home-based HIV counseling and testing (HBCT) study that was conducted in rural and peri-urban communities in KwaZulu-Natal Province, South Africa, between November 2011 and June 2012²¹. The data were obtained from the Human Sciences Research Council (HSRC) of South Africa²¹. This study was conducted to provide a better understanding of the complexity, severity, and prevalence of non-communicable disease (NCDs) in a community known to have one of the highest rates of HIV incidence and prevalence in the world²¹.

Home-based HIV counseling and testing is a cross-sectional, single-site study in South Africa that aims to increase engagement in HIV care by integrating NCDs screening with community-based HIV testing²². A random sampling approach was used, where 587 participants over the age of 18 were selected from 50,000 people living in the Mpumuza suburb²¹. Anthropometric and biological measures were collected in the survey with the purpose of establishing the prevalence of a range of NCDs and associated risk factors. Eligible individuals participated in a face-to-face interview, physical, psychological and clinical examinations. Persons younger than 18 years living in Mpumuza and all household members not previously enrolled, and members unable to give written consent were excluded from the study. Mobile phones were used for data collection to increase efficiency in data capture and analysis²¹.

In our study, we used a subset data consisting of 76 individuals who self-reported the number of cigarettes smoked, fruit and vegetable consumption. We use the dataset to illustrate the multivariate method in modeling the amount of association between body mass index (BMI) and three exposures (smoking, fruit, and vegetable intakes). BMI was measured in kg/m², while smoking was measured as the average number of cigarettes smoked per day. Initially, fruit and vegetable intakes were measured in terms of the number of servings consumed per day. It is often assumed that a standard portion of fruit/vegetable weighs about 80g⁵. Therefore, for this study, we converted the number of servings to grams per day (g/day) by multiplying the reported number of servings by 80g. The subset data has the following three properties that make it suitable for use in this work: (1) measurement error in the recorded number of cigarettes smoked due to possible misreporting, (2) measurement error in fruits and vegetable consumption due to recall bias, and conversion of the number of servings of fruits and vegetables into grams, and (3) the measurement error in the three exposures is often correlated, for instance, smokers are likely to over-report fruit and vegetable intakes due to their beneficial effects, and to under-report the number of cigarettes they smoke due to the associated harmful effects. Epidemiologically, BMI is used as a risk factor of a health outcome. However, in this study, we model BMI as an outcome as in other several studies, for instance,^23–26. The subset data is only used to illustrate the method and not to draw inference.

Ethical statement

Ethics approval was granted by both HSRC Research Ethics Committee (REC: 1/26/05/11) and the University of Washington Institutional Review Board (48733). Informed written consent was obtained from each participant in the study. Participants were provided with written information on the study (including the study’s background and objectives) and their rights regarding participation and withdrawal at any time.

A measurement error model for the data

An interest in epidemiological study could be to investigate the association between BMI and three exposures namely: fruit, vegetable and smoking using the multiple linear regression

where Y denotes the BMI, β₀ is the intercept, β_X1, β_X2 and β_X3 are the coefficient parameters for the true long-term fruit (X₁), vegetable (X₂) and cigarette (X₃) intake respectively and ϵ is the random error term. In this study, we use vegetable intake and cigarette smoking as confounders and assume that the main interest is in estimating β_X1. In practice, the true intakes are unobservable and, therefore, the intakes recorded in self-reported questionnaires are used. Let W₁, W₂ and W₃ denote the measured versions of X₁, X₂ and X₃, respectively. The use of W_p’s in place of X_p’s, (p = 1, 2, 3), in Equation (1) yields biased estimates ${\widehat{\beta }}_{W_1}$, ${\widehat{\beta }}_{W_2}$ and ${\widehat{\beta }}_{W_3}$ of β_X1, β_X2 and β_X3 respectively. Let ${\widehat{\mathit{\beta }}}_{\mathit{W}}={({\widehat{\beta }}_{W_1},{\widehat{\beta }}_{W_2},{\widehat{\beta }}_{W_3})}^{\text{T}}.$

We assumed that the observed exposures are related to the true exposures with additive measurement error as

where ϵ_W = (ϵ_W1, ϵ_W2, ϵ_W3)^⊤, ϵ_W ∼ N(0, Σ_ϵW); W = (W₁, W₂, W₃)^⊤; α₀ = (α₀₁, α₀₂, α₀₃)^⊤, α₁ = (α₁₁, α₁₂, α₁₃)^⊤; with the terms in α₀ and α₁ quantifying the constant bias and the proportional scaling bias respectively; ϵ_W is a random error term, ϵ_Wi is assumed to be independent of the true exposure X_i and the systematic bias components, α_0i and α_1i.

Bias adjustment methods

A univariate method. In a univariate case, the bias in the association between an outcome and an exposure is adjusted by dividing the unadjusted association estimate by the attenuation factor¹⁹. Attenuation factor (λ_i) is defined as λ_i = var(X_i)/var(W_i), i.e., the ratio of the variance of the true exposure to the variance of the observed exposure, also referred to as reliability ratio. This method ignores correlations between the errors and also the correlation between the true exposures.

Multivariate method. We propose and describe a general approach for handling p-exposures (p≥2) measured with correlated errors. For simplicity and without loss of generality, we assume that W_i is measured without systematic bias (i.e., α_0i = 0, α_1i = 1 in Equation 2). For multiple exposures measured with correlated errors, the adjusted association estimates can be obtained by pre-multiplying the unadjusted association estimates by the inverse of the transpose of attenuation-contamination matrix as

where ${\widehat{\mathit{\beta }}}_{\mathit{X}}^{\mathit{\ast }}$ and ${\widehat{\mathit{\beta }}}_{\mathit{W}}^{\mathit{\ast }}$ denotes vectors of true and biased coefficients for the p-exposures respectively and Λ_p denotes a p × p attenuation-contamination matrix^19,27. The off-diagonal elements in Λ are known as contamination factors while the diagonal elements are called attenuation factors¹⁴. Noteworthy, the attenuation factor quantifies the bias in the association between an outcome and an exposure. In contrast, the contamination factor quantifies the effect of measurement error in one exposure variable on the other exposure variable’s estimate. ${\widehat{\mathit{\beta }}}_{\mathit{W}}^{\mathit{\ast }}$ in Equation (3) can be obtained from the observed questionnaire data.

In the multiple exposures case, the estimate of attenuation-contamination matrix ${\widehat{\text{Λ}}}_p$ is defined as

where ${\widehat{\text{Σ}}}_X^{*}$ is the estimate of covariance matrix of the true exposures, ${\widehat{\text{Σ}}}_W^{*-1}$ is the inverse of the estimate of covariance matrix of the measured exposures, ${\widehat{\sigma }}_{X_i}^2$ is the variance estimate of X_i (i = 1, 2, ..., p) ; ${\widehat{\sigma }}_{X_i{X}_j}$ (j = 1, 2, ..., p; i ≠ j) denotes the covariance estimate between the true exposures; ${\sigma }_{W_i}^2$ is the variance estimate of W_i; ${\widehat{\sigma }}_{W_i{W}_j}$ (i ≠ j) is the covariance estimate between the observed exposures.

The elements of the variance-covariance matrix of the observed exposures, ${\Sigma }_W^{\ast }$, are estimated from the observed data. The variances of the true exposures, ${\sigma }_{X_i}^2$ ’s, can be estimated using validity coefficients for the questionnaire. According to Kipnis et al.⁶, the validity coefficient is given by:

where W_i is assumed to be the measured with error term only and ϵ_Wi is assumed to be independent of X_i. From Equation (5), we estimate the variance of the true exposures as

by incorporating external validation information on ρ_WiXi. To obtain covariances between the true exposures, one of the following two approaches is used: (i) if external information about the correlation between true exposures (i.e. ${\widehat{\rho }}_{X_i{X}_j}$ ) is available, we obtain covariances between true exposures as follows:

where ${\widehat{\sigma }}_{X_i}$ are obtained as shown in Equation (6); (ii) if we can obtain prior information about the correlation between the errors in the observed exposures, ${\widehat{\rho }}_{{ϵ}_{W_i}{ϵ}_{W_j}},$ we can solve for ${\widehat{\sigma }}_{X_i{X}_j}$ by decomposing the covariance of observed exposures into unknown covariance between true exposures and unknown covariance between errors as follows:

where X_i and ϵ_Wj, X_j and ϵ_Wi are assumed to be uncorrelated.

From Equation (2) and Equation (6), the estimate of the error variance ${\widehat{\sigma }}_{{ϵ}_{W_i}}^2$ is

See Appendix B of the extended data²⁸ for the proof.

From Equation (8)–Equation (9), the covariances between the true exposures are given by

Using the observed data and external information, we can determine all the terms required to estimate the attenuation-contamination matrix, Λ, as shown in Equation (4) and adjust for the bias in the association between the exposures measured with error and the outcome using Equation (3).

Illustration of the multivariate method using the study data

We illustrate a method that accounts for uncertainty in the validity measures attributable to heterogeneity in the study populations and in parameter estimation. The proposed Bayesian method applies Markov Chain Monte Carlo (MCMC) estimation approach to combine observed self-reported data and external validation data in adjusting for measurement error in three exposures measured with correlated errors. MCMC is a class of algorithms that samples from the posterior distributions by traversing the parameter space²⁹. The posterior distribution is obtained by updating the prior distribution with observed data. The steps for implementing the trivariate method are described below.

We first obtained external information on validity coefficients and generated validity coefficients for use by interpreting the lower and upper limits obtained from the literature as the 95% credible intervals (CIs) of the distribution of possible values respectively. Due to the skewed distribution of validity coefficients, Fisher’s transformation was used to generate the validity coefficients as explained in the next section.

Second, for the observed exposures, we estimated the posterior distribution of the covariance matrix (Σ_W). The exposures were assumed to follow a multivariate normal distribution with mean and covariance, i.e., W ∼ N₃(µ_W, Σ_W). We assumed a weakly informative multivariate normal prior for µ_W as µ_{W prior} ∼ N₃(0,10⁶ I₃), where I₃ is a 3 × 3 identity matrix. In a multivariate normal distribution, Σ_W must satisfy two conditions: (1) be positive definite (i.e. W^TΣ_WW > 0, for all W) and (2) be a symmetric matrix. The semi-conjugate prior distribution for Σ_W, which has these two properties, is the inverse-Wishart distribution²⁹. To minimize the influence of the prior information on the estimate of Σ_W, we considered weakly informative inverse-wishart prior as Σ_{W prior} ∼ IW(I₃, v), where v = 3 is the degrees of freedom.

Third, using the validity coefficients generated from the external data and the posterior distribution of covariance matrix for observed exposures, we estimated the variance of true intakes, ${\widehat{\sigma }}_{X_i}^2$ (i = 1, 2, 3), using the relationship given in Equation (6) so that

The covariances between true intakes (${\widehat{\sigma }}_{X_i{X}_j}$ ; j = 1, 2, 3) were estimated as,

by incorporating external validation information on correlation between the errors (ρ_{ϵWi ϵWj}). We generated the correlation between errors from a plausible range guided by correlation in the observed data and prior expert information on the most likely sign of the correlation between the exposures, as described in the next section.

Having obtained the covariance matrices of the true and observed exposures, we estimated the attenuation-contamination matrix (Λ₃) from their joint distribution as

where ${\widehat{\text{Σ}}}_X$ is the estimate of covariance matrix of the three true exposures, ${\widehat{\text{Σ}}}_W^{-1}$ is the inverse of the estimate of covariance matrix of the three measured-with-error exposures, ${\widehat{\sigma }}_{X_i}^2$ is the variance estimate of X_i (i = 1, 2, 3); ${\widehat{\sigma }}_{X_i{X}_j}$ (i ≠ j) denotes the covariance estimate between the true exposures; ${\widehat{\sigma }}_{W_i}^2$ is the variance estimate of W_i; ${\widehat{\sigma }}_{W_i{W}_j}$ (i ≠ j) is the covariance estimate between the observed exposures.

Lastly, we fitted a Bayesian multiple linear regression model (hereafter, naive method) to obtain the posterior distributions of the unadjusted coefficient estimates ${\widehat{\mathit{\beta }}}_{\mathit{W}}={({\widehat{\beta }}_{W_1},{\widehat{\beta }}_{W_2},{\widehat{\beta }}_{W_3})}^{\text{T}}.$ In the naive model, we assumed weakly informative normal independent priors by choosing a very small precision (large variance) for the unadjusted coefficient estimates as β_{Wi prior} ∼ N(0, 10⁶). The adjusted coefficient estimates ${\widehat{\mathit{\beta }}}_{\mathit{X}}$ were then obtained from the joint posterior distribution of ${\widehat{\text{Λ}}}_3$ and ${\widehat{\mathit{\beta }}}_{\mathit{W}}$ as

Software implementation of the trivariate method

We implemented the trivariate method in R version 3.6.3 using rjags (version 4-10), coda (version 0.19-3), MCMCpack (version 1.4-9), and mvtnorm (version1.1-1) packages. To facilitate Bayesian estimation of the covariance matrix of the observed exposures (Σ_W), rjags package was used to provide an interface from R to the JAGS library³⁰. JAGS is a gibbs sampler that uses MCMC to draw dependent samples from the posterior distribution of the parameters³¹. The Bayesian estimation of Σ_W proceeded in the following steps: (1) defining a model for Σ_W under Bayesian inference using gibbs sampling (BUGS) algorithm in a stand alone file, (2) reading the model file using the jags.model function, (3) updating the model using the update method for jags objects and (4) extracting the posterior samples of the model using the coda.samples function from the coda package.

MCMCregress function from the MCMCpack package was used to generate a posterior density sample from the naive linear regression model³². MCMC convergence diagnostics of all the model parameters was done using trace plots and autocorrelation (ACF) plots from the coda package³³. See extended data: Appendix C²⁸ for convergence diagnostics results. For each model, the burn-in iterations were set to 2,000 and 10,000 MCMC iterations were run after the burn-in iterations. Every first sample value was kept in the MCMC simulations by using a thinning interval of 1. When compiling a JAGS model, an initial sampling step may be needed during which the samplers learn their behaviour to maximize their performance³⁴. Therefore, the number of iterations for adaptation in the the jags model was set to 500. The results were presented in terms of density plots, posterior mean and median. We compared the results obtained under naive, univariate, and trivariate methods. The R code used for analysis is presented in the extended data²⁸.

External information on the validity coefficient and error correlations for the study data

External information on the validity coefficient and error correlations for fruit, vegetable, and cigarette information was obtained from the literature. According to Kaaks et al.¹, the validity coefficient of self-reported fruit intake ranges from 0.33 to 0.79, while that of vegetable intake ranged from 0.30 to 0.60. A meta-analysis study on the validity of questionnaires assessing fruit and vegetable consumption by Collese et al.² reported validity coefficients of 0.26 for vegetables and 0.49 for fruits. Other similar validation studies reported validity coefficients in the aforementioned ranges for fruits and vegetables^3,4,35. Therefore, based on these information we considered a range of 0.3 to 0.8 for fruits and a range of 0.25 and 0.7 for vegetables.

In the Scottish Heart Health Study of 2,849 men and 2,900 women³⁶, the correlation between the self-reported number of cigarettes and biochemical measures was reported between 0.67 and 0.72. In a study on the validation of self-reported smoking by analysis of hair for nicotine and cotinine³⁷, the validity coefficient between the number of cigarettes smoked per day and nicotine/cotinine levels in hair and plasma was found to be between 0.48 and 0.63, while the correlation between the average number of cigarettes smoked and carboxyhemoglobin was 0.70. In a follow-up study to examine the relationships among self-reported cigarette consumption, exhaled carbon monoxide, and urinary cotinine/creatinine ratio in pregnant women³⁸, a validity coefficient in the range of 0.61 to 0.70 was reported. A study by Stram et al.³⁹ found the correlation between the self-reported number of cigarettes smoked and the true lung dose to be between 0.40 and 0.70, and this range was consistent with the findings from the previously discussed related validation studies. Based on this information, we considered a validity coefficient range of 0.40 and 0.70.

We generated the correlation between errors from plausible ranges that were determined based on the correlation in the observed data and the most probable sign of the correlation among fruits, vegetables, and cigarettes as explained below:

a. Since the correlation coefficient between fruit and vegetable intake in the observed data was positive, we also assumed the error correlation between fruit and vegetables to be mostly positive;

b. An investigation on the correlation coefficient between cigarette smoking and fruits/vegetable intake in the observed data showed a negative correlation coefficient. Based on this and the fact that persons who tend to overstate fruit and vegetable consumption are likely to understate the number of cigarettes smoked, we assumed the error correlation to be mostly negative.

We obtained the upper limits of error correlations by assuming that the error covariance equals the covariance in the observed data and set the lower limit of the error correlation to zero, based on the assumption that the covariance in the observed data equals the covariance between the true intakes¹⁴.

Estimating the distribution of ρ_WiXi

Using the range of plausible values obtained from external validation information, we generated the validity coefficients using the Fisher-Z transformation method by assuming that the reported lower and upper limits are 0.05 and 0.95 quantiles of the uncertainty distribution, respectively. Fisher Z-transformation is a commonly used method to transform the sampling distribution of correlation coefficients to become approximately normally distributed^40,41. The procedure is as outlined below:

Sensitivity analysis

We investigated how varying the level of uncertainty assumed for the limits of the validity coefficients reported from literature affected the estimates for fruit, vegetable, and the average number of cigarettes smoked. We also investigated how the estimates varied with the magnitude of the correlation between errors in fruit and vegetable intake, fruit and cigarette smoking, and vegetable and cigarette smoking. This helps determine the estimates’ sensitivity to various magnitudes of CI and the correlation between errors when using the multivariate method.

Source data

Data used in this study are made available to the researcher upon registration and agreeing to the terms and conditions of use in the HSRC web site at http://curation.hsrc.ac.za/ Dataset-565-datafiles.phtml.

Extended data

Figshare: A Method to Adjust for Measurement Error in Multiple Exposures Measured with Correlated Error in the Absence of Internal Validation Study-Supplementary materials. https://doi.org/10.6084/m9.figshare.13147970.v2²⁸

The file shows the validity coefficient derivation, Proof for the estimate of error variance, R code for implementing the methods and convergence diagnostics results (i.e. Trace plots and ACF plots for the standard deviation and naive regression coefficient estimates of the fruits, vegetables and average number of cigarettes smoked, with explanation) and the sensitivity analysis results (supporting Tables) for varying the magnitude of error correlation between the exposures.

The extended data are available under the terms of the Creative Commons Zero (CC0) license.