To investigate the dynamics of Chlamydia muridarum (C.muridarum)colonization in the GI tract, 20 mice were divided into five groups (four mice per group) and inoculated intragastrically with live C. muridarum. GI tissues were collected at five time points: day 3, day 7, day 14, day 28, and day 35 post-infection, with one group sacrificed at each time point. The GI tract was segmented into the stomach, small intestine (SI)—comprising the duodenum, jejunum, and ileum—and large intestine (LI), which includes the cecum, colon, and rectum. All gut contents were removed, and the tissues were homogenized to evaluate the infection burden.

At day 3 post-infection, Chlamydia was detected across all segments of the GI tract. Over the following weeks, the pathogen burden in the stomach and small intestine progressively declined, eventually clearing entirely on day 35, while the burden in the large intestine steadily increased. By day 35, Chlamydia was exclusively detected in the large intestine, indicating its role as the primary site for long-term colonization ( Figure 1A ).

To further study the impact of Chlamydia infection on gut epithelial cells, three additional groups of mice (three per group) were inoculated with live Chlamydia muridarum (Chlamydia), heat-inactivated Chlamydia (De_Chlamydia), or SPG (Naive control). These groups were used for subsequent RNA sequencing analysis. The infection burden was monitored by recovering Chlamydia inclusions from rectal swabs. Long-term colonization was observed only in the group inoculated with live Chlamydia, with the pathogen detectable in rectal swabs up to day 30 post-infection. No Chlamydia was recovered from the rectal swabs of mice in the heat-inactivated Chlamydia or SPG control groups, confirming the necessity of viable Chlamydia for establishing persistent colonization ( Figure 1B ).

To explore the impact of Chlamydia infection on specific cell populations within the mouse colon, we conducted a comprehensive scRNA-seq analysis of colon tissues from mice infected with live Chlamydia, heat-inactivated Chlamydia, and uninfected controls on day 35 post infection. After initial data quality control, we obtained transcriptomic profiles from a total of 59,197 cells, which we normalized and subjected to principal component analysis and UMAP clustering. This analysis identified 27 distinct cell clusters (C0-26), which were annotated using known cellular markers ( Figure 2A ; Supplementary Figure S1A ).

Cells were primarily categorized into three major compartments within the intestinal mucosa: epithelial cells, stromal cells, and immune cells ( Supplementary Figure S1B ). Epithelial cells, comprising 73.6% of the sampled cells and predominantly expressing Epcam, play crucial roles in nutrient digestion and absorption, pathogen defense, and hormone secretion (Huang et al., 2018). The underlying stromal layer’s fibroblasts, expressing collagen markers (Col1a1 or Col11a1), provide structural and nutritional support to the epithelial cells (Wang et al., 2022). Immune cells in the gut, constituting 22.4% of the population and marked by Cd45 (Ptprc), are vital for maintaining mucosal homeostasis in response to microbial infections and physical or chemical stimuli (Zhu et al., 2011). This compartment includes a variety of immune cells such as T cells, B cells, plasma cells, macrophages, and mast cells. Additionally, endothelial cells expressing Pecam1 and Cdh5 and smooth muscle cells (SMC) expressing Actg2 were identified among the stromal cells ( Figure 2B ; Supplementary Figure S1C ). Interestingly, cluster C17, characterized by its unique expression of pancreatic enzyme-related genes (e.g., Cpa1, Cpa2, Ctrb1), appeared to be derived from pancreatic acinar cells and was only found in a single sample ( Supplementary Figure S1D ), underscoring its non-intestinal origin. Consequently, this cluster was excluded from further analyses to ensure a more accurate and reliable understanding of intestinal cell characteristics.

Further analysis based on the top 50 markers for each cell type revealed significant functional enrichment for immune and inflammatory response pathways in immune cells, while epithelial cells showed activation of antimicrobial and defensive responses. Stromal cells, including endothelial cells involved in angiogenesis, fibroblasts related to extracellular matrix and wound healing, and smooth muscle cells associated with muscle contraction, indicated distinct functional engagements ( Figure 2C ).

Particular attention was paid to the expression of tight junction (TJ) related markers in epithelial cells. Post-infection, a slight decrease in the expression of certain TJ genes, notably claudins (Cldns), occludin (Ocln), and junctional proteins (Tjps), was observed ( Figure 2D ). This highlights a possible disruption in the intestinal epithelial barrier, crucial for maintaining mucosal integrity (Lee, 2015).

