Male ApoE^-⁣/- (8 weeks old; 20–24 g) and wildtype C57BL/6 (8 weeks old, 22–25 g) mice were purchased from SiPeiFu BioTech Co., Ltd. (Beijing, China) and housed in an environment with a 12-h light cycle, controlled temperature (21–24 ℃), and free access to water and food. After ApoE^-⁣/- mice were fed an HFD (BEIJING KEAO XIELI Co., Ltd., # D12079B, Beijing, China) containing 21% fat, 50% carbohydrate, and 20% protein for 14 or 18 weeks, a modified mixture of short-term stimuli at the end of the experiment was used to activate unstable plaques, as reported by Wang et al. [19]. The mice undergoing short-term co-stimulation were divided into Model B-14wks (n = 12) and Model B-18wks (n = 12) groups. Specifically, after 14 or 18 weeks of HFD, mice in Model B groups were treated with short-term stimuli, including an iced water bath (0 °C) for 5–10 min on day 1 and 2, intraperitoneal injection with 1 µg/kg/day of LPS on day 1 and 2, and then intraperitoneal injection with 8 µg/kg of phenylephrine on day 3 (Fig. 1A). The ApoE^-⁣/- mice fed with HFD without co-stimulation were divided into “Model A-14wks” (n = 10) and “Model A-18wks” (n = 9). C57BL/6 mice (n = 10) served as the control group and were fed a standard chow diet (3% calories from fat). Body weight and food intake were recorded weekly for all mice.

Fig. 1.

Body weight and biochemical parameters of mice. (A) Experimental procedure. (B) Body weight of mice at 14 and 18 weeks. (C) Serum interleukin (IL)-1$\beta$ levels at 18 weeks. (D) Total cholesterol (TC), triglyceride (TG), low-density lipoprotein (LDL), and high-density lipoprotein (HDL) levels at 14 and 18 weeks. ^#⁢#p 0.01 vs. the control group at the same period. **p 0.01 vs. the corresponding model group at 14 weeks. ^&p 0.05 vs. Model A at the same period. ^&&p 0.01 vs. Model A at the same period.

After 12 h of fasting, the mice were sacrificed under 5% isoflurane inhalation anesthesia and blood samples, the heart’s aortic root, the entire aorta, BAT on the scapular surface, and epididymal white adipose tissue (eWAT) were harvested for follow-up analysis.

Serum total triglyceride (TG), total cholesterol (TC), low-density lipoprotein cholesterol (LDL), and high-density lipoprotein cholesterol (HDL) levels were measured using commercially accessible kits (Nanjing Jiancheng Bioengineering Institute, A110-1-1, A111-1-1, A113-1-1, A112-1-1, Nanjing, Jiangsu, China). Mouse IL-1$\beta$ was determined using enzyme-linked immunosorbent assay (ELISA) kits (DAKEWE Biotech Co., Ltd., 1210122, Shenzhen, China).

Liquid-liquid-mass spectrometry (LC-MS) technology is the foundation of lipidomic research [24, 25]. Lipidomic analysis was analyzed at Novogene Technology Co., Ltd. (Beijing, China). Vortex treatment was performed on a 100 µL plasma sample (IAS) mixed with 750 µL methanol. Subsequently, at room temperature (21–24 ℃), 2.5 mL of methyl tertiary butyl ether (MTBE) was added, and the mixture was shaken for 1 h. Phase separation was induced by adding 0.625 mL of MS-grade water, followed by a 10-min incubation. The samples were then centrifuged for 10 min at 1000 g. The upper phase was collected, while the lower phase underwent a second extraction using 1 mL of a solvent mixture (MTBE:methanol:water = 10:3:2.5, v/v/v). The organic phases from both extractions were combined, dried, and stored at –80 °C for no longer than 3 months. Before analysis, the dried residue was re-dissolved in 100 µL of isopropanol. Quality control (QC) samples were prepared by pooling and equilibrating aliquots from each plasma sample. Plasma lipidomic ultra-high-performance liquid chromatography coupled to tandem mass spectrometry (UHPLC-MS/MS) analysis was conducted using a high-field quadrupole orbitrap (Q Exactive™ HF-X) mass spectrometer system (Thermo Fisher Scientific Inc., Dreieich, Germany), coupled with a Thermo Vanquish™ UHPLC system (Thermo Fisher Scientific Inc., Dreieich, Germany), equipped with a heated electrospray ionization interface (ESI) (UHPLC-Q-Exactive-MS/MS analysis). Chromatographic separation was performed on a Thermo Accucore C30 (2.1 $\times$ 150 mm, 2.6 µm) chip (Thermo Fisher Scientific Inc., Dreieich, Germany). Data were captured and analyzed using Xcalibur 3.1 (Thermo Fisher Scientific Inc., Dreieich, Germany). The column temperature was set at 40 °C, the flow rate was 0.35 mL/min, and the injection volume was 5 µL. The mobile phase consisted of phases A (acetonitrile/H₂O, 6/4) and B (isopropanol/acetonitrile, 1/9) containing 0.1% formic acid and 10 mM ammonium acetate in a 20 min gradient detection as follows: 30% B, initial; 30% B, 2 min; 43% B, 5 min; 55% B, 5.1 min; 70% B, 11 min; 99% B, 16 min; 30% B, 18.1 min [26]. Analyses were performed in both positive and negative ionic modes. The operating parameters were optimized as follows: sheath gas, 20 arbitrary units; sweep gas, 1 arbitrary unit; auxiliary gas rate, 5 for positive mode and 7 for negative mode; spray voltage: 3 kV; capillary temperature: 350 ℃; heater temperature: 400 ℃; S-Lens Radio Frequency (RF) level: 50; automatic gain control target: 10⁶; normalized collision energy: 25, 30 for positive mode and 20, 24, and 28 for negative mode; injection time: 100 ms; isolation window: 1 m/z; automatic gain control target (MS₂): 10⁵; dynamic exclusion: 15 s, and a mass range recorded at 114–1700 m/z. The raw data were imported into the Compound Discoverer 3.01 (CD3.1, Thermo Fisher Inc., Dreieich, Germany) for preliminary screening. The parameters were set: Rt deviation of 0.2 min; m/z deviation of 5 ppm; signal strength deviation of 30%; a signal-to-noise ratio of 3, minimum signal intensity of 10⁵; and peak extraction of information such as addition ions. Subsequently, the peak area was quantified, and the target ions were integrated. The molecular ion peaks and fragmentation ions predicted the molecular formula. This was then compared with Lipidmaps (https://lipidmaps.org) and Lipid Blast databases (https://fiehnlab.ucdavis.edu/projects/LipidBlast). Blank samples were used to remove background interference. Quantitative results were obtained after normalizing data.

Fecal samples for 16S rDNA gene sequencing were analyzed at Novogene Technology Co., Ltd. (Beijing, China). Community DNA fragments underwent paired-end sequencing using the Illumina HiSeq PE250 platform (NovaSeq6000, San Diego, CA, USA). The bioinformatics pipeline for microbiome analysis was implemented in QIIME 2 (https://qiime2.org). Amplification of the V3–V4 regions of the 16S rRNA gene was achieved with the specific primers F515 (5^′-CACGGTCGKCGGCGCCATT-3^′) and R806 (5^′-GGACTACHVGGGTWTCTAAT-3^′), using the gg_13 database (https://www.arb-silva.de/) and DADA2 (https://benjjneb.github.io/dada2/) as the primary denoising method. Experimental data reliability was evaluated through rarefaction and species accumulation curves. Species diversity within samples was assessed with $\alpha$-diversity indices (Chao1 and Shannon), while $\beta$-diversity measures, including Principal Coordinate Analysis (PCoA) were used to examine inter-sample variation. A series of statistical analyses (T-test, MetagenomeSeq, and linear discriminant analysis effect size [LEfSe]) were performed to determine the difference in community structure. Differential metabolites were defined as those meeting the criteria of variable importance in projection (VIP) $>$1, p 0.05, and fold change (FC) $\ge$2 or FC $\le$0.5.

Frozen adipose tissue samples were homogenized for 30 min with radioimmunoprecipitation buffer (1% phenylmethylsulfonyl fluoride (PMSF)). In addition, protein expression levels of PPAR$\gamma$ (1:1000; Proteintech, #16643-1-AP, Wuhan, Hubei, China) and glyceraldehyde 3-phosphate dehydrogenase (GAPDH) (1:2000; Service, #GB12002, Wuhan, China) were determined using western blotting methods.

Image-Pro Plus6.0 software (Media Cybernetics, Rockville, MD, USA) was used to analyze histopathological images. Statistical analyses were conducted with GraphPad Prism 9.0 (GraphPad Software, San Diego, CA, USA), and the results were expressed as mean $\pm$ SD. To assess statistical significance, a two-way analysis of variance (ANOVA) and an unpaired two-tailed Student’s t-test were utilized, with p 0.05 considered significant. Partial least squares discriminant analysis (PLS-DA) was performed using SIMCA 14.1 software (Sartorius Stedim Data Analytics, Uppsala, Sweden), and hierarchical clustering analysis (HCA) was performed using R-language heat map packages (https://cran.r-project.org). Differential metabolites were identified based on VIP $>$1, p 0.05, and a FC $\ge$2 or $\le$0.5. Pathway analyses of significantly different metabolites were performed using the MetaboAnalyst 5.0 (https://www.metaboanalyst.ca), and the pathways having impact values $>$–0.05 and –log(p) $>$ –2 were analyzed.

The weights of the mice and serum biochemical parameters were analyzed to determine the effect of the combined intervention and HFD on serum lipid levels. In contrast to the weight of mice in the control group at 14 weeks, the weight of mice in the Model B-14wks group was reduced (Fig. 1B) (p 0.01), and the serum TC and TG levels were significantly higher than those in the Model A-14wks group (Fig. 1D) (p 0.01). As expected, the serum TC, TG, and LDL levels were significantly increased in the Model A-18wks and Model B-18wks groups compared to those in the Model A-14wks and Model B-14wks groups and the HDL levels were decreased (Fig. 1D) (p 0.01).

Notably, at 18 weeks, the serum IL-1$\beta$ level was noticeably higher in the two Model groups than in the control group, while the IL-1$\beta$ level was significantly higher in the Model B-18wks group than in the Model A-18wks group (Fig. 1C) (p 0.01). The results showed that consuming a HFD for a longer period increased blood lipid and IL-1$\beta$ levels; co-stimulation did not affect the serum lipid levels but lowered the mouse body weight at 14 weeks and further exacerbated inflammation at 18 weeks, suggesting that the effect of cold stimulation in a mouse model demonstrates time-dependent characteristics. At week 14, cold stimulation effectively promoted heat production in the adipose tissue of mice, thereby increasing overall metabolism, which may have subsequently led to weight loss. However, by week 18, this heat-producing capacity of adipose tissue decreased significantly, appearing to lose its previous metabolic function. This variation suggests that adipose tissue exhibits distinct responses to cold stimulation during different developmental stages, potentially impacting weight management and energy metabolism in mice.

The results of Oil Red O staining of the en face of the whole aorta are shown in Fig. 2A. With short-term co-stimulation, the ratio of plaque area to the aorta (the positive region of Oil Red O staining) was significantly higher in the two model groups at 18 weeks than in the corresponding groups at 14 weeks (p 0.05, p 0.01), and the rise in the ratio in the Model B-18wks group was more significant than that in the Model B-14wks group (p 0.01). Oil O staining of the aortic sinus slices revealed a similar result (Fig. 2B), showing a larger ratio of Oil O staining area to plaque area in Model B mice than in Model A mice during the same period (p 0.01).

Fig. 2.

Oil red O staining of the face of the whole aorta and cross-sections of the aortic sinus. (A) Oil red O staining of the whole aorta. (B) Oil red O staining of the aortic sinus. Scale bar: 500 µm. **p 0.01, vs. the corresponding model group at 14 weeks, ^&p 0.05, ^&&p 0.01 vs. Model A at the same period, n = 6. HFD, high-fat diet.

H&E staining of the transverse aortic sinus was utilized to assess the cross-sectional vessel area, the plaque-to-vessel area ratio, and the necrotic core-to-plaque area ratio. The collagen content was evaluated using Masson’s trichrome staining. A damaged fibrous cap and erythrocyte accumulation layer indicated previous rupture and healing processes [11]. As shown in Fig. 3A,B, when ApoE^-⁣/- mice were fed the HFD and co-stimulated, a substantially greater plaque-to-vascular area ratio was observed in the Model B group than in the Model A group at the same period (p 0.05). Accordingly, the percentage of necrotic cores was the largest, or they had even ruptured in Model B-18wks mice, and the cross-sectional vessel area was significantly larger in model mice at 18 weeks than in those at 14 weeks (p 0.05) (Supplementary Fig. 1). In addition, Masson staining (Fig. 3C) showed that the collagen fibers in the Model B-14wks, Model B-18wks, and Model A-18wks groups were fewer and sparser than those in the Model A-14wks group (p 0.01).

Fig. 3.

Plaque characteristics of ApoE^-⁣/- mice after combined stimulation. (A) H&E staining of the aortic sinus of Model A and B mice fed with an HFD after 14 or 18 weeks. The box shows a substantially greater plaque-to-vascular area ratio was observed in the Model B group than in the Model A group at the same period (p 0.05). (B) Local swelling showing the necrotic cores (NC) or plaque rupture. The red arrows point to the fibrous cap. (C) Masson’s trichrome staining. (D) $\alpha$-SMA immunostaining. (E) F4/80 immunostaining. Scale bar: 200 µm (A), 100 µm (B), 200 µm (C), 100 µm (D), 100 µm (E). (F) Vulnerable plaque index = (macrophage area + lipid area)/(smooth muscle cells (SMC) area + collagen area). Model A mice were only fed an HFD. Model B mice were treated with an HFD plus combined stimulation, including ice water (0 ℃, not submerged) for 5–10 min/time/day and intraperitoneal injection of 1.0 µg/kg/day lipopolysaccharide (LPS) for 2 days, followed by 8.0 µg/kg phenylephrine on day. **p 0.01, *p 0.05 vs. the corresponding model group at 14 weeks, ^&p 0.05, ^&&p 0.01 vs. Model A at same period, n = 6.

To quantify the contents of macrophages and smooth muscle cells in aortic plaques, immunohistochemical staining was performed using F4/80 and $\alpha$-SMA as specific markers, respectively. The results revealed no change in the percentage of F4/80-positive area to plaque area between the groups and a significant decrease in the rate of $\alpha$-SMA (Fig. 3D,E). The vulnerability index of plaques in Model B-18wks mice was the highest and was 1.5-fold higher than that in Model B-14wks and Model A-18wks mice and 3.6-fold higher than that in Model A-14wks mice (p 0.05, p 0.01), as shown in Fig. 3F.

These results suggest that a prolonged HFD promotes abnormal arterial vessel dilation and the rupture of elastic fibers. Furthermore, short-term co-stimulation accelerated lipid deposition onto the vessels, forming sparser, more extensive, and vulnerable plaques in the short term. The plaque phenotype observed in Model B-18wks mice was similar to that of unstable human plaques [27], including a huge plaque size, necrotic core, thin fibrous cap, decrease in entirety collagen, and positive remodeling, that is, a significant increase in vessel area. As the HFD time increased, these characteristics became more pronounced.

In the case of an excess energy supply, WAT causes adipocyte hypertrophy (increase in volume) or hyperplasia (increase in number). During adipose tissue remodeling, hyperplasia of WAT is generally more advantageous than adipocyte hypertrophy. This is because small adipocytes exhibit increased insulin sensitivity compared with large adipocytes, whereas large adipocytes often indicate metabolic derangements [28, 29]. In this study, the adipocyte phenotype in mice was assessed (Fig. 4). With co-stimulation, the size of eWAT adipocytes in Model B-18wks mice was severely increased compared to that of Model B-14wks mice (Fig. 4A), and the BAT adipocytes in Model B mice showed widespread hypertrophy both at 14 and 18 weeks in comparison to those of Model A at the same period (Fig. 4B). In addition, the ratio of both eWAT and BAT to body weight (mg/g) of model mice at 18 weeks was higher than that of the corresponding mice at 14 weeks (p 0.01, or p 0.05) (Fig. 4C).

Fig. 4.

Combined treatment altered adipose tissue phenotype. (A) H&E staining of white epididymal tissue (eWAT), n = 4–6. (B) H&E staining of brown adipose tissue (BAT), n = 4–6. (C) Ratio (mg/g) of BAT and eWAT to body weight, n = 10. (D) Peroxisome proliferator-activated receptor $\gamma$ (PPAR$\gamma$) protein of BAT and eWAT of mice subjected to western blotting, glyceraldehyde 3-phosphate dehydrogenase (GAPDH) was used as a loading control, n = 3. ^#⁢#p 0.01, ^#p 0.05 vs. Control group at same period, *p 0.05, **p 0.01 vs. the corresponding model group at 14 weeks, ^&p 0.05, ^&&p 0.01 vs. Model A at the same period. Scale bar: 100 µm.

PPAR$\gamma$, a transcription factor, is highly expressed in adipocytes and is the critical regulator of adipogenesis. It also regulates lipid and glucose homeostasis. Several studies have demonstrated the inhibitory role of PPAR$\gamma$ in atherosclerosis [30]. The results showed that co-stimulation effectively reduced PPAR$\gamma$ expression in eWAT and BAT at 18 weeks compared to controls. In contrast, PPAR$\gamma$ expression in the eWAT of Model A-18wks mice significantly increased (Fig. 4D). These findings suggest that co-stimulation did not affect total body weight, but it was more likely to cause further adipocyte hypertrophy than hyperplasia, especially in eWAT, by attenuating PPAR$\gamma$ expression in HFD-fed mice.

Plasma samples from each group at 18 weeks were analyzed using UHPLC-Q-Exactive-MS/MS under both electrospray lonization (ESI) positive and negative ion modes. To maintain data reliability, QC samples were included in the sequence, which was run after every seven injections. The Pearson correlation coefficient between QC samples was between 0.99 and 1, indicating that the method was stable and appropriate for these samples throughout the testing process. Typical base peak ion chromatograms are shown in Supplementary Fig. 2. To minimize false positives, we further screened the data according toSkotland et al. [31]. A total of 988 lipid molecules was identified and analyzed, and these were mainly grouped into 16 subclasses, including fatty acyl (FA), triacylglycerols, PC, phosphatidylethanolamine, phosphatidylserine, phosphatidylglycerol, phosphatidic acid, phosphatidylinositol, SM, Cer, and HexCer. Detailed information about the identified lipids is presented in Supplementary Fig. 3.

Heatmaps of HCA and PLS-DA plots were used to visualize the cluster overview of all plasma samples from the Model A-18wks, Model B-18wks, and control groups having significantly different lipid compositions (Fig. 5A and Supplementary Fig. 4). The volcano plot analysis results revealed that the plasma lipid compositions of Models B and A at 18 weeks were notably altered compared to those of the control group; that is, more than 400 molecules were upregulated, and 30 molecules were downregulated in both groups (Supplementary Fig. 5). The apparent separation between Model A-18wks and Model B-18wks is shown in Fig. 5B, with good model quality (R2: 0.69 and Q2: –0.78 in positive ion modes; R2: 0.76 and Q2: –0.68 in negative ion modes, respectively; p 0.05). The above results suggest that the lipid metabolism disorders were more severe following combined stimulation with HFD than with HFD alone.

Fig. 5.

Co-stimulation altered plasma lipidomic profiles. (A) Clustering heat maps of hierarchical clustering analysis (HCA) between three groups. (B) Partial least squares discriminant analysis (PLS-DA) score plot of Model A-18wks vs. Model B-18wks group. (C) Ratio of the peak area of lipid subclasses to the plasma lipid (%). (D) Potential lipid markers associated with plaque vulnerability were assessed using the thresholds of variable importance in projection (VIP) $>$1.4 compared with the control group, and fold change (FC) $>$1.5 or FC 0.5 compared with the Model A-18wks group, area $>$10⁷, and p 0.05, n = 6. (E) Pathway analysis of significantly different metabolites, impact value $>$0.05 and –log(p) $>$ 2. (F) Several cardiovascular risk markers reported in some clinical studies were evaluated, ^#⁢#p 0.01, ^#p 0.05 vs. Control group, ^&p 0.05, ^&&p 0.01 vs. Model at the same period.

In line with recent studies [7, 8, 9, 10], the statistical significance of the plasma lipid variables was determined using PLS-DA based on the criteria of VIP $>$1, FC $>$1.5, FC 0.67, and p 0.05. The most striking changes occurred in specific lipid subclasses (Fig. 5C). Compared to the control, the ratios (%) of the peak area of Cer, SM, and HexCer to the total plasma lipids significantly increased in both model groups and PC, FA (–54%), and phosphatidylinositol (–29%) were reduced in the Model B-18wks group. Only a slight increase in the HexCer (+34%) ratio was observed in the Model B-18wks group compared with that in the Model A-18wks group (p 0.05, p 0.01), indicating that upregulation of Cer, SM, and HexCer levels was closely related to the atherosclerotic characteristics caused by HFD (Fig. 5D). In contrast, a further increase in HexCer levels may be more closely related to the vulnerability of plaques induced by co-stimulation in HFD-fed mice. Cer, SM, and HexCer belong to sphingolipid subclasses. Numerous studies have demonstrated that sphingolipids participate in atherosclerosis development and progression [32]. Ceramides have also been linked in recent research to adipocyte malfunction that underlies metabolic disorders [33]. According to the enrichment analysis, this study revealed that co-stimulation mainly impacted the sphingolipid metabolic pathway (Fig. 5E) by regulating Cer, SM, and HexCer metabolism.

To determine the potential lipid markers linked to AS plaque vulnerability, individual molecules showing significant differences were analyzed. Compared with those in the control and Model A-18wks group, 19 lipid components were significantly elevated in the Model B-18wks group (p 0.01 or p 0.05), including PC (18:0e/20:1, 16:1/18:1, 18:0e/20:3, 22:3e/18:0, 22:5e/20:5), SM (d14:1/18:0, 14:1/26:0, d18:1/22:0), Hexosylceramide-N-Dihydrosphingosine (HexCer-NDS) (d20:0/12:1, d18:0/16:0, d28:0/12:1), Hexosylceramide-N-Sphingosine (HexCer-NS) (d18:2/16:0, d18:1/22:0, d18:2/22:0, d28:1/12:1, d20:2/20:0, d22:1/12:1), and glycerophosphoethanolamine (PE) (16:0/20:0, 20:0/18:1) (p 0.01 or p 0.05), as shown in Fig. 5D. The results suggest that these lipid components may influence plaque stability.

Several clinical cohort trials [9, 13] have confirmed that multiple individual lipid components, including Cer (d18:1/16:0, d18:1/18:0, and d18:1/24:1), HexCer (d18:1, d18:2, C16–C18, C20–C26), SM-16, and PC (16:0/16:0, 16:0/18:1, 18:0/16:1), are positively associated with the risk of cardiovascular disease (CVD) outcomes, while lysophosphatidylcholine (LPC) (16:0, 18:0, 18:1, 18:2, and 20:4) is inversely correlated [34]. Additionally, Chen et al. [35] analyzed plasma lipidomics in C57BL/6J mice on a normal diet, C57BL/6J mice, and ApoE^-⁣/- mice on a 16-week HFD; after controlling for the effects of the HFD, they found that reductions in C16 sphingosine and PC (17:1/22:6) levels and increases in PC (18:2/20:4), PC (16:0/16:0), PC (18:0/16:1), SM (d18:1/24:1), SM (d16:0/28:5), and SM (d18:1/16:0) levels were positively associated with AS. Therefore, we further assessed these reported cardiovascular risk markers. Notably, 17 lipid species were consistent with the clinical studies, given that Cer-NS (d18:1/16:0, d18:1/24:1, and d18:1/24:0), HexCer-NS (d18:1/16:0), SM (d14:1/16:0 and d26:3/16:0), and PC (16:0/16:0, 16:0/18:1, 18:0/18:1, 18:1/18:1, and 18:4e/18:1) levels were significantly higher, while LPC (18:2, 18:0, 16:0, 20:4) and PC (18:2/20:4) levels were significantly lower in the two model groups than in the control group (p 0.01 or p 0.05). Moreover, only PCs (18:0/18:1) were significantly higher in the Model B-18wks group than in the Model A-18wks group (p 0.05 or p 0.01) (Fig. 5F). These observations suggested that co-stimulation further exacerbated lipid metabolism disorders in mice fed an HFD and that PC (16:0/16:0, 18:2/20:4, 18:1/18:1), HexCer (C12:1, C16:0), Cer (d18:1/16:0), SM (C16:0), and SM containing long chains of saturated fatty acids may be markers of the risk of vulnerable plaque formation. The sphingolipid metabolic pathway is crucial in the progression and formation of atherosclerotic plaques.

Using 16S rRNA sequencing of the colon contents, we examined the effects of long-term HFD and short-term co-stimulation on the intestinal microbiota of ApoE^-⁣/- mice. The results are shown in Fig. 6A, where the operational taxonomic unit rarefaction curve is seen to gradually flatten as the sequencing volume increases, indicating adequate sequencing data. Compared with the control group, the results of $\alpha$-diversity (Fig. 6B) showed increases in the richness (chao1 index) and diversity (Shannon index) of microbial communities in the Model B-18wks group (Shannon: p = 0.061). PCoA and orthogonal PLS-DA were performed to observe the $\beta$-diversity of the microbial communities (Fig. 6C,D). The results demonstrated a partial overlap between the control and Model A-18wks groups. However, the microbial communities showed clear segregation between the Model A-18wks and Model B-18wks groups, demonstrating that co-stimulation notably altered the gut microbiota composition.

Fig. 6.

Combined stimulation increases intestinal microflora disorders. (A) Operational taxonomic unit rarefaction curve. (B) Chao1 and Shannon indices reflect the $\alpha$-diversity of the gut microbiota. Data are presented as mean $\pm$ SD. (C) Principal Coordinate Analysis (PCoA) analysis reflects $\beta$-diversity. (D) Orthogonal PLS-DA analysis between the Model A and Model B groups.

Subsequently, the community composition and species abundance were examined across various taxonomic levels. At the phylum level (Fig. 7A,B), the Model B-18wks group had significantly greater levels of dangerous bacteria (p_Firmicutes) and significantly lower relative abundances of beneficial bacteria (p_Bacteroidota) compared to the Model A-18wks groups (p 0.05 or p 0.01), resulting in a considerably higher ratio of p_Firmiciutes to p_Bacteroidoa (F/B) in the Model B-18wks group than in the other two groups (p 0.05). Additionally, differences in microbiota composition were visualized using LEfSe and clustering heatmap analysis. Fig. 7C,D show the relative abundance of microbiota clustering heatmaps of the TOP10 at the family level and TOP30 at the genus level, respectively. These heatmaps reveal that the changes at the family and genus levels were more significant in the Model B-18wks group than in the Model A-18wks group. Based on similar trends in both model groups, compared with the control group, the representative differential bacteria (p 0.05) in the Model B-18wks group are shown in Fig. 7E. Of these, the abundances of f_Muribaculaceae, f_Lachnospiraceae, g_Family_XIII_UCG-001 [14, 36], which promotes the formation of short-chain fatty acids (SCFAs), and g_Muribaculum [37], which maintains intestinal homeostasis, were significantly reduced. The Model B-18wks group showed significant increases in f_Erysipelotrichaceae, g_Faecalibaculum [14], f_Peptostreptococcaceae, f_Bifidobacteriaceae [14], g_Romboutsia [38], g_Clostridium_sensu_stricto_1 [14] abundance. These harmful bacteria are positively correlated with high-fat and high-sugar diets, obesity, dyslipidemia, or enteritis. With linear discriminant analysis (LDA) $>$3.5, the LEfSe analysis results from the phylum to the species level are shown in Fig. 7F,G. Compared with the control group, both model groups showed an increased abundance of g_Faecalibaculum and a decreased abundance of beneficial bacteria such as g_Ileibacterium, s_Ileibacterium_valens [39, 40], which are involved in bile acid metabolism. Thus, co-stimulation with HFD further enhanced the abundance of harmful bacteria associated with abnormal glycolipid metabolism and lowered the abundance of beneficial bacteria promoting SCFA and comparing bile acid metabolism to HFD alone. The stability of atherosclerotic plaques may be affected by the changes in their abundance.

Fig. 7.

Combined stimulation alters the dominant microflora composition from the phylum to the species level. (A) Stacked chart of microbial flora at the phylum level. (B) Relative abundance of Firmicutes and Bacteroidota and the ratio of Firmiciutes to Bacteroidota (F/B). Data are presented as mean $\pm$ SD. (C) HCA heatmap of microbial flora at the family level. (D) HCA heatmap of microbial flora at the genus level. (E) Representative differential microflora at the family and genus levels. ^#p 0.05 vs. control group, ^#⁢#p 0.01 vs. control group, ^&&p 0.01, ^&p 0.05 vs. Model A-18wks group, n = 6. (F) Linear discriminant analysis effect size (LEfSe) from the phylum to species level between the control and Model A-18wks groups. Linear discriminant analysis (LDA) $>$3.5. (G) LEfSe from the phylum to species level between the control and Model B-18wks groups. LDA $>$3.5. The red boxes in (F,G) highlight key bacterial taxa with significant differences in both two model groups, compared with the control group.