Antibiotic resistance profiles of gut microbiota across various primate species in Guangxi

Introduction Understanding the gut microbiota and antibiotic resistance gene (ARG) profiles in non-human primates (NHPs) is crucial for evaluating their potential impact on human health and the environment. Methods In this study, we performed metagenomic analysis of 203 primate fecal samples, including nine NHP species and humans, to comprehensively characterize their gut microbiota and ARGs. Results Our study reveals the prevailing phyla in primates as Firmicutes, Bacteroidetes, Euryarchaeota, and Proteobacteria. The captive NHPs exhibited higher ARG abundance compared to their wild counterparts, with tetracycline and beta-lactam resistance genes prevailing. Notably, ARG subtypes in Trachypithecus leucocephalus (T. leucocephalus) residing in karst limestone habitats displayed a more dispersed distribution compared to other species. Interestingly, ARG profiles of NHPs clustered based on geographic location and captivity status. Co-occurrence network analysis revealed intricate correlations between ARG subtypes and bacterial taxa. Procrustes analysis unveiled a significant correlation between ARGs and microbial phylogenetic community structure. Taxonomic composition analysis further highlighted differences in microbial abundance among NHPs and humans. Discussion Our study underscores the impact of lifestyle and geographical location on NHP gut microbiota and ARGs, providing essential insights into the potential risks posed by NHPs to antibiotic resistance dissemination. This comprehensive analysis enhances our understanding of the interplay between NHPs and the gut resistome, offering a critical reference for future research on antibiotic resistance and host-microbe interactions.


Introduction
Antibiotic resistance is a significant global threat (Bolan et al., 2012;Brinkac et al., 2017;Qiao et al., 2018;Gil-Gil et al., 2021;Mustafa et al., 2021).More than 100,000 tons of antibiotics are produced worldwide each year and are extensively used in aquaculture, healthcare, and clinical laboratories, leading to widespread drug resistance in many bacteria and the emergence of ARGs (van Boeckel et al., 2015;Zhang et al., 2015;Xu et al., 2017;Sun et al., 2019;Wang et al., 2019;Qiu et al., 2020).The World Bank has warned that if we do not take strong actions, antimicrobial resistance is expected to kill 10 million people by 2050, with accumulative economic losses reaching 100 trillion dollars (He et al., 2019;Qian et al., 2021).What's worse, ARGs are widely distributed in sewage, soil, surface water, medical waste, and animal feces, and they can be transmitted to humans through agriculture products, water, or air in both captive and wild environments, which will increase the difficulty of dealing with antibiotic resistance (Chi et al., 2020;Yang et al., 2022).From the perspective of "One Health, " an effective way to accelerate harmless and resourceful treatment of ARGs is crucial (Su et al., 2023).
NHPs are genetically and evolutional relatives of humans, and research has shown that the gut microbiome of them is biologically similar to humans (Yildirim et al., 2010).It was once believed that humans evolved from apes, so NHP would be an excellent model for studying the spread and prevention of ARGs.At present, there are many kinds of NHPs lives in Guangxi, China, including Trachypithecus leucocephalus (T.leucocephalus), Trachypithecus francoisi (T.francoisi), Macaca multta (M.multta), etc., which is reliable resources for the study of non-human primates.Studies on them have focused on genome adaptation, conservation, and behavior (Huang and Li, 2005;Li et al., 2020;Liu Z. J. et al., 2020;Fan and Pedersen, 2021;Que et al., 2021), but little is about the diversity and composition of ARGs.Identification of ARGs in NHPs can be used as a proactive method for monitoring antibiotic resistance and as a control measure against infection (Schuchat et al., 2001;Roca et al., 2015;Liu et al., 2016).A study showed that seasonal variations in the gut microbiota of freerange white-faced capuchins, but their ARGs composition remained unknown (Orkin et al., 2019).The composition and abundance of ARGs carried by the gut microbiota of cynomolgus monkeys were significantly affected by dietary changes (Yan et al., 2022).Some antibiotics, such as beta-lactam, aminoglycoside, macrolides, and bacitracin in captive pandas were significantly higher than those in wild pandas (Guo et al., 2019).Nevertheless, few studies have been conducted on the gut microbiota of multiple NHPs, and their gut microbiota and composition of the ARGs they carried are unknown.Therefore, there is an urgent need to fully understand the gut microbiota of NHPs and the composition of ARGs carried by NHPs, especially in the endangered primate T. leucocephalus.
Karst limestone is one of the typical geographical features of Guangxi, and it is the habitat of T. leucocephalus.Whether karst limestone has an impact on the gut microbiota and ARGs of primates, including T. leucocephalus, remains to be further investigated.Our previous findings showed that the gut microbiota of T. leucocephalus is mainly related to heat tolerance and cellulose degradation.Xiang et al. (2020) found that the diversity, and composition of microbiota in karst water in southwest China were closely related to ARGs and antibiotic usage.Stange and Tiehm (2020) detected beta-lactams, aminoglycosides, tetracyclines, macrolides, and beta-lactams in karst lava springs in Germany, the frequency of ermB (42.1%), tet (C) (40.8%), sul2 (39.5%), and sul1 (36.8%) was higher.With increasing detection levels of ARGs detected in the environment, we speculate that the karst limestone in Guangxi has an impact on the gut microbiota and ARGs of local primates, and they play an important role in the development and human-to-human transmission of ARGs (Allen et al., 2010).
To uncover these, 150 fecal samples of 9 NHPs were collected in Guangxi.At the same time, metagenomic data of Chinese population was downloaded from the database.In this study, high-throughput sequencing and metagenomic analysis of ARGs in NHPs in Guangxi was conducted to explore the driving factors of lifestyle, geographical location and host species on ARGs, and to study the symbiotic network between ARG subtypes and related microbiota.This study will provide a comprehensive understanding of the distribution of ARGs in the NHPs gut microbiota in Guangxi and reveal their potential risks to human health.S1; Supplementary Figure S1).All the samples were collected during the spring of 2018 to mitigate potential seasonal variations in microbial content.

