Methods to estimate underlying blood pressure: The Atherosclerosis Risk in Communities (ARIC) Study

Antihypertensive medications complicate studies of blood pressure (BP) natural history; BP if untreated (“underlying BP”) needs to be estimated. Our objectives were to compare validity of five missing data imputation methods to estimate underlying BP and longitudinal associations of underlying BP and age. We simulated BP treatment in untreated hypertensive participants from Atherosclerosis Risk in Communities (ARIC) in visits 1–5 (1987–2013) using matched treated hypertensive participants. The underlying BP was imputed: #1, set as missing; #2, add 10 mmHg for systolic, 5 mmHg for diastolic; #3, add medication class-specific constant; #4, truncated normal regression; and #5, truncated normal regression including prior visit data. Longitudinal associations were estimated using linear mixed models of imputed underlying BP for simulated treated and measured BP for untreated participants. Method 3 was the best-performing for systolic BP; lowest relative bias (5.3% for intercept at age 50, 0% for age coefficient) and average deviation from expected (0.04 to -1.79). Method 2 performed best for diastolic BP; lowest relative bias (0.6% intercept at age 50, 33.3% age <60, 9.1% age 60+) and average deviation (-1.25 to -1.68). Methods 4 and 5 were comparable or slightly inferior. In conclusion, constant addition methods yielded valid and precise underlying BP and longitudinal associations.

Increasing use of antihypertensive medications is certainly important for preventing complications of hypertension but makes studies of blood pressure (BP) natural history difficult to investigate. For instance, setting BP to missing for all subjects receiving antihypertensive treatment can substantially reduce the sample size and power to detect effects of interest. Alternatively, including treated BP measures while investigating BP natural history can bias results towards subset of non-hypertensive subjects. In order to prevent these biases, several methods have been proposed to impute the "underlying BP", BP that would have been measured if the individual was not taking antihypertensive medications. Adding a constant BP to the measured BP of the treated individual has been most commonly used, especially in genetic epidemiology [8][9][10]. More recent methodology incorporates regression-based imputation [8,11,12]. However most of these studies focus on evaluating in cross-sectionally [8,11], and few have investigated longitudinal settings [12], a time frame more relevant to natural history. Indeed, multiple BP measurements in longitudinal studies may amplify biases of not accounting for treatment. Alternatively, prior BP under no antihypertensive medications in longitudinal settings may be informative to estimate underlying BP in case an individual subsequently starts taking antihypertensive medications.
Therefore, the objective of this study is to comprehensively compare existing methods to impute underlying BP among individuals taking antihypertensive medications in a longitudinal setting using data from a large community-based cohort, the Atherosclerosis Risk in Communities (ARIC) Study, with repeated BP measurements over 26 years of follow-up.

Study population
The ARIC Study recruited 15,792 individuals aged 45-64 years from four communities in Forsyth County, NC; Jackson, MS; northwest suburbs of Minneapolis, MN; and Washington County, MD. To date, ARIC has completed the baseline in 1987-89 and four follow-up visits in 1990-92, 1993-95, 1996-98 and 2011-13, respectively. Detailed study design and methods have been previously described [13,14]. All visits included interviews, anthropometric measurements, and BP measurements under a common protocol.
BP was measured via a random zero sphygmomanometer in visits 1-4 and an automatic sphygmomanometer (OMRON HEM-907 XL) in visit 5 by a certified trained technician in accordance with ARIC manual of procedures [13,15]. Briefly, participants were asked not to smoke, eat, perform physical exertion or experience cold temperatures 30 minutes prior to measurements. Participants were also asked to stay in a sitting position for at least 5 minutes prior to measurements. Each participant contributed three BP measurements per visit for visits 1-3 and 5 and the reported BP was an average of the second and third measurements [13][14][15]. For visit 4, BP was measured twice and reported as the average [13,14].
For this study, participants with measurements of SBP and DBP at 2 or more visits and data on antihypertensive medication use and class of medication were included. Other requisite variables were age, sex, race, and body mass index (BMI). Detailed inclusion criteria for visits 1-5 are presented in S2

