The role of gut microbiome and its interaction with arsenic exposure in carotid intima-media thickness in a Bangladesh population.

BACKGROUND: Emerging data suggest that inorganic arsenic exposure and gut microbiome are associated with the risk of cardiovascular disease. The gut microbiome may modify disease risk associated with arsenic exposure. Our aim was to examine the inter-relationships between arsenic exposure, the gut microbiome, and carotid intima-media thickness (IMT)-a surrogate marker for atherosclerosis.
METHODS: We recruited 250 participants from the Health Effects of Arsenic Longitudinal Study in Bangladesh, measured IMT and collected fecal samples in year 2015-2016. 16S rRNA gene sequencing was conducted on microbial DNA extracted from the fecal samples. Arsenic exposure was measured using data on arsenic concentration in drinking water wells over time to derive a time-weighted water arsenic index. Multivariable linear regression models were used to test the inter-relationships between arsenic exposure, relative abundance of selected bacterial taxa from phylum to genus levels, and IMT.
RESULTS: We identified nominally significant associations between arsenic exposure, measured using either time-weighted water arsenic or urinary arsenic, and the relative abundances of several bacterial taxa from the phylum Tenericutes, Proteobacteria, and Firmicutes. However, none of the associations retained significance after correction for multiple testing. The relative abundances of the family Aeromonadaceae and genus Citrobacter were significantly associated with IMT after correction for multiple testing (P-value = 0.02 and 0.03, respectively). Every 1% increase in the relative abundance of Aeromonadaceae and Citrobacter was related to an 18.2-μm (95% CI: 7.8, 28.5) and 97.3-μm (95% CI: 42.3, 152.3) difference in IMT, respectively. These two taxa were also the only selected family and genus using the LASSO variable selection method. There was a significant interaction between Citrobacter and time-weighted water arsenic in IMT (P for interaction = 0.04).
CONCLUSIONS: Our findings suggest a role of Citrobacter in the development of atherosclerosis, especially among individuals with higher levels of arsenic exposure.

Introduction

Millions of people in the U.S. and worldwide are exposed to elevated levels of inorganic arsenic from drinking contaminated water. Chronic arsenic exposure has been associated with increased risks of preclinical and clinical endpoints of CVD (Chen et al., 2011; Chen et al., 2013a; Chen et al., 2013b; Chen et al., 2013c; Jiang et al., 2015; Moon et al., 2012; Moon et al., 2013). However, studies that address inter-individual variability in the association between arsenic exposure and CVD-related endpoints are needed to identify susceptibility. The assessment of susceptibility can stress the potential importance of the cardiovascular effects of arsenic exposure, as the effects may be stronger in population subgroups with susceptibility. Evidence suggests that arsenic can induce its effects on CVD through its induction of oxidative stress, inflammatory responses, and endothelial dysfunction (Bunderson et al., 2002; Kumagai and Pi, 2004; Simeonova and Luster, 2004).

The gut microbiome comprises the collective genome of approximately 100 trillion microbes residing in the gastrointestinal tract (gut microbiota), which is involved in normal physiological functions, such as digestion and metabolism (Ley et al., 2006). Alterations in the composition of the gut microbiome have been implicated in the pathogenesis of several risk factors of cardiovascular disease (CVD) ranging from obesity (Ley et al., 2006; Tilg and Kaser, 2011) to insulin resistance (Pedersen et al., 2016) and type 2 diabetes (Forslund et al., 2015; Qin et al., 2012). However, the evidence is not well-established and additional population-based studies are needed.

The gut microbiome may also interact with other external factors, such as diet, lifestyle, and environmental exposures (Marques et al., 2018; Singh et al., 2017; Wong, 2014), in affecting the cardiovascular health. For instance, several studies have suggested that the gut microbiome metabolizes dietary choline, phosphatidylcholine, and l-carnitine to produce trimethylamine-N-oxide (TMAO)—a metabolite that has been strongly linked to CVD risk (Koeth et al., 2013; Tang et al., 2013; Wang et al., 2011). Recent experimental studies suggest that the gut microbiome may contribute to presystemic biotransformation of ingested arsenic (Hall et al., 1997; Kubachka et al., 2009; Kuroda et al., 2004; Pinyayev et al., 2011; Rowland and Davies, 1981; Van de Wiele et al., 2010). Arsenic exposure may affect microbes in the gut and the associated metabolic activity (Dheer et al., 2015; Lu et al., 2014). However, very few epidemiologic studies have been conducted on these associations (Dong et al., 2017). Epidemiologic studies on the interaction between the gut microbiome and environmental exposures such as arsenic exposure in disease risk are also lacking.