Sample collection and public datasets acquisition
To obtain wild fecal samples, we strategically placed food on the ground and concealed ourselves in vegetation.Upon defecation, we used latex gloves to select fresh samples based on visual cues related to sample integrity and color.For captive specimens, healthy NHPs were randomly chosen, and their fecal samples were promptly collected in an aseptic manner post-defecation.Using a spoon, feces from the center portion were collected and deposited into sterile containers, which were then transported to the laboratory in a low-temperature box accompanied by dry ice.Subsequently, samples were immediately frozen and stored at −80°C for subsequent processing.DNA was extracted within 48 h following samples delivery.The research adhered to protocols established by the China Wildlife Conservation Association and the legal requirements of China.
For comparative analysis, metagenomic data (PRJEB21612) of 53 health humans in China were downloaded from the European Nucleotide Archive Database.1These datasets were incorporated into our study to enhance the breadth of our analysis.DNA extraction was conducted from 200 mg frozen fecal samples using the Qiagen DNA extraction kit (Qiagen, Germany) according to the manufacturer's instructions.DNA concentration and purity were assessed using the Qubit ®2.0 Flurometer (Life Technologies, CA, USA) and NanoPhotometer (IMPLEN, CA, USA), respectively.DNA integrity was evaluated through 1% agarose gel electrophoresis.Metagenomic DNA was fragmented using ultrasonic Processor, yielding an average size of approximately 350 bp.A paired-end library was generated using the TruSeq DNA Sample Preparation Guide (Illumina,15,026,486 Rev.C).High-throughput paired-end sequencing was performed on the Illumina HiSeq X-ten platform (Illumina, USA) utilizing the index PE150 sequencing strategy.Subsequently, sequences with low-quality bases (Q ≤ 19) at the terminal regions were trimmed to ensure high-quality data.Reads containing substantial "N" (≥5% of the read), adapter contamination, or those falling below a specified length (>150-bp) were filtered out.Following these steps, a set of clean data possessing high quality were obtained for further analysis.