Simulation study
We generated hypothetical studies where the goal was to evaluate the natural history of underlying BP but where some patients were treated for hypertension. The hypothetical study population consisted of all non-hypertensive participants (SBP <140 mmHg, DBP <90 mmHg, not taking antihypertensive medications) and untreated hypertensive participants (not taking antihypertensive medications despite either SBP !140 mmHg or DBP !90 mmHg) at each visit. We will refer to the BP among the untreated hypertensive participants as "measured untreated BP". Then, we simulated treatment for the untreated hypertensive participants by assigning a hypothetical treated BP ("simulated treated BP") as the measured BP from a matched treated hypertensive participant (taking antihypertensive medications). Matching was performed using strong predictors of BP, including age ± 5 years, sex, race and BMI ± 5 kg/m 2 . Using antihypertensive medication literature to guide us [6][7][8]10,16,17], we also included in the matching criteria that the simulated treated BP to be within a conservative interval of 20 mmHg for SBP and 10 mmHg for DBP of the measured untreated BP. If more than one treated hypertensive subject satisfied the above matching criteria, then one of the subjects that met the criteria was chosen at random.
We subsequently applied five methods to account for the treatment of the hypertensive patients within our hypothetical study. First, we set the underlying BP value to missing for all participants in the hypothetical study that were taking antihypertensive medications (Method 1). The remaining four methods imputed the underlying BP using the simulated treated BP. Method 2 (simple constant addition) added a constant of 10 mmHg to the simulated treated SBP and 5 mmHg to the simulated treated DBP [8,10,12] and method 3 (medication class-specific constant addition) was to add a constant to the simulated treated SBP and DBP depending on the individuals' medication class. The constants were derived from a summary report of the effects of antihypertensive medication classes [16]. For six representative medication classes [3,16], we applied the following expected treatment effect in mmHg for SBP/DBP as constants: 14.3/10.4 if European American and 6.8/6.6 if African Americans for angiotensin-converting enzyme (ACE) inhibitors, 15.5/11.7 for alpha blockers, 14.8/12.2 for beta blockers, 15.3/10.5 for calcium channel blockers (including dihydropyridine and non-dihydropyridine), 15.5/9.0 for diuretics (including thiazide-like and loop), and 14.8/10.5 for miscellaneous BP-lowering medications [16]. Race specific medication constants were used for ACE inhibitors due to the clear difference by race unlike the other five medication classes (S1 Table) [16]. If the participant was taking combination therapy, the medication with the largest expected monotherapy effect was chosen as the primary medication and the non-primary medications had a percentage of their monotherapy effect listed above [16]. Detailed description of the estimated treatment effects is provided elsewhere [16], and a summary can be found in S1 Table. For example, if a participant was taking diuretics and beta blockers, the total expected effect is 17.1 mmHg SBP/13.1 mmHg DBP; the primary medication effect for diuretics is 15.5 mmHg SBP/ 9.0 mmHg DBP and the non-primary medication effect for beta-blockers is 1.6 mmHg SBP (14.8x11%)/4.1 mmHg DBP (12.2x34%).
Methods 4 and 5 imputed the underlying BP using a regression-based approach. Method 4 (truncated normal regression) used a censored normal distribution to generate the imputed underlying BP by assuming that the distribution of underlying BP is normally distributed and that the underlying BP is higher than the measured BP under treatment [8]. The censored normal model was fit using a priori selected predictors of underlying BP from previous literature [1][2][3] and empirically included if P<0.05. The SBP model included untreated hypertensive participants' age, age 2 , race, sex, BMI, height, height 2 , and interaction terms for age and age 2 with the sex and height variables. The DBP model included the above variables and instead of age 2 , a spline term for age at 60 years and its respective interaction terms with sex, height, and height 2 . Detailed methodology of censored normal regression has been previously published [8]. Briefly, the underlying BP for a treated participant was imputed by randomly sampling from a normal distribution with mean defined as the model predicted mean and standard deviation defined as the model mean squared error where the imputed underlying BP was greater than the simulated treated BP. This random sampling was repeated 15 times and the final imputed underlying BP was an average of the 15 values in order to incorporate variability. Method 5 (truncated normal regression with prior visit BP and antihypertensive treatment) was the same as method 4 but included two extra predictors, measured BP and the use of antihypertensive medications at the prior visit, in the censored normal regression model.
We simulated 1000 hypothetical studies and applied each of the five methods above to each simulated study (S1 Fig). As sensitivity analyses, we modified methods 4 and 5 and used 1 instead of 15 repeated sampling of matching in order to incorporate BP variability. All analysis was performed using R version 3.2.0 [18].