Carotid artery intima-media thickness (IMT) is a non-invasive measure of subclinical atherosclerosis (de Groot et al., 2004; Stein et al., 2008) and predictive of incident CVD (Chambless et al., 2000; O'Leary et al., 1999). Studies using IMT as the endpoint therefore may reveal risk factors for atherosclerosis and CVD. Instead of the binary nature of cardiovascular events, IMT describes the preclinical and clinical degree of the atherogenic state on a continuous scale. This is of relevance both from a biological perspective to investigate the etiology of the longterm process of atherogenesis and in the context of primary prevention. In 2000, we established a large prospective cohort study in Bangladesh to evaluate the health effects of arsenic exposure. Here, we hypothesized that the association between arsenic exposure and IMT differs by gut microbiome composition, and individuals with certain gut microbiome composition may be more susceptible to effects of arsenic exposure on IMT. We conducted a study of 250 participants recruited from the parent cohort who underwent IMT measurement and provided fecal samples for 16S rRNA gene sequencing. We assessed the association between arsenic exposure and the gut microbiome, between the gut microbiome and IMT, as well as the interaction between arsenic exposure and the gut microbiome in IMT.

Materials and methods

Study participants

The parent study, the Health Effects of Arsenic Longitudinal Study (HEALS), is an ongoing, prospective cohort study in Araihazar, Bangladesh. The selection of cohort participants, study design, and methodologies have been described in detail elsewhere (Ahsan et al., 2006). In short, between October 2000 and May 2002, we recruited 11,746 married adults (original cohort) aged 18–75 years who were residents of the study area for at least 5 years and primarily consumed drinking water from a local tube well. At baseline, trained physicians conducted in-person interviews using a standardized questionnaire and clinical examinations, and collected spot urine samples from participants in their homes using structured protocols. Subsequently, an additional 8287 participants were recruited between 2006 and 2008 as the expansion cohort following the same methodologies. Overall participation rate was 97%. The cohort is being actively followed up biennially with similar in-person visits (Ahsan et al., 2006). A field clinic was set up exclusively for the cohort participants to receive medical diagnoses and treatments and facilitate the follow-up (Ahsan et al., 2006) and ancillary studies. The study procedures were approved by the Ethical Committee of the Bangladesh Medical Research Council and the Institutional Review Boards of Columbia University and the University of Chicago. Informed consent was obtained from all the study participants.

For the present study, recruitment, fecal sample collection, and IMT measurements were conducted between February 2015 and November 2016. We randomly selected 400 participants aged 25–50 years, free from diabetes and CVD, and lived in the 6 nearby villages surrounding the clinic (15-minute by walking or rickshaw to the clinic). These participants did not differ significantly from the remaining HEALS participants regarding demographics and lifestyle factors (data not shown). A total of 300 participants were recruited among 328 participants who fulfilled the eligibility requirements including absence of antibiotics use in past month and willing to come to the clinic to provide stool samples and undergo IMT measurement. Out of the 300 participants, 250 provided fecal samples.

Measurement of carotid IMT

We have extensive experience with IMT measurements as we have investigated a set of risk factors of atherosclerosis in the HEALS (Chen et al., 2013a; McClintock et al., 2014; Wu et al., 2016; Wu et al., 2014). For the present study, the IMT measurements were taken at the field clinic following the same protocols as previously described (Chen et al., 2013a). We confined IMT measurements to the common carotid artery (CCA), as recommended (Stein et al., 2008). Briefly, the IMT was measured in the near and far walls of the CCA at both sides of the neck using a SonoSite MicroMaxx ultrasound machine (SonoSite, Inc., Bothell, WA) equipped with an L38e/10-5-MHz transducer by a single physician trained and certified according to a specific protocol (Desvarieux et al., 2005). The IMT measurements were analyzed offline with Matlab (Mathworks, Natick, MA). Similar to previous studies (Holewijn et al., 2009; Polak et al., 2015; Takasu et al., 2010), the primary outcome variable was defined as the mean of all four CCA IMT values.

Arsenic exposure measurements

All participants had data on water and urinary arsenic from previous visits in the parent HEALS. Prior to subject recruitment, water samples from all wells (n = 10,971) in the study area were tested for total water arsenic concentration, using high-resolution inductively-coupled plasma mass spectrometry with a detection limit of < 0.2 μg/L. Information on duration of well use and source of exposure was also collected at baseline. All study participants used the tested wells as their exclusive source of drinking water, and well arsenic level does not fluctuate appreciably over time (Ahsan et al., 2006). At each follow-up, we collected data on current well use and whether the participants switched to a different well. We derived a time-weighted well arsenic concentration (TWA) as a function of drinking durations and well arsenic concentrations (TWA, μg/L) =Σ C/ ΣT, where C denotes the well arsenic concentration at ith visit and T denotes the drinking duration for respective well. Spot urine samples were collected at baseline and at all follow-up visits. Urinary total arsenic concentration was measured by graphite furnace atomic absorption, using the Analyst 600 graphite furnace system (Perkin-Elmer, Waltham, MA) with a detection limit of 2 μg/L (Nixon et al., 1991). Urinary creatinine was measured by a colorimetric method based on the Jaffé reaction (Slot, 1965), and creatinine-adjusted urinary total arsenic concentration was subsequently expressed as micrograms per gram creatinine. We derived a time-weighted urinary arsenic concentration using the same procedure as that of well arsenic concentration.

Fecal specimen collection

Because the participants were not accustomed to collect, transport, and store biological samples, especially fecal samples, on their own, we asked the participants to come to the clinic for sample collection within 2 days post IMT measurement. When participants arrived at the clinic, a senior research officer provided an empty ThermoFisher Scientific vial (Waltham, MA) to each participant for fecal collection. Stool from a single bowel movement was collected. The samples were stored in −20 °C freezer immediately and kept frozen at the clinic until processing at the field laboratory for DNA extraction.

DNA extraction and 16S rRNA gene sequencing

Total DNA was extracted from fecal samples using the MOBIO PowerSoil DNA Isolation Kit (MO BIO Laboratories, Carlsbad, CA) at the field laboratory. Polymerase chain reaction (PCR) amplification and sequencing were performed as described previously (Yang et al., 2009). Briefly, the 16S rRNA gene was amplified targeting the hypervariable V3-V4 region using universal primer set 347F/803R. PCR products were purified using Agencourt AMPure XP (Beckman Coulter Life Sciences, IN) and quantified using the Agilent 4200 TapeStation (Agilent Technologies, CA). Amplicon libraries were pooled at equimolar concentrations and sequenced using the Illumina MiSeq 600-cycle (2 × 300 bp) reagent kit (Illumina, Inc., San Diego, CA). One sample failed PCR amplification and was not further sequenced.

Bioinformatics and quality control

Sequencing data were processed using QIIME (Quantitative Insights Into Microbial Ecology) 1.8.0 (Caporaso et al., 2010). Sequencing reads were demultiplexed and low-quality reads with a quality score of < 25 were filtered out. Chimeric sequences were removed using Chimer-aSlayer (Haas et al., 2011). The pre-processed sequences were then clustered into operational taxonomic units (OTUs) at 97% identity using UCLUST (Edgar, 2010) against the Greengenes database 13.8. The most abundant sequence in each OTU was selected as a representative and taxonomy was then assigned for each OTU from phylum to genus levels using the Greengenes database 13.8 and Ribosomal Database Project (RDP) Classifier 2.2 (Wang et al., 2007). Several alpha diversity metrics at the OTU level (Lozupone and Knight, 2008; McMurdie and Holmes, 2013) were estimated based on rarefied sequence count (10,000 sequences per sample) using QIIME.

Statistical analyses

We first calculated descriptive statistics for demographic, lifestyle, and arsenic exposure variables, as well as alpha diversity metrics, mean and standard deviation for continuous variables and distribution (%) for categorical variables, by tertiles of IMT levels. We used Chi-square and Kruskal-Wallis tests to detect group differences in categorical and continuous variables, respectively. Microbial composition was expressed as relative abundance of bacterial taxa at the levels from phylum to genus by dividing the number of reads from the respective taxa by the total number of reads from all taxa within each individual. We calculated nonparametric Spearman's partial correlations of alpha diversity metrics with time-weighted water arsenic or urinary arsenic and IMT, as well as of individual bacterial genera with age (years), body mass index (BMI, kg/m2), education (years), systolic blood pressure (SBP), and waist-to-hip ratio (WHR). Heatmaps were generated to illustrate the latter correlations using the heatmap.2 function from the gplots package in R (version 3.4.0; R Core Team, Vienna, Austria). Multivariable linear regression models were used to examine the associations of sex and smoking status with individual bacterial genera adjusting for sex, age, BMI, smoking status (never and ever), and education. All P values were corrected for multiple testing using the Benjamini-Hochberg procedure (Benjamini and Hochberg, 1995) to control false discovery rate (FDR) ≤ 5%.

We assessed whether overall microbial composition (β diversity) differed by arsenic exposure or IMT using the microbiome regression-based kernel association test (MiRKAT) (Zhao et al., 2015). MiRKAT can be used to test the association between a microbial community and a continuous/binary phenotype via a kernel metric, where the kernel can be constructed using phylogenetic or non-phylogenetic distance matrices. Unweighted and weighted UniFrac distance (Lozupone and Knight, 2005; Lozupone et al., 2007) and Bray-Curtis dissimilarity (Bray and Curtis, 1957) were used to measure the microbial composition profiles for the MiRKAT test. Since each distance metric has strength in capturing a specific association pattern, we also conducted an omnibus test to simultaneously consider the aforementioned distance matrices (Chen et al., 2018). We conducted multivariable linear regression analyses to test the associations between arsenic exposure, defined by time-weighted water arsenic or urinary arsenic (independent variable), and each of the taxa (dependent variable), adjusting for sex, age, BMI, smoking status, and education. We also used multivariable linear regression models to test differences in IMT (dependent variable) in relation to each of the taxa (independent variable), adjusting for betel quid chewing (betel quid is an addictive substance that is wrapped in a piper-betel leaf and chewed widely in South Asian populations and has been associated with IMT in our population (McClintock et al., 2014); never and ever) and SBP in addition to the covariates mentioned above. All P values were adjusted for multiple testing per taxonomic level. We further performed taxa selection using the least absolute shrinkage and selection operator (LASSO) linear regression (Tibshirani, 1996) implemented in R ‘glmnet’ package. This method uses an L1 penalty to shrink the regression coefficients for unimportant taxa to zero, resulting in a parsimonious model where only the taxa with the strongest associations with the outcome will be selected. Taxa were selected from the LASSO linear regression models with the optimal value of lambda from leave-one-out cross-validation at each taxonomic level. Covariates included in the traditional linear regression model were controlled in the LASSO taxa selection process.

For bacterial taxa that were significantly associated with IMT or arsenic exposure after controlling for multiple testing, we tested their interactions with arsenic exposure in IMT. The relative abundance of each bacterial taxon was tested separately in the interaction models: Y = α + β + β + β * F, where Y is IMT, As is time-weighted water arsenic or urinary arsenic as a continuous variable, F is the relative abundance of a given bacterial taxon, and the term As * F denotes the cross-product of arsenic exposure and the relative abundance of the taxon. We estimated regression coefficients and their 95% confidence intervals (CIs) for the difference in IMT in relation to 1) a 1-standard deviation (SD) increase in time-weighted water arsenic or urinary arsenic (βA), in the absence of bacterial taxa, 2) a 1% increase in the relative abundance of bacterial taxa in the absence of arsenic exposure (βF), and 3) synergistic effect of the joint presence of a 1-SD increase in arsenic exposure and a 1% increase in the relative abundance of bacterial taxa beyond the sum of their individual effects (βAF). The significance of the interaction on the additive level was judged based on the P value for coefficient associated with the cross-product (βAF). Potential confounders included sex, age, BMI, smoking status, education, betel quid chewing, and SBP. All statistical analyses were conducted using SAS (version 9.4; SAS Institute, Inc., Cary, NC).