Moreover, the interaction among various cell populations in healthy intestinal mucosa contributes to the functionality of the gut barrier. Misinterpretation of molecular signals or improper cell interactions might lead to barrier impairment and disease progression (Di Tommaso et al., 2021). Post-infection, besides the intrinsic changes within epithelial cells, the quantity of interaction signals between epithelial and other immune cell groups also showed significant alteration ( Figure 2E ). Notably, the interactions involving mast cells and macrophages, crucial for immune surveillance and pathogen phagocytosis (Abraham and St John, 2010), were diminished. This reflects the crucial role of altered immune-epithelial interactions in maintaining intestinal homeostasis following pathogen infection.

Overall, these results provide a comprehensive panorama of the unique cellular composition and alterations in the mouse colon following Chlamydia infection, illustrating how various cell groups cooperatively regulate epithelial barrier homeostasis, innate immunity, tissue repair, and disease progression.

Intestinal epithelial cells serve as the primary barrier against pathogen invasion, with each cell lineage playing a distinct role in pathogen resistance (Iftekhar and Sigal, 2021). To further delineate this, we isolated all epithelial cells from the colon and performed a detailed clustering, identifying nine major epithelial cell types. These include absorptive enterocytes (EC), transient-amplifying cells (TA), mucus-secreting goblet cells (GC), hormone-producing enteroendocrine cells (EEC), antimicrobial peptide-secreting Paneth cells (PC), pathogen-sensing tuft cells (TC), antigen-transporting microfold-like cells (M), enterocyte progenitors (EP), and intestinal stem cells (ISC) ( Figure 3A ). ECs, as the predominant cell type in the colon, are further categorized into proximal (EC_Proximal) and distal cells (EC_Distal) based on region-specific gene expression, which reflects their functional diversity ( Figure 3B ). For instance, proximal cells are specialized in absorbing iron and nutrients, whereas distal cells predominantly absorb bile acids and vitamin B12. Notably, the colonic crypts are rich in GC and TA cells located at the base, with TA cells being rapidly proliferating progeny of stem cells that give rise to ECs or EPs, and subsequently differentiate into mature ECs, GCs, and other cell types. This proliferative and differentiative capacity is vital for maintaining the integrity and renewal of the colonic mucosal layer.

Post-Chlamydia infection, a dramatic reduction in the number of TA cells was detected, suggesting a potential compromise in the epithelial renewal and repair capabilities of the gut (Sanman et al., 2021). Conversely, an increase in EC numbers may represent a compensatory mechanism where, under infection or damage, ECs proliferate rapidly to replace damaged cells, thereby maintaining the integrity and functionality of the intestinal barrier (Chelakkot et al., 2018). Importantly, we noted an increase in GC cells following heat-inactivated Chlamydia infection but a significant decrease in non-inactivated Chlamydia infections ( Figure 3C ). Considering the primary function of GC cells in secreting mucus and antimicrobial substances to prevent pathogen invasion and clear infected cells (Pelaseyed et al., 2014), this phenomenon suggests that while the heat-inactivated Chlamydia may still act as an antigenic stimulus provoking a protective response, the virulent Chlamydia strains might directly suppress GC proliferation.

Further analysis of TJ gene expression averages across different epithelial cell types revealed a notably lower TJ score in GC cells ( Figure 3D ), suggesting that Chlamydia might initially invade GC cells, thereby gradually compromising the integrity of the intestinal barrier. Subsequent analysis of differential gene expression within GC cells identified four expression patterns associated with Chlamydia infection, with patterns 1 and 3 composed of downregulated genes, and patterns 2 and 4 of upregulated genes ( Figure 3E ). GO enrichment analysis indicated significant enrichment of the “response to bacterium” pathway in patterns 1 and 3 ( Figure 3F ), implying that these downregulated genes may play a proactive role in defending against Chlamydia infection. Focusing on antimicrobial and defense-related downregulated genes, angiogenin 4 (Ang4) emerged as the most significantly downregulated gene in the Chlamydia vs. Control comparison. Ang4, known for its antimicrobial properties and crucial role in maintaining intestinal epithelial homeostasis, appears to be strongly suppressed by Chlamydia to evade host defenses. However, its expression showed an upward trend in the De_Chlamydia group, reflecting a possible adaptive response of host cells to inactivated pathogens. Moreover, a consistent downregulation of antimicrobial genes such as Reg3b and Reg3g from the inactivated to the virulent Chlamydia infection stages indicates a weakened defensive capability of the host against pathogen invasion. Notably, the specific downregulation of the Ido1 gene in GC cells post-infection, which encodes indoleamine 2,3-dioxygenase involved in tryptophan catabolism and immune regulation (Seo and Kwon, 2023), ranked among the top three significantly downregulated genes in both the De_Chlamydia and Chlamydia groups ( Figure 3G ; Supplementary Figure S2 ), pinpointing it as a potent target following Chlamydia infection in GC cells.