Validity of imputation methods
The validity of the imputation methods was assessed several ways. First, we compared the distribution of the imputed underlying BP based on methods 2 through 5 to the measured untreated BP. The average difference between the measured untreated BP and the imputed underlying BP was calculated for each visit. In addition, we assessed the difference in the empirical distributions for each of the 1000 simulations using the Anderson-Darling test. The utility of the Anderson-Darling test may be limited in conditions with large sample size and deviances from a normal distribution. We interpreted our results within the context of these limitations.
Subsequently, we obtained the true association between underlying BP and age using the measured BP values for the hypothetical study population; all non-hypertensive participants (SBP <140 mmHg, DBP <90 mmHg, not taking antihypertensive medications) and untreated hypertensive participants; by fitting a linear mixed model including fixed effects for age (spline at 60 years for DBP, as aforementioned) and random effects at the patient level for both the intercept (centered at 50 years) and age term(s). The same linear mixed models were fit to the 1000 hypothetical studies applying methods 1 through 5. The relative bias of each method was defined as the relative difference in the average regression coefficient compared to the true regression coefficient and was calculated for the model intercept and age term(s).
Our a priori hypothesis was that method 5 would outperform the other methods because the method was regression-based and incorporated the tracking of prior BP. All analysis was performed using R version 3.2.0 [18].

Matching and simulation
In brief, 10,283 participants were included in the hypothetical study population, of whom at visit 1, 9,122 were non-hypertensive and 1,161 were untreated hypertensive and were assigned simulated treatment. Detailed exclusion criteria (S2 Table) and demographic characteristics of participants at all other visits (S3-S6 Tables) are provided as supporting information. Results are primarily presented below using visit 2 since visit 1 was not informative for imputation methods using prior visit data (method 5) and since results were consistent across the remaining visits with the largest sample size at visit 2.
A total of 1,078 untreated hypertensive participants and 4,720 treated hypertensive participants had complete data at visit 2 (Fig 1). Table 1 shows the demographic and clinical characteristics of 1,078 untreated hypertensive participants compared to the 4,720 treated hypertensive participants at visit 2. The untreated and treated hypertensive participants were of comparable age and BMI. In general, the untreated hypertensive participants were more likely to be male, White and educated at least high school level. Compared to the treated hypertensive participants, comorbid conditions were less prevalent in untreated hypertensive participants, including kidney dysfunction (estimated glomerular filtration rate [eGFR] <60 ml/ min/1.73m 2 ), diabetes, prevalent coronary heart disease (CHD), parental history of CHD, and history of heart failure. All of the untreated hypertensive participants were matched to a treated hypertensive participant according to the criteria previously described (Fig 1). The distribution of simulated treated BP was shifted towards lower BP compared to the distribution of measured untreated BP, as expected ( Table 2).
Imputed underlying BP Table 2 shows the relative accuracy of methods 2 through 5 at imputing the underlying SBP and DBP (Table 2). In general, the best-performing method was the same for each visit for SBP and DBP, respectively. Method 3 (class-specific constant addition) was the most accurate for  2, method 5 followed by methods 4 and 3 had the lowest mean difference between measured untreated and imputed underlying DBP. The distributions of the imputed underlying BP were significantly different than the distributions of the measured untreated BP (Table 2), although the impact of the large sample size on the statistical significance cannot be discounted. While the mean difference between the imputed underlying and underlying BP gives an estimate of the average accuracy, the relative accuracy of the imputed underlying BP differed by the underlying BP. Using one round of SBP imputation at visit 2 as an example (Fig 2), the imputed underlying SBP tended to underestimate the underlying BP at the outer ranges of the measured untreated SBP, especially the upper ranges, resulting in V-shaped Bland Altman plots. This V-shaped pattern was apparent in methods 2-4. The V-shape was partly due to the matching criteria, which enforced a lower and upper bound in the simulated treated BP, and partly due to the applicability of each method's assumptions at the extremes of the simulated treated BP. The precision of the estimation was also comparable across methods although the regression-based methods tended to be more precise as expected (Fig 2). Although methods 2-5 demonstrated wide 95% limits of agreement (34.65 to 51.99 mmHg), method 5 (truncated normal regression with prior visit BP and antihypertensive treatment) yielded the most precise estimation followed by method 4 (truncated normal regression). The pattern observed for this