Results

Sequencing summary

We obtained ~8.6 million quality-filtered sequencing reads from the 249 samples, with an average of 34,520 reads per sample. These reads were clustered into 108,987 OTUs, with a mean of 3297 (SD 1295) OTUs and a range of 38–8909 OTUs per sample (Supplemental Table 1). The OTUs were assigned to 19 phyla, 39 classes, 78 orders, 142 families, and 273 genera. After removing taxa without a name assigned at each taxonomic level and the very rare taxa (relative abundance < 0.01%), a total of 8 phyla, 16 classes, 24 orders, 39 families, and 54 genera were included in final analyses. A list of bacterial taxa with a mean relative abundance of > 1% at the levels of phylum and genus was presented in Supplemental Table 2.

Characteristics of the study subjects

The study population consists of 59% women and 41% men, with a mean age of 48.6 (SD 7.9) years (Table 1). The participants were lean with an average BMI of 21.5 kg/m2 and had an average of 2.4 years of formal education. The proportion of ever smokers was much higher in men (76.5%) than in women (13.6%). The prevalence of betel quid chewing in the overall population was 50.6%. Carotid IMT was measured on average 11.6 years after baseline. The mean and median of IMT was 819.5 μm and 804.2 μm, respectively, with a range of 653.9 μm to 1206.7 μm. The participants were exposed to moderate levels of well water arsenic (78.2 μg/L), with mean urinary arsenic of 273.1 μg/g creatinine. Male sex, increasing age, ever smoking among both men and women, betel quid chewing, elevated blood pressure, and increasing WHR were associated with higher levels of carotid IMT. There was no association of BMI and educational attainment with IMT. There was also no evidence of an association of arsenic exposure and alpha diversity metrics with IMT (Table 1).