Profiling of ARGs and microorganism composition
We quantified the profiles of ARGs in the 203 metagenomic datasets using the ARGs-OAP pipeline (version 2.0) (Feng et al., 2018;Yin et al., 2018;Yang et al., 2019;Qian et al., 2021). 2 The datasets were aligned to a structured ARG (SARG) database following established protocols outlined in prior publications.Specifically, we combined the metadata file of samples with the clean reads files in a single directory within our local Linux system.The initial Perl script was executed to identify ARG-like sequences and 16S ribosomal RNA (rRNA) gene sequences, with parameters set at an alignment length of 25 amino acids, 80% sequence similarity, and an e-value of 1 × 10 −7 (Feng et al., 2018;Yin et al., 2018).Subsequently, a candidate ARG sequences file and an online metadata file were matched against the integrated structure to align with SARG reference databases, encompassing 24 types and 1,244 subtypes of ARGs, using the second Perl script within our local Linux system.ARG abundances are normalized with respect to total 16S reads (unit: copies of ARG per copy of 16S rRNA), the number of cells (unit: ARG copies per cell), and ppm (reads carrying ARGs per million reads), a unit normalized by dividing the total reads carrying ARGs by a million read counts (Zheng et al., 2022).This process yielded the ARG profile abundance at different levels, type, subtype level, and gene (reference sequence).
For taxonomic composition analysis, metagenomes data underwent taxonomic classification and relative abundance quantification (in terms of the number of cells, not just fraction of reads) through HUMAnN2 pipeline (v0.11.1) integrated with the MetaPhlAn database of marker genes mpa_v20_m200 with settings (Li et al., 2015;Gacesa et al., 2022).Differential

Procrustes analysis
We conducted correlation analyses between ARGs-ARGs and ARGs-genus taxa using 'vegan' package in R (v.4.3.1)(Yang et al., 2022).For comparisons of M 2 values with 999 permutations, Spearman's correlation coefficients were calculated using the protest function (Qian et al., 2021).The M 2 value, representing the sum of the squared distances between matched sample pairs, was computed.To obtain the p-value, the original M 2 value was compared to the simulated distribution of M 2 values.When the p-value was below 0.05, a significant correlation between antibiotic resistance and bacterial species was established (Feng et al., 2018).

Statistical and network analysis
Principal coordination analysis (PCoA) and principal component analysis (PCA) were conducted base on the ARG subtype and taxonomic abundance tables.A heat map illustrating the abundance of ARG types and subtypes was generated using Hemi1.0.Column diagram depicting the abundance of ARG types were drawn using GraphPad prism 9.0.Petal figure showing specific genes was drawn using Venn (Yang et al., 2022). 3Simpson index, Shannon index and Chao1 index (measuring within-sample diversity) were calculated based on ARG relative abundance profiles using the PAST 4.09, and visualization was accomplished using origin (2022) software.
In total, 107 opportunistic pathogens were detected across all species through comparison with a previously published pathogen list (Li et al., 2015).Among these, 10 pathogens were shared across all ten primates (Figure 2C).Escherichia coli and Bacteroides eggerthii prevailed most frequently, with Homo sapiens (59), N. coucang (55) and N. pygmaeus (45) hosting the broadest spectrum of pathogens.Clinically pertinent pathogens were identified in NHPs, including Clostridium perfringens, Klebsiella pneumoniae, and Clostridium difficile, which were also found in the migratory birds.

Comparative analysis of resistome between NPHs and humans
Previous studies have underscored the close resemblance between the gut microbiomes of NHPs and humans, suggesting the potential Our analysis revealed the existence of 16 ARGs shared between both humans and NHPs.Moreover, rifamycin resistance genes were exclusively identified in NHPs, albeit at a limited frequency.Among these genes, tetracycline stood out as the most prevalent and dominant ARGs in metagenomic resistome profiles of both humans and NHPs, followed by MLS, beta-lactam, and multidrug (Figure 4A).Diversity index analysis demonstrated that human exhibited relatively higher Simpson, Shannon, and Chao1 indices, compared to NHPs.Notably, captive Papio displayed the highest diversity indexes.PCA plots based on ARGs relative abundance indicated that nearly all NHPs samples (particularly those from T. francoisi, T. cristatus, N. pygmaeus, and M. mulatta) were closely situated to human samples  -D,F).However, some T. leucocephalus and Lemur catta samples diverged from human samples in the PCA plot (Supplementary Figures S10A,E).A considerable abundance of ARGs was identified in samples from the Wuzhou Langur Breeding and Research Center and Nanning Zoon's T. francoisi (Supplementary Figure S9).consequently, it can be inferred that lifestyles and geography exert an influence on similarities and abundance of ARGs among the NHPs.Thus, NHPs resistomes represent a potentially crucial reservoir for studying antibiotic resistance within human gut microbiota.This warrants attention as captive NHPs frequently interact with humans, raising concerns for potential public health implications.Furthermore, it is noteworthy that both wild T. leucocephalus and M. mulatta harbor numerous ARGs that align with those found in human beings.