Longitudinal associations of imputed underlying BP with age
The results of the longitudinal analysis are presented in Table 3 for SBP and DBP. Although estimations of intercept and age are important for validity, we highlighted the age associations given the longitudinal nature of our study. The linear association between underlying SBP and age was underestimated (by 28%) using method 1, where we set the underlying SPB missing for patients who were treated for hypertension. In comparison, methods 2-4 tended to have little to no bias (<6%), with lowest bias for method 3. For DBP, the coefficients measuring the true association between underlying DBP and age were small (0.03 and -0.11 before and after age 60 years, respectively) and thus even small differences in the average age coefficients across the imputation methods translated to larger relative bias compared to the SBP analysis. Nevertheless, method 1 yielded the highest bias (200% and 72.7%, respectively). Method 2 yielded the least biased estimates for both coefficients for age (relative bias 33.3% and 9.1%, respectively). Method 5 also yielded relatively little bias (relative bias 33.3% and 45.5%, respectively).

Discussion
This longitudinal simulation study demonstrated the need for imputation methods to account for the effect of antihypertensive medications, especially in cohorts drawn from the general population where medication use is common. As anticipated, setting the treated BP measurements as missing was especially biased. In addition, constant addition methods (methods 2 and 3) tended to be superior to the other methods in both estimating underlying BP and in the longitudinal analysis of underlying BP and its change over age. Although the truncated normal regression methods (methods 4 and 5) were comparable, unexpectedly they did not outperform methods 2 and 3 overall. Therefore, a hybrid of the method 3 for SBP and method 2 for DBP may be an option for epidemiological studies focusing on both SBP and DBP.
Our findings are consistent with previous cross-sectional literature. Tobin et al. reported that the constant addition is better than or equivalent to the truncated normal regression in cross-sectional setting [8] and surprisingly, we confirmed this in a longitudinal setting. There may be several reasons behind the validity of constant addition methods. Firstly since the definition and management of hypertension are well established, treatment itself and the medication classes may be more relevant to estimating underlying BP than individual characteristics such as age and sex. Second, previous work to account for medication use addressed the issue of underlying BP as a missing data problem [8,11,12,19]. In theory, imputation methods such as truncated normal regression that are based on variables in the dataset assume that the dependent variable (i.e. underlying BP) is missing at random (MAR) [8,11,12,19]. However since underlying BP can be missing not at random (MNAR) or can be due to variables that were unobserved, not measured well or not included in the model, regression based methods may not perform as well as expected [8,11,12,19]. This could also contribute to why ignoring the missingness of underlying BP (method 1) would not be appropriate, as reported previously [8,11,12,19]. For other traits less variable over time such as cholesterol, regression-based methods may outperform constant addition methods [20,21]. This may even be true for BP if measured more often. Nevertheless, constant addition methods are simpler compared to truncated normal regression and therefore advantageous in the context of large-scale epidemiological investigations. Among constant methods, method 3 for SBP and method 2 for DBP performed well. The discrepancy could be explained in part due to how DBP vs. SBP responds to medication. In general SBP is well studied since it is a more established risk factor of future cardiovascular disease and other complications [2,3]. Thus, medications may be investigated for their effect on SBP more than their effect on DBP. Also just as the relationship between age and DBP is more complex [2], the relationship between DBP and medication class may be more difficult to summarize as a series of constants. In addition, the range of DBP is smaller than SBP and thus even slight measurement error in the studies measuring medication effects can be more influential for DBP. It is also important to note the large relative bias estimates for DBP is in part due to the small coefficients in the linear mixed models. Finally, we cannot discount that DBP compared to SBP may be less influenced by the medication class or combination of medications used for treatment. Thus adding 5mmHg may better capture the effects of antihypertensive medication use for DBP. This is the first study to our knowledge to comprehensively assess the existing imputation methods of underlying BP in a longitudinal setting. While our imputation methods relied on previously published methods [8,12,22], we expanded on these methods by incorporating the heterogeneity of medication classes in method 3 and the tracking of BP over time in method 5. Also unlike previous studies, we used data from an existing cohort rather than simulate data based on a known model [8,12,22]. This allowed for a more empiric assessment of the imputation methods. Our empiric dataset was robust for addressing antihypertensive use; it was derived from a large, multi-center cohort with a large proportion of treated and untreated hypertensive participants. The ARIC cohort also included standard and reliable measurement of BP and covariates over nearly three decades.
There were some limitations to our study. Firstly, since the number of matches considerably decreased, we did not take into account comorbid conditions such as kidney dysfunction and diabetes during matching despite these altering the indication and the effects of antihypertensive medications. For instance, the guidelines for BP target levels are different among individuals with and without comorbidities such as diabetes. Among individuals with diabetes, the SBP target levels are <130-140 mmHg and the DBP target levels are <80 mmHg [23]. Nonetheless, the relative performance of the imputation methods was consistent with the main findings even when we included kidney function and diabetes in the matching criteria as a sensitivity analysis (S7 Table). Secondly, our assessment of BP relied on an average of two measurements on a single day in several years, although in clinical practice hypertension would be determined based on BP measurements at multiple occasions over several months. Thirdly, a constant threshold (SBP !140 mmHg or DBP !90 mmHg) was used in all five ARIC visits for the simulation even though ARIC visits occurred from 1987-2013 and thresholds for diagnosing and treating hypertension changed overtime (i.e., 160 mmHg SBP/90 mmHg DBP before 1992-3 and 140/90 mmHg since then) [24,25]. Since a constant cutoff over all visits irrespective of comorbidities was used, bias may have been introduced in classifying participants as untreated hypertensive and in the matching process. Finally, the adherence and dose of antihypertensive medications were not measured. Both factors could lead to heterogeneity in the effect of antihypertensive medications on BP and could affect the observed longitudinal associations.