Table 1

Distribution of population characteristics and alpha diversity metrics by tertiles of carotid IMT[a].

	Overall(n = 249)	Tertiles of carotid IMT[b]			P-value[c]
		1 (n = 83)	2 (n = 83)	3 (n = 83)
IMT (μm)	819.5 ± 98.5	730.5 ± 24.3	801.5 ± 18.3	926.6 ± 91.9	–
Male (%)	41.0	27.7	41.0	54.2	< 0.01
Age (years)	48.6 ± 7.9	43.9 ± 7.2	49.0 ± 7.7	52.9 ± 6.2	< 0.01
Body mass index (kg/m2)	21.5 ± 4.1	21.4 ± 4.2	21.6 ± 4.2	21.6 ± 3.9	0.81
Waist-to-hip ratio	0.91 ± 0.08	0.90 ± 0.07	0.90 ± 0.07	0.94 ± 0.08	< 0.01
Ever smoking (%)
Men	76.5	65.2	73.5	84.4	0.07
Women	13.6	6.7	14.3	23.7	0.02
Education (years)[d]	2.4 ± 3.3	2.6 ± 3.3	2.2 ± 3.3	2.2 ± 3.4	0.52
Betel quid chewing (%)	50.6	34.9	57.8	59.0	< 0.01
Systolic blood pressure (mm Hg)	112.0 ± 17.4	105.5 ± 14.6	110.9 ± 16.7	119.5 ± 18.1	< 0.01
Diastolic blood pressure (mm Hg)	71.0 ± 11.5	67.4 ± 10.6	71.2 ± 9.8	74.5 ± 12.8	< 0.01
Time-weighted well water arsenic (μg/L)	78.2 ± 75.7	77.6 ± 79.9	85.6 ± 76.3	71.4 ± 70.7	0.46
Time-weighted urinary arsenic (μg/g creatinine)	273.1 ± 193.8	268.6 ± 184.9	303.1 ± 205.7	247.6 ± 188.4	0.08
Time between baseline and IMT measurement (years)	11.6 ± 2.8	11.2 ± 2.9	11.6 ± 2.7	11.9 ± 2.7	0.44
Observed OTUs	3297 ± 1295	3416 ± 1389	3237 ± 1376	3237 ± 1108	0.67
Chao1 richness estimator	9943 ± 3587	10,404 ± 3883	9695 ± 3663	9729 ± 3178	0.42
ACE	10,889 ± 3928	11,422 ± 4301	10,642 ± 4015	10,601 ± 3409	0.43
Shannon diversity index	4.5 ± 0.5	4.6 ± 0.6	4.4 ± 0.5	4.5 ± 0.6	0.15
Simpson's index	0.92 ± 0.04	0.93 ± 0.04	0.92 ± 0.03	0.92 ± 0.04	0.25
PD_whole_tree	219.0 ± 71.4	226.2 ± 75.8	215.5 ± 76.8	215.4 ± 60.8	0.60

