2.1.

Data Acquisition

We consider two datasets containing both DTI and T1-weighted images: the Baltimore Longitudinal Study of Aging (BLSA)^25,26 and the Vanderbilt Memory and Aging Project (VMAP).²⁷ VMAP data were collected by Vanderbilt Memory and Alzheimer’s Center Investigators at Vanderbilt University Medical Center. BLSA DTI scans were acquired on a 3T Phillips scanner in 32 directions at a $b$-value of $700\mathrm{s}/{\mathrm{mm}}^2$ with voxel dimensions of $2.2\times 2.2\times 2.2{\mathrm{mm}}^3$ that were resampled to $0.8125\times 0.8125\times 2.2{\mathrm{mm}}^3$. All VMAP DTI scans considered were acquired on a 3T Phillips scanner with an 8ch SENSE head coil in 32 directions at a $b$-value of $1000\mathrm{s}/{\mathrm{mm}}^2$ with voxel dimensions of $2\times 2\times 2{\mathrm{mm}}^3$. For both sites, T1 scans were acquired in the same scanning session as the DTI images. BLSA T1 scans were acquired with voxel dimensions of $1.2\times 1\times 1{\mathrm{mm}}^3$ and VMAP T1 scans were acquired with voxel dimensions of $1\times 1\times 1{\mathrm{mm}}^3$.

2.2.

Silver Standard Cohort

We consider both cognitively unimpaired and mild cognitive impairment (MCI) participants across VMAP and BLSA, matched by cognitive status, sex, and age within 4 years while keeping track of APOE2 positivity, APOE4 positivity, race/ethnicity, and years of education. The final size of our silver standard cohort is $N=358$, with each site contributing an equal number of participants (Table 1).

Table 1

Demographic information for each dataset in the silver standard cohort.

2.3.

Data Processing and Pre-processing