Co-occurrence network analysis of ARGs in NHPs
Procrustes analysis was employed to elucidate correlations between ARG subtypes and bacteria taxa at the genus level.Deviation square M 2 yielded a value of 0.9374, with a corresponding p value 0.001 for ARGs and ARGs (Figure 5A).Additionally, the deviation square M 2 for ARGs and microbiota amounted to 0.9708, with p value of 0.001 (Figure 5B).
To unravel the co-occurrence patterns of ARGs within the ten NHPs, a co-occurrence network was constructed based on robust and statistically significant correlations among different NPHs found in individuals (|ρ| > 0.6, p < 0.01).The co-occurrence network diagrams for ARGs contained 12 antibiotic types encompassing 79 resistance genes (Figure 5C).Positive correlations emerged among similar ARGs indicative of the shared selective pressure imposed by specific antibiotics.For example, tetracycline ARGs tetO, tetM, tet32, tetX2, tetX and tetX3 exhibited such correlations.Similar patterns were observed for tetracycline-related ARGs vanR and vanS.Only one negative correlation was observed between different ARGs, such as bacA (bacitracin ARG) and vatB (MLS ARG).This could be attributed to competition among host bacteria or inhibition-related dynamics (Yan et al., 2022).
Illustrated in Figure 5D, the network comprised 61 nodes and 84 edges (with a modularity index of 0.602), including Bifidobacterium, Clostridium, Subdoligranulum, and Ruminococcaceae_noname.Furthermore, 38 species levels were linked to 10 resistance gene types, encompassing tetracycline, vancomycin, aminoglycoside, bacitracin, multidrug, MLS, beta-lactamase, beta-lactam, trimethoprim and chloramphenicol.Notably, Treponema exhibited a positive correlation with vatB (MLS ARG), indicating that its potential role as a host for ARGs.Concurrently, certain species displayed negative correlation with ARGs, exemplified by the Spirochaetes phylum, Treponema, and bacA (bacitracin ARG).Thus, intricate interaction between various ARGs and microbes were evident, underscoring the intimate connection between NHP gut microbes and ARGs.Besides, based on these results we can tract the hosts for specific ARG subtypes (Li et al., 2023).