Data were missing on body mass index for 1 subject.

Tertile 1: 653.9–767.2 μm; tertile 2: 767.3–833.1 μm; tertile 3: 833.2–1206.7 μm.

P values were computed with the Chi-square or Kruskal-Wallis test.

Assessed at baseline.

Associations of microbial composition with population characteristics

The Spearman correlation showed that each characteristic, including age, BMI, education, SBP, and WHR was associated with individual genera in a distinct pattern (Fig. 1). Within the Firmicutes phylum, BMI was positively associated with the relative abundance of the genus Acidaminococcus from the family Veillonellaceae and inversely related to the relative abundance of genus Oscillospira from the family Ruminococcaceae (adjusted P-value = 0.0003 and 0.005, respectively). Acidaminococcus and Oscillospira were present among 55.4% and 98.8% participants, respectively. We did not detect any other significant association between the individual genera and the characteristics. Multiple linear regression showed that the relative abundance of genus Oscillospira was significantly higher in women (adjusted P-value = 0.04) than in men (Table 2).

Fig. 1.

Heatmap of Spearman's partial correlations summarizing the associations of selected demographic, lifestyle, and anthropometric characteristics (x-axis) with relative abundance of individual bacterial genera (y-axis). The coefficients were adjusted for sex, age, BMI, smoking, and education. Color-coding adjacent to the vertical dendrogram corresponds to phylum membership of the taxa for that row of the heatmap. BMI = body mass index, SBP = systolic blood pressure, WHR = waist-to-hip ratio.

Table 2

Associations of sex and smoking status with relative abundance of selected genera.

Taxonomy	β[a] (95% CI; relativeabundance, %)	P-value	Adjusted P-value[b]
Sex
Genus Oscillospira	0.13 (0.05, 0.20)	0.0007	0.04
Genus Desulfovibrio	0.06 (0.02, 0.10)	0.004	0.11
Genus Prevotella	−6.8 (−12.5, −1.2)	0.02	0.33
Smoking
Genus Catenibacterium	0.58 (0.14, 1.02)	0.01	0.21
Genus Slackia	0.01 (< 0.01, 0.02)	0.01	0.21
Genus Collinsella	0.25 (0.06, 0.44)	0.01	0.21
Genus Dialister	0.70 (0.05, 1.34)	0.03	0.35
Genus [Eubacterium]	0.28 (0.01, 0.56)	0.05	0.35
Genus Peptococcus	0.01 (< 0.01, 0.03)	0.05	0.35

Coefficient from linear regression model indicates difference in relative abundance of genera in women relative to men and in ever smokers relative to never smokers, adjusting for sex, age, BMI, smoking, and education (sex or smoking was not adjusted for in its respective main analyses).

P values were corrected for multiple testing using the Benjamini-Hochberg procedure to control false discovery rate (FDR) ≤ 5%.

Associations between microbial composition and arsenic exposure

Although there was a tendency of an inverse association between time-weighted urinary arsenic and microbial richness and evenness, none of the associations attained statistical significance (Supplemental Table 3). The associations between arsenic exposure and overall microbial diversity were not statistically significant (Supplemental Table 4). At the nominal level, time-weighted water arsenic was significantly positively associated with the relative abundance of the class RF3 and order ML615J-28 within this class in the Tenericutes phylum (Table 3); time-weighted urinary arsenic was also significantly positively associated with the class Epsilonproteobacteria and order Campylobacterales within this class in the phylum Proteobacteria, genus Anaerostipes from the family Lachnospiraceae, and genus Faecalibacterium from the family Ruminococcaceae in the Firmicutes phylum; for instance, every 1-SD increase in time-weighted urinary arsenic (193.8 μg/g creatinine) was related to a 0.26% increase in the relative abundance of the genus Faecalibacterium. However, none of the associations survived correction for multiple testing.

Table 3

Associations of time-weighted water arsenic and urinary arsenic with relative abundance of individual bacterial taxa.