Conclusions
In summary, our longitudinal simulation study demonstrated that imputation methods can account for treatment in large-scale epidemiological studies of BP natural history. Although both constant addition methods and regression methods are useful, constant addition methods may outperform regression methods in estimating and modeling underlying BP, especially SBP. Further studies are required to assess generalizability and robustness of our findings accounting for possible modifiers of underlying BP such as dose, adherence, and cause of combination therapy.
Supporting information S1 Fig. Summary of simulation study: Example participant. The figure illustrates an example participant ('Participant A') with elevated blood pressure (BP) who does not take antihypertensive medications (untreated hypertensive) and thus was eligible for the simulation study. Participant A's BP (systolic blood pressure 153/diastolic blood pressure 93) was measured without antihypertensive treatment (measured untreated BP). Participant A was matched to a treated hypertensive participant with similar age, sex, race and body mass index (Participant X). Participant X's BP (150/90) was used as the simulated treated BP (BP that would have been measured for Participant A had he/she been taking antihypertensive medications). Using the simulated treated BP, the following imputation methods were conducted: #1, set as missing; #2, single constant addition (systolic BP 10/diastolic BP 5); #3, class-specific constant addition; #4, truncated normal regression; #5, truncated normal regression with prior visit BP and antihypertensive treatment. :935-942. The expected effect of each medication class was estimated and listed overall and by race for ACE inhibitors. If the medication was used as part of combination therapy, the average effects of the non-primary medication were reduced. This weighted effect was expected to be influenced by whether a diuretic was part of the combination therapy. Abbreviations: BP, blood pressure; SBP, systolic blood pressure; DBP, diastolic blood pressure; ACE, angiotensin-converting enzyme. (DOCX) S2 Table. Summary of participants and inclusion criteria at visits 1-5. At each visit, participants were included in the analysis if they did not have missing age, sex, race, body mass index (BMI), antihypertensive medication use and antihypertensive medication class. Participants were further excluded if they had blood pressure (BP) measured at only 1 visit. Of the remaining participants, they were classified according to medication status (treated hypertensive, untreated) and BP levels among untreated (non-hypertensive, untreated hypertensive). Untreated hypertensive participants were included in the simulation study if they had at least 1 match among the treated hypertensive according to age, sex, race and BMI.  Table. Accuracy of imputed underlying BP to measured untreated BP matched on diabetes and kidney function. The relative bias (%) was calculated as the average relative difference in regression coefficient from each simulation (100 x [β imputation method −β measured untreated BP ] / β measured untreated BP ). Controls were matched on diabetes and kidney dysfunction status. Abbreviations: SBP, systolic blood pressure; Var, variance; DBP, diastolic blood pressure.