### Study population

The NAS is an ongoing longitudinal study of aging launched in 1963 by the Veterans Administration [10]. 2,280 community-dwelling healthy men living in the Greater Boston area were enrolled between 1963 and 1968. Every 3–5 years, after an overnight fast and abstention from smoking, participants visited the clinic for an extensive physical examination, laboratory tests, blood collection, and a self-administered questionnaire on medical history, food intake, medication usage, alcohol consumption, smoking history, and other factors that could affect health.

### Outcome data

The cardiac rhythm was measured for approximately 7 minutes. From November 2000 to December 2009, square root of the mean squared differences of successive normal-to-normal intervals (rMSSD), 5-min standard deviation (SD) of normal-to-normal intervals (SDNN), HF, LF, and LF/HF measurements were obtained during each participant’s regularly scheduled visit [11, 12]. To date, HRV measurements have been made for 700 participants, in some cases on multiple occasions (214 participants with one, 216 with two, 259 with three, and 11 with four measurements). A total of 1467 valid measurements were included in our analyses. All participants gave written consent and this study had Institutional Review Board approval.

### Exposure and weather data

We measured concentrations of ambient BC continuously at the Harvard University Countway Library stationary ambient monitoring site that is located < 1 km from the clinical laboratory where subjects were examined. Hourly ambient carbon monoxide (CO), nitrogen monoxide (NO), and nitrogen dioxide (NO_{2}) concentrations were obtained from the Massachusetts Department of Environmental Protection (MDEP) for monitoring stations within 20 km of the clinical laboratory, including those in Boston (Bremen, Kenmore Square, North-End, and Roxbury), Waltham and Lynn MA. The pollutant exposure was averaged across all the MDEP monitoring stations; data from the MDEP monitoring stations were used to calculate the mean hourly concentrations of each pollutant. We obtained data on dew point temperature and ambient temperature from the first order National Weather Service station at Boston Logan airport. We calculated apparent temperature (ATemp), a human discomfort index [13], as: ATemp = −2.653 + (0.994 × AT) + (0.0153 × DPT^{2}), where DPT is the dew point temperature in Celsius and AT is the ambient temperature in Celsius.

### Statistical analysis

We examined the association between latent traffic pollution and latent parasympathetic tone and between latent traffic pollution and sympathetic tone marker, LF/HF, using SEMs in a Bayesian framework that account for repeated measures. SEMs consist of two key components; a structural model part which shows the association dependencies between two latent constructs, between a latent construct and a measured variable or between two measured variables, and a measurement model part which shows the relation between latent constructs and their indicators/markers [9, 14]. HF, SDNN, rMSSD, and LF/HF were log transformed to satisfy the standard SEM assumption that all continuous variables are approximately normally distributed. Descriptive statistics and correlations between HRV markers, traffic-related pollutants and covariates were calculated.

#### Structural equation model with repeated measures in Bayesian framework

We examined associations between traffic-related pollution and markers of parasympathetic tone and sympathetic tone using SEMs, parametric Directed Acyclic Graphs (DAGs) that contain specified paths connecting observed and latent variables [9, 14]. SEMs in Bayesian framework consist of (1) distributional assumptions on the latent variables, (2) prior distributions for the parameters, and (3) probability models for the observed variables given latent variables and model parameters. SEMs also specify linear models relating the unknown means of the latent traffic pollution, sympathetic tone marker LF/HF and latent parasympathetic tone variables to fixed covariates (e.g. confounding factors), other latent variables, and a subject level random intercept to accommodate the repeated measures design of the study. For the effect of traffic pollution on both LF/HF and the latent parasympathetic tone, exposure periods were assessed using 4-hr, 24-hr, 48-hr and 72-hr moving average pre-visit, given results from previous studies that show effects at similar exposure windows [4, 11].

For the measurement part of the SEM, we specified measurement models describing the relationships between observed measures and latent variables. The relationship between the traffic-related observed measures (BC, CO, NO, and NO_{2}) and the latent exposure variable *Traffic* was first modeled (Figure 1).

For moving average concentration pre-visit measurement of participants *i* = 1, …*I*, of traffic-related pollutants *j* = 1, …*J*, on participant’s scheduled time visits *t* = 1, …*T*,

{X}_{\mathit{ijt}}~N\phantom{\rule{0.25em}{0ex}}\left({\mu}_{\mathit{ijt}}{}^{X},{\tau}_{j}{}^{X}\right)

(1)

where *X*
_{
ijt
} is the normally distributed traffic-related exposure variable for moving average concentration pre-visit measurement of participant *i* of traffic-related pollutant *j* on participant’s scheduled time visit *t*, with unknown mean *μ*
_{
ijt
}
^{X} and precision (1/variance) *τ*
_{
j
}
^{X}. The mean exposure of traffic-related pollutant *j*, *μ*
_{
ijt
}
^{X}, was related to the latent exposure variable *Traffic*
_{
it
} through a linear model,