Taxonomy	Mean (range) relative abundance,%	Standard deviation of relative abundance,%	β[a] (95% CI; relative abundance, %)	P-value	Adjusted P-value[b]
Time-weighted water arsenic
Class RF3	0.02 (0–1.3)	0.09	0.01 (< 0.01, 0.02)	0.04	0.60
Order ML615J-28	0.02 (0–1.3)	0.09	0.01 (< 0.01, 0.02)	0.04	0.90
Time-weighted urinary arsenic
Class Epsilonproteobacteria	0.06 (0–1.7)	0.22	0.03 (< 0.01, 0.06)	0.04	0.49
Order Campylobacterales	0.06 (0–1.7)	0.22	0.03 (< 0.01, 0.06)	0.04	0.63
Genus Anaerostipes	0.03 (0–0.4)	0.04	0.01 (< 0.01, 0.01)	0.04	0.82
Genus Faecalibacterium	2.9 (0–13.6)	1.86	0.26 (0.03, 0.50)	0.03	0.82

Coefficient from linear regression model indicates difference in relative abundance of individual bacterial taxa in relation to per 1-SD increase in time-weighted water arsenic (75.7 μg/L) or urinary arsenic (193.8 μg/g creatinine), adjusting for sex, age, BMI, smoking, and education.

P values were corrected for multiple testing using the Benjamini-Hochberg procedure to control false discovery rate (FDR) ≤ 5%.

Associations between microbial composition and IMT

Shannon diversity index and Simpson's index were inversely related to IMT; however, none of the associations reached statistical significance (Supplemental Table 3). The association between overall microbial diversity and IMT was not statistically significant (Supplemental Table 4). The relative abundance of 8 taxa were nominally associated with IMT, including the order Enterobacteriales, family Aeromonadaceae, family Enterobacteriaceae, and genus Citrobacter within the phylum Proteobacteria, and genus Butyricimonas from the family [Odoribacteraceae] of phylum Bacteroidetes (Table 4). The associations of Aeromonadaceae and Citrobacter with IMT retained significance after correction for multiple comparisons (adjusted P-value = 0.02 and 0.03, respectively). The prevalence of Aeromonadaceae and Citrobacter was 70.3% and 71.9%, respectively. Every 1% increase in the relative abundance of Aeromonadaceae and Citrobacter was related to an 18.2-μm (95% CI: 7.8, 28.5) and 97.3-μm (95% CI: 42.3, 152.3) increase in IMT, respectively. Since the average relative abundance was low for these taxa, we also expressed the results using per 1-SD increase in the relative abundance. Every 1-SD increase in the relative abundance of Aeromonadaceae (1.0%) and Citrobacter (0.19%) was related to an 18.3-μm (95% CI: 7.9, 28.7) and 18.1-μm (95% CI: 7.9, 28.3) increase in IMT, respectively.

Table 4

Associations between relative abundance of individual bacterial taxa and carotid IMT.

Taxonomy	Mean (range) relative abundance, %	β[a] (95% CI; carotid IMT, μm)	P-value	Adjusted P-value[b]
Order Enterobacteriales	3.7 (0–40.1)	1.9 (0.06, 3.7)	0.04	0.52
Family Aeromonadaceae	0.09 (0–15.8)	18.2 (7.8, 28.5)	0.0006	0.02
Family Enterobacteriaceae	3.7 (0–40.1)	1.9 (0.06, 3.7)	0.04	0.81
Genus Butyricimonas	0.02 (0–1.0)	130.6 (6.0, 255.3)	0.04	0.79
Genus Citrobacter	0.05 (0–2.4)	97.3 (42.3, 152.3)	0.0006	0.03

Coefficient from linear regression model indicates difference in carotid IMT in relation to per 1% increase in relative abundance of individual bacterial taxa, adjusting for sex, age, BMI, smoking, education, betel quid chewing, and SBP.

P values were corrected for multiple testing using the Benjamini-Hochberg procedure to control false discovery rate (FDR) ≤ 5%.

Of the 8 phyla, 16 classes, 24 orders, 39 families, and 54 genera, the LASSO selected five taxa that were related to IMT, including phylum Proteobacteria, order Enterobacteriales and RF39, family Aeromonadaceae, and genus Citrobacter. Of note, except the order RF39, all other taxa were significantly associated with IMT in the traditional linear regression model (Table 4).

Interaction between microbial composition and arsenic exposure in IMT

Since no bacterial taxa were significantly related to arsenic exposure after correcting the P-values for multiple comparisons (Table 3), we conducted interaction analyses between arsenic exposure and the bacterial taxa that had a significant main effect on IMT including Aeromonadaceae and Citrobacter, to reduce multiple testing (Table 5). Only Citrobacter had a significant interaction with time-weighted water arsenic in IMT (P-value = 0.04). IMT difference was −6.1 μm (95% CI: −17.1, 4.9) for every 1-SD difference of 75.7 μg/L in time-weighted water arsenic in the absence of Citrobacter and was −34.9 μm (95% CI: −172.3, 103.4) for every 1% increase in the relative abundance of the genus in the absence of water arsenic. There was a difference of 121.7 μm (95% CI: 4.8, 238.5) in IMT for the joint presence of a higher level of water arsenic and the genus beyond the sum of their individual effects.

Table 5

Interaction between time-weighted water arsenic and relative abundance of individual bacterial taxa in carotid IMT[a].