Discussion
NHPs, situated proximately to humans within the animal kingdom, share a substantial degree of genetic and physiological resemblance (Li et al., 2018;Yan et al., 2022).Exploring and contrasting the gut microbiota and antibiotic resistome within captive and wild primate populations can yield valuable insights into human interaction and captivity on NHPs (Campbell et al., 2020).Consequently, a comprehensive understanding of the gut microbiota and ARGs across diverse NHPs and their distinctions from humans is of notable import.Nevertheless, the existing research in this domain remains limited.In this study, we observed a greater similarity in gut microbial composition particularly between captive NHPs (N.coucang, N. pygmaeus) and humans.Our investigation revealed 107 opportunistic pathogens and 580 distinct ARGs, unveiling active correlation between the ARGs to ARGs, as well as ARGs and microbiota at genus level, across ten primate species.
We established that Firmicutes, Bacteroidetes, Euryarchaeota and Proteobacteria emerged as the dominant phyla in NHPs.Notably, the gut microbiome clustering exhibited species affiliations within specific geographic region.Human gut microbial composition exhibited a closer affinity to that of N. coucang and N. pygmaeus compared to other NHPs.As such, we infer that NHPs captivity may foster similarities in gut microbiomes akin to humans.Previous study reconstructed over 1,000 previously uncharacterized microbial species from 6 available NHP metagenomic cohorts, almost 90% of which associated with NHPs, suggesting that captivity possibly humanizes the primate intestinal microbiome (Manara et al., 2019).However, Li et al. (2018) identified in total 1.9 million non-redundant bacterial genes from twenty cynomolgus macaques of which only 39.49% is present in the human.We realize that discernible differences persist in the gut microbiota of captive NHPs and modern humans following diverse diets (Yang et al., 2022).This observation aligns with our identification of 10 species (Ruminococcus_bromii, Roseburia_ inulinivorans, Roseburia_intestinalis, Faecalibacterium_prausnitzii, Eubacterium_rectale, Eubacterium_eligens, Bifidobacterium_longum, Bacteroides_vulgatus, Bacteroides_stercoris) that exhibited significantly greater abundance within the human gut, as determined through LEfSe analysis.Furthermore, our investigation unveiled a notable trend of tetracycline and beta-lactam resistance among captive primates.Remarkably, tetracycline resistance was also discerned in T. leucocephalus, albert at a proportion lower than that in M. mulatta.Tetracycline and beta-lactam, recognized as broad-spectrum antibiotics, find primary use as first-line treatments for non-communicable diseases resulting from bacterial infections (Fair and Tor, 2014;Hu et al., 2021).Their antibacterial spectrum encompasses Gram-positive and Gram-negative bacteria, chlamydia, mycoplasma, rickettsia, spirochetes, mycobacteria, and some protozoa.Consequently, tetracycline and beta-lactam are extensively employed in veterinary clinics for dogs, cats, cattle, sheep, pigs, turkey, and chickens.This widespread application potentially underpins the elevated prevalence of tetracycline and beta-lactam resistance among captive NHPs.Among the 17 identified tetracycline-related ARGs, tetX stands out as an enzyme modification gene, capable of inactivating antibiotic enzymes through enzymatic hydrolysis.Notable, tetQ, tetW, tetO, tet32, and tet44 genes, which are also prevalent resistance genes in pigs and cattle, manifest substantial abundance (Lim et al., 2020).For instance, tet32, initially identified in Clostridium difficileassociated human colonic anaerobe K10, finds wide distribution in bovine rumen and pig feces (Cao et al., 2020).The identification of numerous tetracycline-related ARGs underscores the necessity for heightened resistance surveillance and judicious antibiotic usage.Strategies such as alternative antibiotic selection or modifications in the living and dietary conditions of NHPs warrant consideration as control measures.
Co-occurrence network analysis offers a high-resolution depiction of correlations between ARGs and microbiota; however, it does not elucidate definitive relationships.Nonrandom symbiosis models can be used to elucidate potential hosts of ARGs within animal feces and environment samples.Essentially, this nonrandom symbiosis pattern suggests that specific microbiota may harbor distinct ARGs.In our study, we postulated that when the abundance trends of ARGs align significantly with coexisting microbiota across different species, a nonrandom co-occurrence pattern potentially exists between ARGs and microbiota.This led us to infer that Clostridium (with 11 ARG subtypes), Bifidobacterium (with 10 ARG subtypes), and T4likevirus (with 7 ARG subtypes) might emerge as primary potential hosts for tetracycline and MLS.Given that these ARG-carrying bacteria may pose substantial risks to animal health, and potentially transmit to humans, it becomes imperative to implement stringent control measures to mitigate these concerns and curtail the risk of therapeutic failure.Several potential ARGs hosts were confirmed through experimental validation, such as mdtO in Escherichia, underscoring the value of co-occurrence network analysis as a robust analytical tool capable of unearthing fresh insights into ARGs and their potential hosts within intricate environmental systems.Nevertheless, these potential symbiosis relationships necessitate further exploration through alternate methods, such as metagenomic assembly and microbiological cultivation.
Compared to human, we observed a resemblance in human gut ARGs to those present in captive NHPs, with numerous overlapping ARG detected in both populations.An investigation involving chimpanzees, gorillas, and humans with different lifestyles reported a segregation of ARGs in wild chimpanzees and gorillas, while captive individuals exhibited a clustering closer to humans (Campbell et al., 2020).Our result further revealed an elevated ARGs abundance in captive primates compared to their wild counterparts.Concurrently, primate species within the same geographic region, such as N. coucang and N. pygmaeus, exhibited similarities in their ARG profiles, whereas significant differences in ARGs abundance were evident in T. francoisi across distinct regions.Remarkably, the detected ARG in humans were also observed in NHPs, underscoring the notable impact of lifestyle and geographic location on gut microbiota-associated ARGs.Moreover, the dispersion of ARGs in T. leucocephalus, with certain samples approximating the human pattern, potentially signifies interactions with humans, possibly through contact with agricultural products.In totality, our findings suggest that lifestyle and geographical location exert substantial influence on NHPs ARGs dynamics.
Recognizing areas warranting refinement, we acknowledge the necessity for multivariate analysis to unveil the distinct contributions of lifestyle and geographic factors to the variability of ARG profiles.While our study presents a pioneering exploration of the ARG composition in T. leucocephalus and its comparison to humans, the inclusion of only one captive T. leucocephalus precludes a comprehensive analysis of the specific impact of the karst environment on ARGs.In our analyses, spearman correlation elucidated associations between ARG and genus-level microbiota.However, it is crucial to note that numerous ARGs can transfer among bacteria through self-proliferation and horizontal gene transfer in real-world settings.Therefore, certain genes, particularly those widely distributed across bacterial species, might not be solely discerned through correlation analysis (Liu et al., 2019).Consequently, further investigations are warranted to deepen our insights in subsequent studies.Regardless, our study unveils the gut microbiota and resistome of different primates by metagenomic sequencing, providing essential evidence for NHPs and humans to combat the threat of antibiotic resistance.