All DTI scans are preprocessed using v 1.0.8 of the PreQual preprocessing pipeline²⁸ for denoising and to remove susceptibility-induced and eddy current distortions (Fig. 1). Fractional anisotropy (FA) and MD are calculated from the preprocessed data. The EVE Type-I, EVE Type-II, and EVE Type-III JHU atlases^29,30 are registered to the diffusion space of each participant using ANTs SyN registration³¹ and FSL’s epi_reg.³² The T1 brain mask used for epi_reg is calculated via SLANT.³³ The epi_reg transform is converted to the same format used by ANTs via the Convert3D³⁴ tool developed by the ITK-SNAP team. The two transformations are then applied in a single registration step. After registration of the atlases to the diffusion space, mean FA, MD, AD, and RD are calculated for each region of the three atlases using the MRtrix3 software.³⁵ We use the regions of interest (ROIs) from the EVE Type-III atlas, as they are all WM regions. The code used for the process is available at: ( https://github.com/MASILab/AtlasToDiffusionReg).

Fig. 1

After registration of the JHU EVE-III Atlas, mean FA values were calculated in all the regions for each participant in the silver standard cohort. A point in the experimental space is “feasible” if the sample size for either site is at least $N=6$, the imbalance level does not result in $N$ for either site exceeding the available number of participants for that site, and if sampling of participants yielded a covariate shift within 1 year of the target age difference between sites. For each feasible point in the experimental space, 10 bootstraps were subsampled from the silver standard cohort, and the FA values for the subsamples were harmonized by ComBat. The resulting parameters were then compared to those from the silver standard to determine reliability of ComBat at that location in the experimental space.

2.4.

ComBat Harmonization

The ComBat model proposed by Johnson et al.¹⁴ assumes that $Y_{isf}$, the original input scalars for feature $f$ extracted from participant scan $i$ that comes from site s are defined as

We use (i) the mean FA values for 112 of the 118 regions of the EVE Type-III atlas for each participant; (ii) covariates of age, sex, cognitive status, race, education, APOE2 carrier status, and APOE4 carrier status, the covariates used in DTI harmonization by Yang et al.;³⁶ and (iii) a covariate indicating the dataset the participant scan came from. Six regions are excluded from consideration because the registration process resulted in regions with zero volume, which could induce a large shift in the mean of the distribution of mean FA values for the cohort. Such shifts could have substantial impact on the ComBat harmonization procedure. The education covariate is a continuous variable indicating years of education, and the APOE2 and APOE4 covariates are categorical variables indicating the presence or absence of the respective APOE allele. The cognitive status covariate indicates whether the participant is cognitively unimpaired or is diagnosed with MCI. The version of ComBat used for this analysis is implemented by Fortin et al. ( https://github.com/Jfortin1/neuroCombat).¹⁵

2.5.

Experimental Search Space

Bayer et al. highlighted a variety of factors that might influence the performance of ComBat.¹⁷ To examine the robustness of ComBat, we use 740 different permutations of 19 levels of total sample size, 10 levels of sample size imbalance, and 6 levels of covariate shift (Fig. 1). Total sample size is the number of participants whose features are input to ComBat for harmonization. In the context of this analysis, we define sample size imbalance to be the ratio $X:10$ of participants from one site relative to the other site, where 10:10 is perfect balance of sample size between sites. We consider levels from $X=1$ to $X=10$. To compare different experimental permutations, we keep VMAP as the single site whose sample size changes with respect to BLSA. For the covariate shift, the value is the difference in the mean age of the participants from a single dataset compared to the other one. Fortin et al. demonstrated ComBat harmonization of datasets with different age ranges, but we evaluate performance at multiple levels to provide a more controlled and comprehensive analysis of changes in performance.³⁷ For all levels of difference, we consider scenarios both when the mean age for VMAP is greater than BLSA and when the mean age of VMAP is less.

Fig. 2

The root mean squared error (RMSE) of standardized ${\beta }_{\mathrm{AGE}}$ estimates for mean FA versus age compared to the silver standard indicate that ComBat is not stable with all experimental permutations considered, as the error increases when the cohort changes to have an average mean age difference between VMAP and BLSA of (a) 0 years, (b) 2 years, (c) 4 years, (d) 6 years, (e) 8 years, and (f) 10 years. The values represent the mean normalized RMSE across EVE Type-III Atlas regions averaged across 10 iterations of each feasible point in the experimental space. For each subplot, total sample size of the cohort is on the $x$-axis and sample size imbalance is on the $y$-axis, where $Y:10$ represents $Y$ participants at VMAP for every 10 at BLSA. Any non-feasible experimental permutations are represented in gray.

2.6.

Bootstrapping Experimental Space

We bootstrap 10 simulations at each of these experimental permutations through sampling from the silver standard cohort without replacement. For each iteration, ComBat harmonization is performed as described in Sec. 2.4 to obtain the harmonized data and estimated parameters ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$ for each region and site combination. The mean FA vs Age regressions are calculated according to Sec. 2.7. Bootstrapping subsamples of the silver standard can be done either with or without replacement. As sampling with replacement can introduce artifacts in the data, such as repeated participants, we choose to sample without replacement. As a consequence, not all permutations in the experimental space defined in Sec. 2.5 are feasible. In addition, the subsamples at some permutations will be more highly correlated than others due to lack of participant variability for a site that fit the experimental criteria. We also do not include any permutations that result in any site having fewer than six participants.

2.7.

Mean FA Versus Age Regression

To estimate the associations of FA with age from the harmonized data, we perform a linear regression

2.8.

Comparison to Silver Standard

To evaluate the robustness of ComBat at each experimental permutation, we compare to the silver standard cohort using three different error metrics: (i) the average root mean square error in normalized effect size for ${\beta }_{\mathrm{AGE}}$ across all regions, (ii) the average root mean squared error of ${\widehat{\gamma }}_{sf}^{*}$ across all regions, and (iii) the average root mean squared error of the log of ${\widehat{\delta }}_{sf}^{*}$ across all regions. The standard error for each regression estimate is also obtained from the linear regression estimation. We normalize ${\beta }_{\mathrm{AGE}}$ for each region by dividing by its respective silver standard standard error value from the regression estimation to compare ${\beta }_{\mathrm{AGE}}$ effect size across all regions for the permutations. Unlike ${\beta }_{\mathrm{AGE}}$, we cannot get silver standard estimates for the standard errors of ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$, so we cannot normalize the values with reference to the silver standard. Additionally, since ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$ are estimated iteratively via the EM algorithm, they are highly dependent on each other, and normalizing them independently may introduce bias. Thus, we leave ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$ unchanged from ComBat for comparisons to the silver standard. We also look at the average root mean squared errors of the differences between ${\widehat{\gamma }}_{\mathrm{BLSA},f}^{*}$ and ${\widehat{\gamma }}_{\mathrm{VMAP},f}^{*}$ and between the log of ${\widehat{\delta }}_{\mathrm{BLSA},f}^{*}$ and the log of ${\widehat{\delta }}_{\mathrm{VMAP},f}^{*}$ across all regions in order to assess the relative scalings and shifts of the feature distributions.

2.9.

Checking Assumptions

To determine if the instability of ComBat is related to the assumptions that the model makes, we assess the following.

Note that the ComBat model for the data in Eq. (1) represents biological covariates as causing linearly independent variation in the data. Further, ComBat assumes that the biological variation in the data is separable from the variation due to site/scanner biases, which requires the biological covariates to not be strongly correlated with the batch variable. We assess (1) and (2) by evaluating the negative log likelihoods of the residual distributions compared to the assumed prior distributions. For the residuals of the model and ${\widehat{\gamma }}_{sf}^{*}$, we use the mean and standard deviations of the empirical distributions to generate a normal distribution fit to the data and assess the average negative log likelihood across all empirical data points compared to the generated normal distribution. For ${\widehat{\delta }}_{sf}^{*}$, we estimate the inverse gamma distribution as

2.10.

Comparison to Other Linear Models

To assess the stability of ComBat compared to other linear models, we also perform the same analysis of root mean squared error in normalized ${\beta }_{\mathrm{AGE}}$ compared to the silver standard on the same data subsets of ComBat experimental runs for an ordinary least square (OLS) and linear mixed effects (LME) model. The OLS is the same as Eq. (6), but without the prior ComBat harmonization of the data. The LME is modeled as

3.1.

Total Sample Size and ${\beta }_{\mathrm{AGE}}$

To evaluate ComBat’s reliability with respect to $N$, we plot root mean squared error of normalized ${\beta }_{\mathrm{AGE}}$ for each ROI averaged across all 10 iterations against the total sample size (Fig. 4). The total sample size is spaced logarithmically, so the expected trend in error as $N$ increases will be a linear decrease with respect to logarithmic increases in sample size. At small sample sizes, the decrease in error does not follow this trend. We only consider experimental permutations that have a zero-year covariate shift and a 10:10 imbalance in sample size. We consider ComBat to be “stable” as long as ${\beta }_{\mathrm{AGE}}$ for 50% of the ROIs are within one standard deviation of the respective silver standard, and “unstable” otherwise. For harmonization of mean FA values, we observe that ComBat becomes unstable at sample sizes of $N<162$. For a more conservative threshold, where all ROIs are within one standard deviation away, we observe instability for $N<252$.

Fig. 4

(a) For decreasing sample size, the expected trend for a well-behaved model is error increasing by a factor of $\sqrt{N}$. Thus, the trend in error in the logarithmic space is linear with increasing $N$. We consider ComBat to be “stable” with respect to $N$ when the RMSE for standardized ${\beta }_{\mathrm{AGE}}$ of mean FA versus age compared to the silver standard is below 1 (blue line) or the error follows the error trend stated. On this criterion, we suggest that ComBat is unstable when $N<162$, as over 50% of the errors of ROIs for $N=144$ are above 1 and the increase in error in not linear in the logarithmic space. (b) RMSE for standardized ${\beta }_{\mathrm{AGE}}$ of mean MD versus age shows ComBat remaining stable for $N>162$, indicating that different DTI scalars have different levels of sensitivity to changes in sample size for ComBat. To observe the effect of only sample size on ComBat, only permutations with an imbalance level of 10:10 and no covariate shift were considered.

3.2.

Sample Size Imbalance and ${\beta }_{\mathrm{AGE}}$

To examine the error in ${\beta }_{\mathrm{AGE}}$ estimation with respect to sample size imbalance, we plot the root mean squared error of ${\beta }_{\mathrm{AGE}}$ for each ROI averaged across all 10 iterations against the sample size imbalance for both sample sizes that would make ComBat stable and unstable (Fig. 5, Fig. S3 in Supplementary Material 1). We only consider experimental permutations that have a zero-year covariate shift. We do not observe sample size imbalance to have an effect on the estimation of ${\beta }_{\mathrm{AGE}}$, as the error fluctuates for all levels of imbalance.

Fig. 5

RMSE for standardized ${\beta }_{\mathrm{AGE}}$ of mean FA versus age compared to the silver standard for sample sizes of (a) $N=306$, (b) $N=288$, (c) $N=270$, (d) $N=252$, (e) $N=234$, (f) $N=216$, (g) $N=198$, (h) $N=180$, (i) $N=162$, (j) $N=144$, (k) $N=126$, and (l) $N=108$. Covariate shift does not seem to have a definitive threshold at which the error in estimation of ${\beta }_{\mathrm{AGE}}$ is much larger compared to the respective experimental run with no covariate shift. For an all-encompassing threshold ambiguous to the size of $N$, we suggest a maximum covariate shift of 2 years between sites because a covariate shift of either 4 or 6 years increases the error in estimation of ${\beta }_{\mathrm{AGE}}$ depending on $N$. Only experimental permutations that have an imbalance ratio of 10:10 were considered.

3.3.

Covariate Shift and ${\beta }_{\mathrm{AGE}}$

To examine the error in ${\beta }_{\mathrm{AGE}}$ estimation with respect to covariate shift, we plot the root mean squared error of ${\beta }_{\mathrm{AGE}}$ for each ROI averaged across all 10 iterations against the covariate shift (Fig. 6, Fig. S4 in Supplementary Material 1). We only consider experimental permutations that have a 10:10 sample size imbalance level. For stable sample sizes, we observe stability along the covariate shift axis at mean age differences of up to 2 to 4 years between sites. As this threshold fluctuates among sample sizes, we consider a conservative threshold at 2 years and a looser threshold at 4 years. The average effect size for each level of covariate shift can be found in Table S1 in Supplementary Material 1.

Fig. 6

Sample size imbalance alone does not substantially affect estimation of ${\beta }_{\mathrm{AGE}}$ for mean FA harmonization; only at smaller $N$ does it appear to have an effect. However, this is likely due to the small sample size at these experimental permutations. Only experimental permutations that have a covariate shift of 0 years were considered. Panels (a)–(l) are the same $N$ as Fig. 5.

3.4.

Comparison to Silver Standard - ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$

To examine the error in harmonization, we perform a comprehensive visualization of root mean squared error of ${\widehat{\gamma }}_{sf}^{*}$ and the log of ${\widehat{\delta }}_{sf}^{*}$ for each site at each feasible experimental parameter (Figs. 7 and 8). Unlike sample size with ${\beta }_{\mathrm{AGE}}$, we do not have expectations of well-behaved models for changes in ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$. As the ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$ estimates are the measures of site bias, we consider ComBat to be unstable with any deviation from the silver standard values, as any error would result either in site bias not being removed completely or removing variability attributed to biological factors. The gradation in color indicates increasing error as we move further away from the silver standard.

Fig. 7

${\widehat{\gamma }}_{sf}^{*}$ and $\log {\widehat{\delta }}_{sf}^{*}$ do not follow the same trend in error as ${\beta }_{\mathrm{AGE}}$ for harmonization of mean FA values. (Top left) RMSE error in ${\widehat{\gamma }}_{sf}^{*}$ estimates for VMAP averaged across ROIs with total sample size on the $y$-axis and sample size imbalance on the $x$-axis at covariate shift levels of (a) 0 years, (b) 2 years, (c) 4 years, (d) 6 years, (e) 8 years, and (f) 10 years. (Top right) Error in $\log {\widehat{\delta }}_{sf}^{*}$ estimates for VMAP averaged across ROIs and plotted in the same order as top left. Bottom left and bottom right plots are the RMSE in ${\widehat{\gamma }}_{sf}^{*}$ and $\log {\widehat{\delta }}_{sf}^{*}$ estimates for BLSA respectively, with slices along the covariate shift axis presented the same as top left. For ComBat to accurately estimate site bias, ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$ should be as close to the silver standard values as possible. Thus, we suggest a maximum covariate shift of 2 years, an imbalance of 9:10 and a total sample size of $N>252$ for stability in both ComBat parameters, as these experimental parameters have a relative error close to zero.

Fig. 8

RMSE error in estimates for (top left) VMAP ${\widehat{\gamma }}_{sf}^{*}$, (top right) VMAP $\log {\widehat{\delta }}_{sf}^{*}$, (bottom left) BLSA ${\widehat{\gamma }}_{sf}^{*}$, and (bottom right) BLSA $\log {\widehat{\delta }}_{sf}^{*}$ for harmonization of mean MD values increase much more quickly than they do for harmonization of FA values. Slices along the covariate shift axis are plotted similarly to Fig. 7.

For ${\widehat{\gamma }}_{sf}^{*}$, in terms of sample size, we observe a threshold of around $N\ge 252$ for ComBat stability, similar to the conservative threshold for sample size in error estimation of ${\beta }_{\mathrm{AGE}}$; along the sample size imbalance axis, we observe a threshold of around a 10 : 9 imbalance ratio between sites; and along the covariate shift axis, we observe a threshold of around four years. For ${\widehat{\delta }}_{sf}^{*}$, in terms of sample size, we observe a threshold of around $N\ge 308$ for ComBat stability; along the sample size imbalance axis, any imbalance in sample size results in ComBat instability; and along the covariate shift axis, we observe a threshold of up to 2 years in mean age difference between sites for ComBat stability.

For the shift differences compared to those of the silver standard, we observe changes in ${\widehat{\gamma }}_{\mathrm{BLSA},f}^{*}-{\widehat{\gamma }}_{\mathrm{VMAP},f}^{*}$ that are less consistently responsive to changes in the experimental parameters (Fig. 9, Fig. S5 in Supplementary Material 1) than for the site-wise root mean squared errors of ${\widehat{\gamma }}_{sf}^{*}$ for both sites. Still, we observe a similar threshold of around $N\ge 252$ for ${\widehat{\gamma }}_{sf}^{*}$ in terms of respective shift stability, a threshold of 9:10 for imbalance ratio between sites, and a covariate shift threshold of four years. For $\log ({\widehat{\delta }}_{\mathrm{BLSA},f}^{*})-\log ({\widehat{\delta }}_{\mathrm{VMAP},f}^{*})$, we similarly observe less consistent responses to changes in the experimental parameters than for the site-wise root mean squared errors. Similar to the site-wise estimation errors, we observe that any imbalance results in ComBat instability and a threshold of 2 years for imbalance ratio and covariate shift, respectively. However, we observe a threshold of $N\ge 252$ for $\log ({\widehat{\delta }}_{sf}^{*})$ in terms of respective scaling stability, which is less restrictive than that of the lone site-wise estimates.

Fig. 9

RMSE for ${\widehat{\gamma }}_{\mathrm{BLSA},f}^{*}-{\widehat{\gamma }}_{\mathrm{VMAP},f}^{*}$ of experimental runs compared to the silver standard averaged across ROIs shows a similar threshold for stability of around $N\ge 252$ as the lone RMSE of ${\widehat{\gamma }}_{sf}^{*}$ for both sites. However, the stability of the difference for $\log {\widehat{\delta }}_{\mathrm{BLSA},f}^{*}-\log {\widehat{\delta }}_{\mathrm{VMAP},f}^{*}$ shows a looser threshold of $N\ge 252$ compared to the lone RMSE of $\log {\widehat{\delta }}_{sf}^{*}$ for both sites for mean FA harmonization.

3.5.

Analysis of Assumptions

We visualize distributions of experimental run residuals for mean FA compared to the silver standard residuals in Fig. 10 and plot the average of the negative log likelihoods of residuals in Fig. 11 (see Figs. S6 and S7 in Supplementary Material 1 as well). We do see an increase in the negative log likelihoods as we decrease sample size, in spite of an increase in error for the estimation of ${\beta }_{\mathrm{AGE}}$. At the smallest sample sizes, where the estimation for ${\beta }_{\mathrm{AGE}}$ has the most error, the negative log likelihood is at its most negative values, indicating the residuals for small sample sizes are, on average, more normally distributed than other experimental bootstraps and the silver standard. The results from the Anderson–Darling test are also in agreement with the negative log likelihood analysis for sample size, as the residuals are normally distributed more often at lower sample sizes, which do not follow the trend of increasing error of ${\beta }_{\mathrm{AGE}}$ with decreasing sample size (Fig. 12). We see a decrease in the negative log likelihoods as we move further from the silver standard along the covariate shift axis. We see a slight increase in the negative log likelihoods as we move further along the sample size imbalance axis, indicating that sample size imbalance may lead to violation of the assumptions of ComBat. However, the Anderson–Darling results suggest that neither imbalance ratio nor covariate shift have a consistent effect on the normality of the residuals.

Fig. 10

The residuals from the ComBat model for the silver standard (blue, $N=358$) do not adhere to the assumption of normality given the heavy left tail. As we step further from the silver standard in terms of sample size, imbalance, and covariate shift, we expect the residuals to become even less normally distributed if the assumption of residual normality directly impacts the error in ${\beta }_{\mathrm{AGE}}$ for experimental runs (blue). Decreasing sample size does not appear to consistently lessen the tail of the residual distributions. The residuals plotted above are for the left genu of the corpus callosum for mean FA harmonization.

Fig. 11

(a)–(f) The average negative log likelihood that the residual distribution follows a normal distribution (with mean and standard deviation estimated from the residual distribution) decreases as we decrease in sample size, and is smallest when the sample size is $<50$, suggesting that the residual distributions are more normal at low sample sizes for mean FA harmonization. This contrasts with the increasing error of ${\beta }_{\mathrm{AGE}}$ with decreasing sample size, suggesting that looking at the distribution of the assumptions alone cannot indicate if ComBat is appropriate for removing site biases of the given input cohort. Thus, we suggest the bootstrapping methodology to determine reliability of ComBat for site bias removal. Difference of average negative log likelihoods for residual distributions between experimental runs and the silver standard negative log likelihoods (averaged across all ROIs). Slices along the covariate shift axis are plotted similarly to Fig. 2.

Fig. 12

The Anderson–Darling test results suggest that the increase in error for ${\beta }_{\mathrm{AGE}}$, and thus the stability of ComBat, cannot be assessed by a decrease in normality of the residuals for either mean FA (top) or mean MD (bottom) harmonization.

We visualize distributions of ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$ obtained from mean FA harmonization as kernel density estimates compared to respective prior distributions in Fig. 13 and plotted the average negative log likelihoods of both parameters in Figs. 14 and 15. Corresponding plots for mean MD harmonization can be found in Figs. S8–S10 in Supplementary Material 1. For ${\widehat{\gamma }}_{sf}^{*}$, total sample size does not appear to have an effect on the shape of the distributions for experimental runs until $N=36$, at which point there is a large increase in the negative log likelihoods. For mean FA, we observe that ${\widehat{\gamma }}_{sf}^{*}$ distributions are not normally distributed until around $N=36$ according to the Anderson–Darling test (Fig. 16), which agrees with the negative log likelihood values. The results from the Anderson–Darling test for mean MD harmonization are plotted in Fig. S11 in Supplementary Material 1. However, sample size imbalance appears to affect the distributions, as VMAP ${\widehat{\gamma }}_{sf}^{*}$ deviate from normal and BLSA ${\widehat{\gamma }}_{sf}^{*}$ tend to normal as we move further from the silver standard along the covariate shift axis. Covariate shift does not appear to affect the distribution of ${\widehat{\gamma }}_{sf}^{*}$.

Fig. 13

Despite changes in the experimental parameters compared to the silver standard, the empirical distributions (dotted blue) of the ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}^{*}$ parameters for both sites in FA harmonization closely follow the respective estimated distributions. This adherence to the assumptions on ComBat does not correlate with the respective increasing error in estimates of ${\beta }_{\mathrm{AGE}}$ for experimental iterations that are further from the silver standard along the experimental axes. ${\widehat{\gamma }}_{sf}^{*}$ is assumed to have a normal distribution per site ${\widehat{\delta }}_{sf}^{*}$ is assumed to have an inverse gamma distribution.

Fig. 14

Similar to the ComBat model residuals, the VMAP average negative log likelihoods of the ${\widehat{\gamma }}_{sf}^{*}$ estimates (top) for a normal distribution and ${\widehat{\delta }}_{sf}^{*}$ estimates (bottom) for an inverse gamma distribution of experimental runs compared to the respective silver standard values do not correlate with the error trends for ${\beta }_{\mathrm{AGE}}$ in mean FA harmonization. This suggests that we cannot look at the distributions of ${\widehat{\gamma }}_{sf}^{*}$ and ${\widehat{\delta }}_{sf}$ to examine whether the input cohort is suitable for ComBat harmonization based on the premise of a violation of ComBat assumptions.

Fig. 15

BLSA average negative log likelihoods for ${\widehat{\gamma }}_{sf}^{*}$ (top) and ${\widehat{\delta }}_{sf}^{*}$ (bottom) compared to the respective silver standard values do not correlate with the error trends for ${\beta }_{\mathrm{AGE}}$ in mean FA harmonization.

Fig. 16

The Anderson–Darling test results suggest that the increase in error for ${\beta }_{\mathrm{AGE}}$, and thus the stability of ComBat, cannot be assessed by a decrease in normality of ${\widehat{\gamma }}_{sf}^{*}$ for mean FA harmonization.

For ${\widehat{\delta }}_{sf}^{*}$, neither covariate shift nor sample size appears to affect the negative log likelihoods of the distributions. For VMAP, imbalance appears to have no effect on ${\widehat{\delta }}_{sf}^{*}$, whereas a greater imbalance makes BLSA ${\widehat{\delta }}_{sf}^{*}$ closer to an inverse gamma distribution, indicating that a larger proportion of the total sample size will make ${\widehat{\delta }}_{sf}^{*}$ distributed more similarly to an inverse gamma distribution. We also calculate the correlation matrix of covariates for the silver standard cohort (Fig. 17) and do not find unexpectedly high correlations between any covariates.

Fig. 17

For the silver standard, the covariates are not highly correlated. We expect a correlation between age and diagnosis, as participants with MCI are more likely to be older. The correlations between race/ethnicity and site are non-zero but are substantially $<1$.

3.6.

Comparison to Other Linear Models

We do not observe asymmetry in the pairwise comparison of model errors for either mean FA or mean MD for ComBat versus OLS and LME (Figs. 18 and 19). Practically speaking, ComBat is similar to both LME and OLS when trying to harmonize DTI datasets with matched pairs.

Fig. 18

Pairwise comparisons for all experimental permutations of ComBat to OLS and LME models for ${\beta }_{\mathrm{AGE}}$ of mean FA suggest that ComBat is very similar to other linear models for pairwise DTI harmonization.

Fig. 19

Pairwise comparisons for all experimental permutations of ComBat to OLS and LME models for ${\beta }_{\mathrm{AGE}}$ of mean MD suggest that ComBat is very similar to other linear models for pairwise DTI harmonization. For permutations with very large error, ComBat may have slightly more stability.