The stability of transcriptional regulatory networks is crucial for maintaining normal function in host cells. When pathogens invade, these networks can flexibly adjust the expression levels of antimicrobial or bactericidal genes to meet various environmental challenges and survival pressures, often under complex and finely tuned transcriptional control (Besic et al., 2020). To investigate potential transcriptional regulatory mechanisms in epithelial cells following Chlamydia infection, we utilized SCENIC to identify specific regulons, including transcription factors (TFs) and their target genes. These regulons likely coordinate cellular physiological and pathological processes and participate in responses to Chlamydia infection. By calculating the Regulon Specificity Score (RSS), we found that different intestinal epithelial cell groups contain specific activated regulons, such as Brca1(+) in TA cells and Foxa3(+) in GC cells. Significant variations were also noted among different regional ECs, such as Dbp(+) in proximal cells and Ikzf2(+) in distal cells, highlighting the diversity and specificity of transcriptional regulation in intestinal epithelial cells ( Figure 4A ; Supplementary Figure S3A ).

Further analysis of their target genes through GO enrichment demonstrated their involvement in dissimilar biological functions, indicating not only their activity in specific cell types or regions but also their functional specialization ( Supplementary Figure S3B ). We particularly focused on GC cells, which in the Chlamydia vs. Control comparison revealed 80 TFs with significantly differential expression—68 upregulated and 12 downregulated ( Figures 3B , 4B ). In contrast, only 27 differentially expressed TFs were detected in the De_Chlamydia vs. Control comparison, suggesting that the impact of inactivated Chlamydia infection on the transcriptional regulatory network stability of GC cells is less severe than that of viable Chlamydia.

Examination of these differentially expressed TFs in the three groups of GC cells showed that downregulated TFs had more pronounced changes in AUC activity post-infection, particularly Creb3 ( Figure 4C ; Supplementary Figure S3D ). As a member of the cAMP response element-binding protein family, Creb3 plays a pivotal role in the unfolded protein response (UPR), a critical component of cellular defense mechanisms often induced by pathogenic invasion that affects protein folding. Disruption in Creb3 transcriptional regulatory mechanisms could be significantly associated with Chlamydia infection (Smith, 2018). To understand whether Creb3 is involved in regulating the expression of host cell genes related to anti-infection, GO analysis revealed that Creb3-regulated target genes are significantly enriched in responses to bacteria and interferon-beta ( Figure 4E ), suggesting Creb3’s vital regulatory role in antimicrobial and antiviral defense processes. KEGG pathway analysis further indicated that Creb3-regulated target genes are primarily associated with endoplasmic reticulum protein processing and viral infection pathways ( Supplementary Figure S3E ), emphasizing Creb3’s critical role in regulating cellular ER stress and UPR, whose mediated transcriptional regulatory network disruption might be a key pathogenic mechanism allowing Chlamydia to invade GC cells.

Additionally, based on the regulatory strength scores of target genes regulated by each TF, we found that Creb3 might strongly regulate the expression of numerous antimicrobial genes in GC cells, such as Ly6a, Lypd8l, Tap1, Ligp1, Lgals9, and specifically the downregulated Ido1 ( Figure 4D ), further showcasing Creb3’s importance in combating Chlamydia infection and host defense. Moreover, we also discovered a significant positive correlation between Creb3 and Ido1 expression in GC cells ( Supplementary Figure S3F ), which is thought to possibly constrain the specific expression of Ido1 in GC cells by inhibiting Creb3 activity post-Chlamydia infection, leading to impaired GC cell barrier function.