Taxonomy	β[b]As (95% CI; carotid IMT, μm)	β[c]t (95% CI; carotid IMT, μm.)	μ[d]int (95% CI; carotid IMT, μm)	P-value
Family Aeromonadaceae	−4.1 (−15.2, 7.1)	−82.1 (−364.5, 200.3)	60.8 (−110.1, 200.3)	0.48
Genus Citrobacter	−6.1 (−17.1, 4.9)	−34.9 (−172.3, 103.4)	121.7 (4.8, 238.5)	0.04

Linear regression models adjusted for sex, age, BMI, smoking, education, betel quid chewing, and SBP.

Coefficient indicates change in carotid IMT associated with per 1-SD increase in time-weighted water arsenic.

Coefficient indicates change in carotid IMT associated with per 1% increase in relative abundance of individual bacterial taxa.

Coefficient indicates change in carotid IMT associated with multiplicative interaction between per 1-SD increase in time-weighted water arsenic and per 1%

Discussion

Although several bacterial taxa were significantly associated with time-weighted water arsenic or urinary arsenic at the nominal level, none of the associations remained significant after correction for multiple testing. To the best of our knowledge, the present study is the first epidemiologic study to investigate the role of the gut microbiome in subclinical atherosclerosis, represented by carotid IMT. We identified a significant association of the family Aeromonadaceae and genus Citrobacter with IMT (adjusted P = 0.02 and 0.03, respectively). Importantly, these two taxa were also selected using the LASSO as important ones contributing to IMT. This is also the first study to evaluate the inter-relationships between arsenic exposure from drinking water, the gut microbiome, and IMT. We found a significant interaction between Citrobacter and time-weighted water arsenic, in IMT (P for interaction = 0.04).

We identified a very diverse gut microbial community as indicated by the considerably higher number of species-level OTUs (a total of 108,987 OTUs) than those reported in other studies, which typically were several hundreds to a few thousands of OTUs (Escobar et al., 2014; Muegge et al., 2011; Nam et al., 2011). This microbial community shared high similarity with Asian Indians and was significantly different from that of populations of developed countries (Arumugam et al., 2011; Bhute et al., 2016; Lin et al., 2013; Nam et al., 2011). For instance, the genus Prevotella was the most prevalent bacteria in the gut microbiome of both the Indian and Bangladeshi populations (including our study population) (Bhute et al., 2016; Lin et al., 2013), while the genus Bacteroides dominated the gut microbiome of American and European populations (Arumugam et al., 2011; Bhute et al., 2016). Interestingly, we found a positive association between the genus Acidaminococcus and BMI, consistent with a large Korean cohort study on obesity (Yun et al., 2017). In addition, we found an inverse association between BMI and the relative abundance of genus Oscillospira, which has been associated with leanness in both children and adults (Escobar et al., 2014; Konikoff and Gophna, 2016; Tims et al., 2013). These data indicate that despite differences in diversity and composition across different populations, there may be common gut microbial determinants of health.

Arsenic exposure, measured using time-weighted water arsenic and urinary arsenic, was positively related to the relative abundance of several bacterial taxa from the phylum Tenericutes, Proteobacteria, and Firmicutes. However, none of the associations retained significance after correction for multiple testing. These bacterial taxa have either beneficial or adverse effects on human health. For instance, the order ML615J-28 in class RF3 was related to leanness and gluten-free diet (Bonder et al., 2016; Goodrich et al., 2014) and the genus Faecalibacterium was related to anti-inflammatory activity within the gut (Sokol et al., 2008). On the other hand, the class Epsilonproteobacteria and order Campylobacterales are pathogenic microorganisms and enriched in inflammatory bowel disease (Rizzatti et al., 2017). The genus Anaerostipes was enriched in subjects with high CVD risk profile and pre-diabetes (Kelly et al., 2016). Given these evidence and that arsenic is known to have anti-microbial properties and was used to treat infectious diseases, future larger studies are warranted to elucidate the role of arsenic exposure in shaping the gut microbiome.

Several studies investigated the impact of arsenic upon the gut microbiome of model organisms. For instance, a significant decrease in the relative abundance of Firmicutes was observed in arsenic-treated mice (Lu et al., 2014). In another mice study, there was a significant time- and dose-dependent increase in the relative abundance of Bacteroidetes and a proportionate decrease in Firmicutes with increasing arsenic exposure (Dheer et al., 2015). A recent study in larval zebrafish demonstrated that arsenic exposure (10-, 50-, 100-μg/L arsenic-treated zebrafish as a pooled group) altered the developing microbiome, resulting in both increases and decreases in abundance of genera belonging to Proteobacteria and Firmicutes, as well as increases in abundance of genera belonging to Actinobacteria and Bacteroidetes (Dahan et al., 2018). Differences in exposure levels and study designs may explain differences between our studies and studies of model organisms. Our study population was exposed to arsenic at low-to-moderate levels, which may not have considerable impact on the composition of the gut microbiome. We had 0.80 power to detect an effect as small as 0.18 standardized effect (standardized to standard deviation of bacterial taxa) in the relative abundance of a given bacterial taxon associated with one SD difference in arsenic exposure variables. For instance, we had 0.80 power to detect a 0.33% difference in the relative abundance of the genus Faecalibacterium in relation to one SD difference in arsenic exposure. The observed estimates (Table 3) were smaller, and we may need a larger sample size to detect an association.