{\mu}_{\mathit{ijt}}{}^{X}={\lambda}_{0j}+{\lambda}_{1j}\mathit{Traffi}{c}_{\mathit{it}}

(2)

where *Traffic*
_{
it
} is the normally distributed latent exposure variable *Traffic* for participant *i* for scheduled time visit *t*, with unknown mean *μ*
_{
it
}
^{Traffic} and precision (1/variance) *τ*
^{Traffic},

\mathit{Traffi}{c}_{\mathit{it}}~N\left({\mu}_{\mathit{it}}{}^{\mathit{Traffic}},{\tau}^{\mathit{Traffic}}\right).

(3)

The mean *Traffic* exposure, *μ*
_{
it
}
^{Traffic}, was related to confounding factors, through a linear model,

\begin{array}{l}{\mu}_{\mathit{it}}^{\mathit{Traffic}}={\phi}_{1}\mathit{Sine}\phantom{\rule{0.12em}{0ex}}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]+{\phi}_{2}\mathit{Cosine}\phantom{\rule{0.12em}{0ex}}[2\mathit{\pi \theta}\left(t\right)/365.24]\\ \phantom{\rule{5.3em}{0ex}}+{\phi}_{3}\mathit{Residual}{\left(\mathit{ATemp}\right)}_{t}+{\phi}_{4}\mathit{Residua}{l}^{2}{\left(\mathit{ATemp}\right)}_{t}\end{array}

(4)

where *θ* (t) is the calendar day of the year the participant visited the clinic, the sine and cosine terms capture seasonal trends in traffic pollution, ATemp is apparent temperature averaged across the same exposure window, and residual (ATemp) is the residual of a model regressing ATemp against sine and cosine terms of the calendar day of the year in the study. This approach was taken because temperature data was highly seasonal, and this avoids substantial collinearity between the seasonal terms and temperature, which would otherwise be present in the model. Furthermore, Eq. (4) looks like it doesn’t vary with participant *i* but it does vary because the scheduled time visits *t* are unique to that participant.

The relationship between the latent outcome variable, and its observed markers was also modeled. The HRV markers HF, SDNN and rMSSD were used as the markers for the latent *Parasympathetic tone* variable. For HRV measurement of participants *i* = 1, …*I*, of HRV parasympathetic tone markers *k* = 1, …*K*, on participant’s scheduled time visits *t* = 1, …*T*,

{Y}_{\mathit{ikt}}~N\left({\mu}_{\mathit{ikt}}^{Y},{\tau}_{k}^{Y}\right)

(5)

where *Y*
_{
ikt
} represents the normally distributed markers of parasympathetic tone for HRV measurement of participant *i* of parasympathetic tone marker *k* on participant’s scheduled time visit *t*, with unknown mean *μ*
_{
ikt
}
^{Y} and precision (1/variance) *τ*
_{
k
}
^{Y}. For parasympathetic tone marker *k*, the mean outcome, *μ*
_{
ikt
}
^{Y}, was related to the latent outcome variable *Parasympathetic tone* through a linear model,

{\mu}_{\mathit{ikt}}^{Y}={\alpha}_{0k}+{\alpha}_{1k}\mathit{Parasympathetic}\phantom{\rule{0.25em}{0ex}}\mathit{ton}{e}_{\mathit{it}}

(6)

where *Parasympathetic tone*
_{
it
} is the normally distributed latent outcome variable *Parasympathetic tone* for scheduled time visit *t* for participant *i*, with unknown mean *μ*
_{
it
}
^{Parasympathetic tone} and precision (1/variance) *τ*
^{Parasympathetic tone},

\begin{array}{l}\mathit{Parasympathetic}\phantom{\rule{0.25em}{0ex}}\mathit{ton}{e}_{\mathit{it}}~N\\ \phantom{\rule{7em}{0ex}}\left({\mu}_{\mathit{it}}^{\mathit{Parasympathetic}\phantom{\rule{0.25em}{0ex}}\mathit{tone}},{\tau}^{\mathit{Parasympathetic}\phantom{\rule{0.25em}{0ex}}\mathit{tone}}\right).\end{array}

(7)

Furthermore, for LF/HF measurement of participants *i* = 1, …*I*, on participant’s scheduled time visits *t* = 1, …*T*,

\mathrm{LF}/{\mathrm{HF}}_{\mathit{it}}~N\left({\mu}_{\mathit{it}}^{\mathrm{LF}/\mathrm{HF}},{\tau}^{\mathrm{LF}/\mathrm{HF}}\right)

(8)