The integrity of intestinal epithelium plays a critical role in inflammatory diseases and intestinal infections. The maintenance of epithelial homeostasis depends not only on the coordinated actions among different epithelial cell types but also relies on the multitude of immune cells residing within the intestine (Chelakkot et al., 2018). During pathogenic invasion, epithelial cells secrete a vast array of signaling molecules to recruit immune cells for immunoregulation, a process vital for clearing pathogens and preventing the spread of infections (Larsen et al., 2020). Moreover, the crosstalk between immune and epithelial cells is crucial for regulating epithelial proliferation, differentiation, and the restoration of barrier functions, as described in Figure 2E .

We constructed interaction networks between epithelial and immune cell groups under different infection statuses ( Figure 5A ). We observed robust maintenance of appropriate interactions among various epithelial cell groups and between these cells and immune cells, demonstrating an efficient cooperative mode to sustain the integrity of the epithelial barrier. Compared to the control group, interactions among epithelial cell groups, especially with the proliferative TA cells, were notably enhanced in the Chlamydia group. This suggests that Chlamydia infection may indirectly stimulate rapid renewal and regeneration of epithelial cells. Conversely, the interactions between epithelial and immune cell groups, particularly with macrophages and mast cells, were significantly reduced. Macrophages, known for their potent phagocytic abilities to ingest and clear pathogens including Chlamydia, and mast cells, which continuously monitor for pathogen invasion and attract additional immune cells to respond, illustrate that Chlamydia infection likely impedes this cellular communication (Wang et al., 2024), making it challenging for the host to recruit macrophages or mast cells to the infection sites for pathogen clearance.

Furthermore, we specifically highlighted the interaction patterns of goblet cells (GC) during Chlamydia infection, emphasizing the differences in information flow between GC cells and other types of epithelial or immune cells ( Figure 5B ; Supplementary Figure S4A ). Post-infection, the Ceacam1-Ceacam1 signaling between GC cells and other epithelial cells was significantly weakened. Similarly, signals related to TJs (such as Ocln, Cldns) also exhibited a similar weakening, potentially indicating a substantial reduction in adhesive capacity among epithelial cells, diminishing intestinal barrier function ( Figure 5B ; Supplementary Figures S4B, C ). In the signaling pathways involving interactions with immune cells, nearly all ligand-receptor pairs were suppressed after Chlamydia infection, potentially preventing GC cells from recruiting immune cells to clear the pathogen. Notably, we observed significant suppression in the Lgals9-P4hb signaling pathway within interactions involving macrophages and mast cells ( Figure 5B ). Lgals9 plays a critical role in resisting microbial invasions and viral infections; its significant downregulation suggests a weakened intestinal defense against Chlamydia infection. Importantly, the substantial downregulation of the receptor P4hb, which encodes a protein involved in preventing endoplasmic reticulum stress (ERS) and incorrect protein folding, was notable. Correlation analysis further confirmed a strong positive relationship between P4hb expression and Creb3 ( Supplementary Figure S4D ), suggesting that P4hb might activate Creb3-mediated transcriptional regulatory networks associated with anti-infection, therefore, the Lgals9-P4hb signaling pathway in intestinal GC might have substantial potential against Chlamydia infection and immune evasion ( Figure 5C ).

High rates of cell turnover and plasticity contribute to the intestine’s natural adaptability, maintaining homeostasis and health through the differentiation of various epithelial cell lineages (Kim et al., 2024). However, the impact of Chlamydia infection on epithelial lineage differentiation and cell fate remains unclear. In this study, we delineated the lineage structure of intestinal epithelial cells and their pseudo-time changes to simulate the differentiation process. Initially, we inferred four potential fate trajectories (lineages 1-4) from multipotent progenitor cells (EP) to fully differentiated cells. Lineages 1–3 respectively showed the differentiation of TA cells into goblet cells (GC), tuft cells (TC), and enteroendocrine cells (EEC). Lineage 4 likely represents the cell fate of enterocytes (EC) ( Figure 6A ). As anticipated, post-Chlamydia infection, a noticeable reduction in TA and GC cell numbers was observed along with a pronounced slowing of differentiation in lineage 1 as pseudo-time progressed ( Figures 6B, C ).