Our findings that every 1-SD increase in the relative abundance of Aeromonadaceae and Citrobacter was related to an 18-μm increase in IMT are novel. Given a strong association between IMT and CVD risk, we estimated that an 18-μm increase in IMT could translate to a 6% increased risk of coronary heart disease (Chambless et al., 1997). Although the estimation may differ by populations and background risk, the data indicate that inter-individual variability in the abundance of certain gut bacterial taxa may explain differences in CVD risk. Species of Citrobacter are regarded as opportunistic pathogens and can lead to a variety of infectious diseases in humans involving the urinary, gastrointestinal, and respiratory tracts. Citrobacteria have also been implicated in infective endocarditis (Dzeing-Ella et al., 2009; Tellez et al., 2000). Among the 11 Citrobacter species, three (Citrobacter freundii, Citrobacter koseri, Citrobacter amalonaticus) are known to be pathogenic in humans, causing most cases of Citrobacter-related infections; however, emerging studies showed that other species can also cause opportunistic infections (Hirai et al., 2016; Lai et al., 2010; Oyeka and Antony, 2017). A previous study observed a positive association of Citrobacter with weight, BMI and pro-inflammatory cytokine IL-1β (Xiao et al., 2014).

The relative abundance of the phylum Proteobacteria, order Enterobacteriales, and the family Enterobacteriaceae to which Citrobacter belongs were also nominally associated with IMT (P = 0.097, 0.04, and 0.04, respectively). The LASSO also selected Proteobacteria, Enterobacteriales, and Citrobacter as important taxa associated with IMT. Members of the Proteobacteria are facultative, anaerobic gram-negative bacteria with an outer membrane mainly composed of lipopoly-saccharides (LPS) (Allcock et al., 2001). Numerous studies endorse the concept that an increased abundance of Proteobacteria in the gut reflects an unstable structure of the gut microbial community (Shin et al., 2015) which has been associated with intestinal inflammation (Morgan et al., 2012) and metabolic disorders (Fei and Zhao, 2013; Larsen et al., 2010). Many studies support the involvement of Enterobacteriaceae in intestinal inflammation (Garrett et al., 2010; Morgan et al., 2012). Interestingly, a recent study found an increased abundance of Enterobacteriaceae in 218 patients with atherosclerotic CVD relative to 187 healthy controls (Jie et al., 2017). The family Aeromonadaceae has been reported to be more abundant in type 2 diabetes patients (Wang et al., 2017). As this family also belongs to Proteobacteria, which shares many biochemical characteristics with Enterobacteriaceae (Igbinosa et al., 2012), the LPS-induced inflammation may contribute to its association with atherosclerosis.

In a previous separate study of 959 individuals, we reported a positive association between arsenic exposure and IMT (Chen et al., 2013a). In the present study, the main effect of arsenic exposure on IMT was not significant, probably due to the limited sample size, but we still observed a significant positive interaction between Citrobacter and arsenic exposure in IMT. We had 0.80 power to detect an effect as small as 0.18 standardized effect, or 17.7 μm in IMT (SD 98.5 μm) associated with one SD difference in the relative abundance of bacterial taxa or arsenic exposure. The effect of Citrobacter was stronger than expected. Our data suggest that individuals with higher level of Citrobacter may be more susceptible to the cardiovascular effect of arsenic exposure. A higher relative abundance of Citrobacter in the gut microbiome may indicate a higher level of vascular inflammation that could render the host more vulnerable to the effects of other exogenous exposures such as arsenic exposure on endothelial dysfunction and atherosclerosis.

This is the largest study of the gut microbiome conducted in South Asians. Since this rural population lacks basic health care services from the existing health care facilities, the study participants were free from several factors, such as frequent use of medications, supplements, or antibiotics, that may disrupt the gut microbiome. Also, the focus on a homogeneous population with respect to other risk factors of IMT enhanced the internal validity of the study (Dawber, 1980; Szklo, 1998; Willett and Colditz, 1998). The fact that this cross-sectional study is not a prospective study is a potential limitation. However, it is not likely that the participation of the study was dependent on both the composition of the gut microbiome and IMT level (a necessary condition for bias). Another limitation is the lack of data on current arsenic exposure. However, the correlation of water arsenic or urinary arsenic between each follow-up was high (0.66–0.83), suggesting that arsenic exposure did not change substantially during follow-up. Finally, 16S rRNA gene sequencing does not generally provide sequencing resolution below the genus level and does not allow for direct functional profiling. Future studies that investigate the functions and underlying pathways of the gut microbiome and atherosclerosis using metagenomic sequencing are needed to shed light on the underlying mechanisms.

In summary, our study did not identify significant associations between arsenic exposure and the composition of the gut microbiome. We found a significant association of the family Aeromonadaceae and the genus Citrobacter with IMT, as well as a significant interaction between Citrobacter and water arsenic in IMT. The data suggest that the gut microbiome may play an important role in the development of atherosclerosis, especially among individuals with higher levels of arsenic exposure. Future prospective studies with metagenomic sequencing are needed to further elucidate the contribution of the gut microbiome itself and its interaction with arsenic exposure to the development of atherosclerosis and other CVD outcomes.