where LF/HF_{
it
} is the normally distributed sympathetic tone marker outcome for participant *i* for scheduled time visit *t*, with unknown mean *μ*
_{
it
}
^{LF/HF} and precision (1/variance) *τ*
^{LF/HF}.

The mean *Parasympathetic tone* outcome, *μ*
_{
it
}
^{Parasympathetic tone} was related to the latent exposure variable *Traffic*
_{
it
}, confounding factors, predictors of the outcome variable, and a subject level random intercept, *a*
_{
i
}, through a linear model,

\begin{array}{l}{\mu}_{\mathit{it}}^{\mathit{Parasympathetic}\phantom{\rule{0.25em}{0ex}}\mathit{tone}}={\gamma}_{1}\mathit{Traffi}{c}_{\mathit{it}}+{\gamma}_{2}{X}_{1\mathit{it}}+\dots +{\gamma}_{m+1\phantom{\rule{0.25em}{0ex}}}{X}_{\mathit{mit}}\\ \phantom{\rule{1em}{0ex}}+{\gamma}_{m+2\phantom{\rule{0.25em}{0ex}}}\mathit{Sine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\phantom{\rule{0.25em}{0ex}}+{\gamma}_{m+3\phantom{\rule{0.25em}{0ex}}}\mathit{Cosine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\\ \phantom{\rule{1em}{0ex}}+{\gamma}_{m+4\phantom{\rule{0.25em}{0ex}}}\mathit{Residual}{\left(\mathit{ATemp}\right)}_{t}+{\gamma}_{m+5\phantom{\rule{0.25em}{0ex}}}\mathit{Residua}{l}^{2}{\left(\mathit{ATemp}\right)}_{t}+{a}_{i}.\end{array}

(9)

In addition, the mean LF/HF outcome, *μ*
_{
it
}
^{LF/HF}, was also related to the latent exposure variable *Traffic*
_{
it
}, predictors of the outcome variable, confounding factors, and a subject level random intercept, *b*
_{
i
}, through a linear model,

\begin{array}{l}{\mu}_{\mathit{it}}^{\mathrm{LF}/\mathrm{HF}}={\delta}_{1}\mathit{Traffi}{c}_{\mathit{it}}+{\delta}_{2}{X}_{1\mathit{it}}+\dots +{\delta}_{m+1\phantom{\rule{0.25em}{0ex}}}{X}_{\mathit{mit}}\\ \phantom{\rule{0.25em}{0ex}}+{\delta}_{m+2\phantom{\rule{0.25em}{0ex}}}\mathit{Sine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\phantom{\rule{0.25em}{0ex}}+{\delta}_{m+3\phantom{\rule{0.25em}{0ex}}}\mathit{Cosine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\\ \phantom{\rule{0.25em}{0ex}}+{\delta}_{m+4\phantom{\rule{0.25em}{0ex}}}\mathit{Residual}{\left(\mathit{ATemp}\right)}_{t}+{\delta}_{m+5\phantom{\rule{0.25em}{0ex}}}\mathit{Residua}{l}^{2}{\left(\mathit{ATemp}\right)}_{t}+{b}_{i}.\end{array}

(10)

For Eqs. (9) and (10), the fixed covariates (*X*
_{
1
}
*, …, X*
_{
m
}) included age, fasting blood glucose, body mass index (BMI), smoking status (ever/never), room temperature, sitting mean arterial blood pressure, ≥ 2 servings of alcohol per day (yes/no), and cardiac medication usage (beta blocker, angiotensin converting enzyme inhibitor, and calcium channel blocker). These covariates were chosen a priori as potentially important predictors of both *Parasympathetic tone* and LF/HF. To account for seasonal variation in *Parasympathetic tone* and LF/HF, a function of calendar date {sine [2πθ(t)/365.24] + cosine [2πθ(t)/365.24], where *θ* (t) is the calendar day of the year the participant visited the clinic} was used. To account for temperature variation in both *Parasympathetic tone* and LF/HF, we regressed ATemp on sine [2πθ(t)/365.24] + cosine [2πθ(t)/365.24] and used the linear [Residual (ATemp)_{t}] and quadratic [Residual^{2} (ATemp)_{t}] residuals of the equation as part of the fixed effects of the model. We also included a subject level random intercept, *a*
_{
i
} for Eq. (9) and *b*
_{
i
} for Eq. (10), to represent subject-specific permanent effects for each participant to accommodate the repeated measures design of the study. Moreover, Eqs. (9) and (10) implicitly depend on participant *i* because the scheduled time visits *t* depend on participant *i.* To evaluate obesity as an effect modifier, subjects were classified into two groups according to BMI (obese: BMI ≥ 30 vs. non-obese: BMI < 30). We also assessed effect modification by a history of diabetes (diabetic vs. non-diabetic) and by the four possible combinations of obesity and history of diabetes (diabetic and obese, diabetic and non-obese, non-diabetic and obese, and non-diabetic and non-obese). We included latent interaction term constructs between the dichotomized modifier variable and traffic pollution in Equations (9) and (10).

For model identifiability, the location of all the latent variables was centered and set to 0 (i.e. *δ*
_{
0
}
*, φ*
_{
0
} and *γ*
_{
0
} = 0). A scale must be specified for each latent variable because latent variables were constructed from multiple measured parameters, which have different scales of variation. We chose HF as the reference scale for *Parasympathetic tone* and BC as the reference scale for *Traffic* pollution. In order to set the reference scale for *Parasympathetic tone* and *Traffic* and for identifiability and interpretation of these latent variables, α_{1}
^{HF} and λ_{1}
^{BC}, were also constrained to be equal to 1 [6]. For the DAG diagram of the SEMs, see Figure 1.

#### Linear mixed models in Bayesian framework

In order to compare the results from the SEM in Bayesian framework with repeated measures to those obtained from a standard univariate analysis, associations between traffic-related air pollutant BC and mean markers of parasympathetic tone (HF) and sympathetic tone (LF/HF) were estimated using linear mixed models with random subject-specific intercepts in Bayesian linear mixed models (BLMMs) framework. Additionally, we compared the results from the above BLMMs to Frequentist framework estimates (Frequentist linear mixed models, LMMs) obtained from LME in R (http://www.r-project.org). The same fixed and random covariates used in SEMs were included in all the LMMs and BLMMs. To allow more direct comparison for the SEMs, BLMMs, and LMMs, effect size estimates were reported for an inter-quartile range (IQR) increase in BC (reference pollutant of latent variable *Traffic* pollution).

We specified a model for the measured health effect outcomes, HF and LF/HF. The means of both HF and LF/HF for participant *i* on participant’s scheduled time visit *t*, *E[HF*
_{
it
}
*]* and *E[LF/HF*
_{
it
}
*]*, respectively, was related to BC, confounding factors and predictors of the health effect outcomes, and a subject level random intercept, *c*
_{
i
} or *d*
_{
i
}, through a linear model,

\begin{array}{l}E\left[H{F}_{\mathit{it}}\right]={\beta}_{0}+{\beta}_{1}B{C}_{\mathit{it}}+{\beta}_{2}{X}_{1\mathit{it}}+\dots +{\beta}_{m+1\phantom{\rule{0.25em}{0ex}}}{X}_{\mathit{mit}}\\ \phantom{\rule{2.75em}{0ex}}+{\beta}_{m+2}\mathit{Sine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\phantom{\rule{0.25em}{0ex}}+{\beta}_{m+3}\mathit{Cosine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\\ \phantom{\rule{1.75em}{0ex}}+{\beta}_{m+4}\mathit{Residual}{\left(\mathit{ATemp}\right)}_{t}+{\beta}_{m+5}\mathit{Residua}{l}^{2}{\left(\mathit{ATemp}\right)}_{t}+{\mathrm{c}}_{i}\end{array}

(11)

\begin{array}{l}E\left[\mathit{LF}/H{F}_{\mathit{it}}\right]={\rho}_{0}+{\rho}_{1}B{C}_{\mathit{it}}+{\rho}_{2}{X}_{1\mathit{it}}+\dots +{\rho}_{m+1\phantom{\rule{0.25em}{0ex}}}{X}_{\mathit{mit}}\\ \phantom{\rule{0.75em}{0ex}}+{\rho}_{m+2}\mathit{Sine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\phantom{\rule{0.25em}{0ex}}+{\rho}_{m+3}\mathit{Cosine}\left[2\mathit{\pi \theta}\left(t\right)/365.24\right]\\ \phantom{\rule{0.25em}{0ex}}+{\rho}_{m+4}\mathit{Residual}{\left(\mathit{ATemp}\right)}_{t}+{\rho}_{m+5}\mathit{Residua}{l}^{2}{\left(\mathit{ATemp}\right)}_{t}+{d}_{i}.\end{array}

(12)

In terms of Bayesian inference, results were quantified and characterized in terms of the strength of the following hypothesized relationships: (1) inverse relationships between BC and HRV marker HF, and between *Traffic* pollution and *Parasympathetic tone*, and (2) positive relationships between BC and the HRV marker LF/HF, and between *Traffic* pollution and sympathetic tone marker, LF/HF. In addition, results were presented as probability statements, providing evidence that coefficients describing *Traffic* pollution (*γ*
_{
1
} and *δ*
_{
1
}) or BC (*β*
_{
1
} and *ρ*
_{
1
}) are either negative or positive depending on the hypothesized relationship. See Additional file 1 for details for implementation particulars and description of prior distributions of parameters, and see Additional file 2 for the WinBUGS code.