We further examined gene expression changes along lineage 1’s trajectory, clustering them into four gene expression modules (Modules 1-4), where Modules 1 and 5 genes exhibited strong expression in the later stages of pseudo-time, crucial for determining the final fate of GC cells. GO analysis revealed that Module 1 genes were enriched in processes related to protein glycosylation, apoptosis, and endoplasmic reticulum stress (ERS) ( Figure 6D ; Supplementary Figure S5A ), highlighting ERS’s significant role in mediating the differentiation process of the GC lineage. Module 5 genes predominantly represented GC-specific genes, such as the antimicrobial peptide Ang4 and the mucin Muc2, which together form a robust mucus barrier.

Finally, we focused on the pseudo-time expression changes of P4hb and Creb3, which are involved in ERS within GC ( Figure 6E ). Results showed that P4hb exhibited strong upregulation in the later stages of the GC lineage, closely associated with GC cell fate. However, in the Chlamydia group, P4hb showed significantly lower expression, suggesting that Chlamydia infection might suppress P4hb, significantly inhibiting GC cell lineage differentiation and leading to a reduction in the number of intestinal GC cells ( Figure 6D ). In contrast, the pseudo-time changes of Creb3 were relatively minor ( Supplementary Figure S5B ), likely more related to the homeostasis and function of fully differentiated GC cells.

C57BL/6J female mice, aged 5–6 weeks, were sourced from Vital River in Beijing and maintained under controlled environmental conditions (temperature: 22°C; light/dark cycle: 12 hours each). Following a one-week acclimatization period in the laboratory, these mice were used for experiments that received approval from the Ethics Committee of the Institute of the 3^rd Xiangya Hospital, Central South University, adhering to the guidelines set by the Chinese Council on Animal Care.

The C. muridarum strains used in this research were derived from the Nigg3 strain (GenBank accession number CP009760.1). We cultured Chlamydia in HeLa cells, followed by purification of elementary bodies (EBs) using established protocols, in short, Chlamydial organisms are amplified in cell culture and purified as EBs using a gradient centrifugation method, which involves preparing a density gradient, sonicating infected cells to release EBs, centrifuging to remove cell debris, loading supernatant onto the gradient, ultracentrifuging to isolate EBs, and then resuspending the purified EBs in SPG buffer for storage at -80°C. These EBs were aliquoted and stored at -80°C. For deactivation, EBs underwent a 30-minute heat treatment at 56°C.

Mice aged 6 to 7 weeks received intragastric inoculations of purified C. muridarum EBs (either live or heat deactivated), with each mouse receiving 2 × 10^5 inclusion-forming units (IFUs). The naïve control group was administered a sucrose-phosphate-glutamic acid (SPG) vehicle.

To monitor the shedding of live organisms, anorectal swabs were collected different time points after inoculation. Each swab was immersed in 0.5 ml of SPG and vortexed with glass beads to extract Chlamydial organisms, which were then titrated on HeLa cell monolayers in duplicate. The infected cultures were processed for immunofluorescence assays as described below. Inclusions were counted in five random fields per coverslip under a fluorescence microscope. For coverslips with fewer than one IFU per field, the entire coverslip was counted. Coverslips exhibiting significant cytotoxicity in HeLa cells were excluded from analysis. The total number of IFUs per swab was calculated based on the average IFUs per view, the area ratio of the view to that of the well, dilution factors, and inoculation volumes. If applicable, mean IFUs per swab were derived from serially diluted and duplicate samples. The total count of IFUs per swab was converted to log10 and used to compute the mean and standard deviation for the mice within each group at each time point, as previously published (Zhou et al., 2021).

For quantitating live organisms from mouse GI segments, each organ or tissue segment was transferred to a tube containing 0.5 to 5 ml SPG depending the sizes of the organs (Zhou et al., 2021). GI tract tissues include esophagus, stomach, small intestine (duodenum, jejunum, and ileum), and large intestine (cecum, colon, and rectum). The tissue segments were homogenized in cold SPG using a 2-ml tissue grinder or an automatic homogenizer. The homogenates were briefly sonicated and spun at 3,000 rpm for 5 min to pellet remaining large debris. After sonication, the supernatants were titrated for live C. muridarum organisms on HeLa cells as described above. The results were expressed as log10 IFU per tissue segment.