Conclusion
Through comprehensive metagenomic analysis of ten primates sequencing datasets, we have established an extensive and diverse repository of primate gut microbiota and ARGs.Our study underscores the potential influence of factors like captivity, geographic location, and species on the composition of gut microbiota and ARGs in NHPs.Notably, captive NHPs exhibited elevated relative abundance of ARGs compared to their wild counterparts.Moreover, our co-occurrence network analysis revealed discernible correlation between ARGs and distinct microbiota components.This research enriches our comprehension of the parallels and distinctions in gut microbiota and ARGs profiles between NHPs and humans, offering a pivotal reference for future investigation in this domain.

FIGURE 2
FIGURE 2 Comparative analysis of gut microbiome among various primate species.(A) Hierarchical classification tree of the ten primate species.(B) PCoA plot based on Bray-Curtis distance of genus-level relative abundance profile.(C) The difference of genus-level relative abundance by LEfSe (p < 0.01 and LDA > 4).(D) Petal diagram showed shared opportunistic pathogens among the ten species.

FIGURE 3
FIGURE 3 ARGs profiles in different primate species.(A) Relative abundance of major ARG types.(B) Abundance of the top 100 predominant ARG subtypes in the 203 primate fecal samples.(C) Comparative analysis of total ARG abundance across different primate species.(D) PCoA based on Bray-Curtis distance of ARG relative abundances.

FIGURE 4
FIGURE 4 Comparative analysis of gut antibiotic resistome among various primate species.(A) Comparison of ARG types and proportion in the nine NHPs and human.(B) Simpson index Comparison of ARGs in NHPs and human.(C) Shannon indexes Comparison of ARGs in NHPs and human.(D) PCA plot illustrating ARG types among NHPs and human based on relative abundance.(E) Chao1 index comparison of ARGs in NHPs and human.

FIGURE 5
FIGURE 5 Co-occurrence networks of ARGs and microbiota taxa.(A) Procrustes analysis results showing correlations between ARGs.(B) Procrustes analysis results demonstrating correlations between ARG composition and microbial taxonomy.(C) Network analysis illustrating co-occurrence patterns among ARG subtypes.Node size indicates connection strength (degree) with strong (Spearman's correlation coefficient ρ > 0.6) and significant (p value <0.01) correlations.(D) Network analysis illustrating co-occurrence patterns between ARG subtypes and microbial taxa.Node size indicates connection strength (degree) with strong (Spearman's correlation coefficient ρ > 0.6) and significant (p value <0.01) correlations.