The scRNA-seq data included 59,197 cells from 3 Control samples, 3 De_Chlamydia-treated samples and 3 Chlamydia-treated mouse large intestine tissue samples. The UMI count matrix was converted into a Seurat object by the R package Seurat (version 4.3.0) (Butler et al., 2018). Cells with UMI numbers <1000 or with detected genes < 500 or with over 20% mitochondrial-derived UMI counts were considered low-quality cells and were removed. Genes detected in less than 5 cells were removed for downstream analyses.

After quality control, the UMI count matrix was log normalized. Then top 2000 variable genes were used to create potential Anchors with FindIntegrationAnchors function of Seurat. Subsequently, IntegrateData function was used to integrate data. To reduce the dimensionality of the scRNA-Seq dataset, principal component analysis (PCA) was performed on an integrated data matrix. With Elbowplot function of Seurat, top 50 PCs were used to perform the downstream analysis. The main cell clusters were identified with the FindClusters function offered by Seurat, with resolution set as default (res = 0.4). Finally, cells were clustered into 10 major cell types. And then they were visualized with UMAP plots. To identify the cell type for each cluster, we detected gene markers for each cell clusters using the “FindAllMarkers” function in Seurat package (v4.3.0) on a natural log scale was at least 0.5 and the difference of percent of detected cells was at least 0.25 and adjusted p-value was less than 0.05, then we annotated cell types using ScType tools (Ianevski et al., 2022). Function enrichment analysis of differentially expressed genes(DEGs) was performed using “clusterProfiler (4.6.2)” package.

Differentially expressed genes (DEGs) were determined with the FindMarkers/FindAllMarkers function from the Seurat package (one-tailed Wilcoxon rank sum test, pvalues adjusted for multiple testing using the Bonferroni correction). For computing DEGs, all genes were probed that the expression difference on a natural log scale was at least 0.5 and the difference of percent of detected cells was at least 0.15 and adjusted p-value was less than 0.05.

The modules of TFs were identified by the SCENIC (Aibar et al., 2017) python workflow (version 0.11.2) using default parameters (http://scenic.aertslab.org). A human/mouse TF gene list was used from the resources of pySCENIC (https://github.com/aertslab/pySCENIC/tree/master/resources). Activated TFs were identified in the AUC matrix, and differentially activated TFs were selected using R package limma (Ritchie et al., 2015) based on the fold change (logFC ≥ 0.5 or ≤ -0.5) and false discovery rate (FDR ≤ 0.05). To identify cluster-specific regulons (especially for analyses with many cell types, where some regulons are common to multiple of them) we used the Regulon Specificity Score (RSS) (Suo et al., 2018).

Cell–cell interactions based on the expression of known ligand–receptor pairs in different cell types were inferred using CellChat (Jin et al., 2021) (v2.1.1). To identify potential cell–cell communication networks perturbed or induced in mice large intestine, we followed the official workflow and loaded the normalized counts into CellChat and applied the preprocessing functions identifyOverExpressedGenes, identifyOverExpressedInteractions and projectData with standard parameters set. As database, we selected the Secreted signaling pathways and used the precompiled protein–protein-interactions as a priori network information. For the main analyses the core functions computeCommunProb, computeCommunProbPathway and aggregateNet were applied using standard parameters and fixed randomization seeds. Finally, to determine the senders and receivers in the network, the function netAnalysis_signallingRole was applied on the netP data slot.

To sort out functional categories of genes, Gene Ontology (GO) terms and KEGG pathways were identified using KOBAS 2.0 (Xie et al., 2011). Hypergeometric test and Benjamini-Hochberg FDR controlling procedure were used to define the enrichment of each term.

We used Slingshot (Street et al., 2018) to infer developmental differentiation trajectories in the scRNA-seq dataset. After dimensionality reduction and clustering, Slingshot can serve as a component in an analysis pipeline by identifying the global lineage structure with a cluster-based minimum spanning tree, and fitting simultaneous principal curves to describe each lineage, and inferring pseudo-time variables.

We performed time-dependent trend analysis of gene expression using the Mfuzz package of R language (Kumar and EF, 2007). First, the average expression level log2 (TPM/10 + 1) of each gene in each stage among single cells was calculated. Then, the ‘timeclust’ function was used to cluster different expression patterns.

Tight junction scores were calculated using the Seurat function AddModuleScore, which analyses the average expression per cell. For statistical significance testing of two independent groups, an unpaired two-tailed Student’s t-test was used. P-values of < 0.05 were considered statistically significant.