Graphical Abstract

Introduction

Ticks are obligate hematophagous arthropods that predominantly depend on human and animal host blood for development. Their importance as vectors of human and animal pathogens is highlighted by the substantial resources channeled towards controlling the several tick-transmitted pathogens. The Gulf Coast tick (Amblyomma maculatum) is one of the important human-biting ticks in the United States, along with Amblyomma americanum, Ixodes scapularis, and Dermacentor variabilis. Amblyomma maculatum (Am. maculatum) is the most geographically widespread tick species in the southeastern United States (Sonenshine, 2018). It is the competent vector of Rickettsia parkeri (R. parkeri), a spotted fever group rickettsial (SFGR) pathogen that causes mild infections in the United States (Paddock and Goddard, 2015).

Recent advances in the control of ticks and tick-transmitted pathogens have led to the discovery of microbial communities that co-exist alongside pathogens (Narasimhan and Fikrig, 2015; Bonnet et al., 2017; Wei et al., 2021; Wu-Chuang et al., 2021; Narasimhan et al., 2021). The functional role of the tick microbiome has only been partly clarified by tick research. Overall, the primary focus of tick microbiome research has been on the microbial profiling of field-collected ticks to understand the basic bacterial profiles/communities present. The tick microbiome is a moving target and investigating its dynamics within the tick vector is a continuous effort (Narasimhan and Fikrig, 2015; Bonnet et al., 2017; Greay et al., 2018; Narasimhan et al., 2021). Ticks specialize in an exclusive diet of vertebrate blood and have evolved intimate interactions with beneficial symbionts that provide essential B vitamins and cofactors deficient or not readily accessed in the blood diet (Adams, 1963; Gottlieb et al., 2015; Hunter et al., 2015; Smith et al., 2015; Bonnet et al., 2017; Duron et al., 2017; Guizzo et al., 2017; Duron et al., 2018; Olivieri et al., 2019; Ben-Yosef et al., 2020). Correlative evidence also suggests that the tick reproductive fitness may benefit from essential nutrients provided by tick endosymbionts, though this is yet to be experimentally proven (Zhong et al., 2007). Additionally, it has been discovered that specific microorganisms within the tick and other blood-feeding arthropods perform immunomodulatory roles. These microorganisms prime the immune system of the host to continually produce low quantities of antimicrobial peptides, and either inhibit or promote the multiplication of other microbes (de la Fuente et al., 2017; Budachetri et al., 2018; Fogaça et al., 2021; Boulanger and Wikel, 2021; Bonnet and Pollet, 2021).

Several studies have reported on the structure and microbial composition of field-collected and lab-maintained tick populations (Budachetri et al., 2014; Budachetri et al., 2017; Karim et al., 2017; Ross et al., 2018; Varela-Stokes et al., 2018; Portillo et al., 2019; Thapa et al., 2019a; Zhang et al., 2019; Adegoke et al., 2020; Gil et al., 2020; Guizzo et al., 2020; Kumar et al., 2022; Maldonado-Ruiz et al., 2021; Rojas-Jaimes et al., 2021). These studies have revealed a variety of factors responsible for the differences among microbial communities between different tick samples, thus contributing extensively to our understanding of how microbial communities residing within the tick contribute to several aspects of tick biology (Ponnusamy et al., 2014; Budachetri and Karim, 2015; Trout Fryxell and DeBruyn, 2016; Budachetri et al., 2018; Thapa et al., 2019b; Adegoke et al., 2020; Brinkerhoff et al., 2020; Kueneman et al., 2021). However, the structure of Am. maculatum microbial communities and how these communities are shaped by blood meal and pathogen interactions have been the focus of two studies (Budachetri et al., 2014; Budachetri et al., 2018). Previous microbiome studies of Am. maculatum identified unique differences in microbial composition and abundance of specific bacterial genera in field-collected and questing tick tissues (Budachetri et al., 2014). The ability of a tick to serve as a competent vector for any pathogen goes beyond the establishment of the pathogen within the tick tissue. It also extends to the facilitation of pathogen transmission at the tick–host interface, which is made possible in part by salivary proteins. This process is described as saliva-assisted transmission (Nuttall and Labuda, 2008; Šimo et al., 2017; Nuttall, 2019; Karim et al., 2021). Whether the tick salivary microbiome interferes or facilitates pathogen colonization and transmission remains an unanswered question. While a study by Budachetri et al. (2014) reported several Rickettsia species in the saliva of Am. maculatum, only R. parkeri was detected in the salivary glands of Am. maculatum. Limited information exists about how microbial community members interact across tick developmental stages and tissues and how these tissues shape or modulate microbial community assemblages during tick physiological changes, such as those following blood-feeding or pathogen interactions. An elegant study by Lejal et al. (2021) identified several positive and negative correlations between members of the field-collected, questing Ixodes ricinus microbiome and how some of these interactions are modulated by tick-borne pathogens. It has been established that each tick tissue (organ) undergoes unique responses to specific physiological changes, but it is not known whether these responses also induce changes in the microbial communities residing within these tissues (Franta et al., 2010; Schwarz et al., 2014; Crispell et al., 2016; Kumar et al., 2016; Budachetri et al., 2017; Starck et al., 2018; Tirloni et al., 2020; Kurokawa et al., 2020). In the current study, a high-throughput sequencing approach was utilized to examine microbiome changes throughout the ontogeny of uninfected and pathogen-infected Am. maculatum ticks. Characterization of unique microbe-microbe interactions (measure of the degree of relationships between genera) at the different tick developmental stages and tissues revealed potential synergistic and antagonistic interactions among the bacterial communities.

Materials and methodsEthics statement

All animal experiments were conducted in strict accordance with the recommendations of the Guide for the Care and Use of Laboratory Animals of the National Institutes of Health, USA. The protocols for the blood-feeding of immature and mature development stages of Am. maculatum were approved by the Institutional Animal Care and Use Committee of the University of Southern Mississippi (protocols #15101501.2 and 15011402).

Tick maintenance

Rickettsia parkeri-infected and uninfected A. maculatum ticks were maintained under laboratory conditions as previously described (Budachetri et al., 2014; Budachetri et al., 2017; Budachetri et al., 2018). The colonies were established from field-collected adult ticks collected in 2011 from the Mississippi Sandhill Crane, National Wildlife Refuge (https://www.fws.gov/refuge/mississippi_sandhill_crane/) using the drag-cloth method. To generate replete adult females, unfed male and female Am. maculatum ticks were infested on sheep and allowed to feed to repletion and drop off. Fully engorged female ticks were then allowed to oviposit in individual snap-cap vials with perforated lids protected by breathable mesh cloths to allow for oxygen exchange while preventing their escape. Unfed larvae that emerged from the eggs were infested on golden Syrian hamsters until fully engorged, which were then allowed to molt into nymphs and then reinfested on golden Syrian hamsters. The fully fed nymphs were subsequently allowed to molt into adult male and female ticks. The presence or absence of R. parkeri at each developmental stage was confirmed following previously established nested PCR, qPCR and immunofluorescence. Clutches of eggs from ovipositioning females were tested for R. parkeri infection and infected females were separated from uninfected female (Budachetri et al., 2014; Budachetri et al., 2018; Guizzo et al., 2022).

Tick life stages and DNA extraction

DNA was isolated from freshly laid eggs as well as unfed and partially fed ticks at various developmental stages: nymphs, male and female adults. To remove surface contaminants and bacterial carryover, all samples from developmental stages were surface-sterilized using 2% sodium hypochlorite, followed by two 5-minute washes in 70% ethanol and two 5-minute washes in sterile deionized water (Kumar et al., 2022). Tissues were dissected from unfed and partially blood-fed adult females and directly stored in tissue lysis ATL buffer (contains sodium dodecyl sulphate as active component) before extracting DNA. Female ticks were used for tissue collection because of their ability to effectively take a blood meal and increase their size throughout the host attachment period compared to males that take minimal blood meal and shows no physical changes. Tick dissection was performed in a biosafety cabinet while changing blades (Personna American Safety Razor, Mexico) between dissections. Saliva was collected on ice by endogenous induction of salivation following injection of pilocarpine, as previously described (Ribeiro et al., 2004; Valenzuela et al., 2000). DNA was extracted using the DNeasy Blood and Tissue Kit following the manufacturer’s protocol (QIAGEN, Germantown, MD, USA) in a biosafety cabinet to eliminate potential carryover of bacteria between samples. Control samples included a blank control and a mock bacterial sample with a known bacterial population (ZymoBIOMICS™ Microbial Community DNA Standard, Zymo Research, Irvine, CA, USA) as part of the DNA-extraction protocol.

Library preparation for Illumina 16S sequencing

Library preparation and indexing were carried out according to Illumina 16S metagenomic sequencing. The V3–V4 region of the 16S rRNA gene was amplified using DNA samples from different life stages and tissues. A two-step PCR amplification was carried out, which included an initial V3–V4 amplification using primers with specific Illumina adapter overhang nucleotide sequences using the following thermocycling conditions: 95°C for 3 min; followed by 25 cycles of 95°C for 30 s, 55°C for 30 s, and 72°C for 30 s; and a final extension step of 72°C for 5 min. The PCR amplicon was purified using AmPure XP beads (Agencourt Bioscience Corporation, Beverly, MA, USA) following the manufacturer’s protocol, and the purified product was eluted in 50 µl of TE buffer (Qiagen cat. no. 19086, Qiagen, Germantown, MD, USA) and analyzed on a 2% gel to confirm a single DNA band with an estimated amplicon size of 550 bp. The second step involved attaching unique indexes, in the form of nucleotide sequences, to the purified V3–V4 PCR product. Commercially available dual indexes (i5 and i7) from the Illumina Nextera index kit V2 (Illumina, San Diego, CA) were used for the second (index) PCR. Briefly, a reaction including each of the forward and reverse Nextera index primers, high fidelity 1X KAPA HiFi HotStart ReadyMix (Kapa Biosystems, Wilmington, MA, USA), and 5 μl of each purified V3–V4 amplicon was set up using the following thermocycling conditions: 95°C for 3 min; followed by 8 cycles of 95°C for 30 s, 55°C for 30 s, and 72°C for 30 s; and a final extension step of 72°C for 5 min. The PCR product was briefly purified using AmPure XP beads (Agencourt Bioscience Corporation, Beverly, MA, USA), following the manufacturers protocol. The purified product was eluted in 25 µl of TE buffer and analyzed on a 2% gel to confirm a single DNA band with an estimated amplicon size of 630bp. The purified indexed PCR library products were quantified using the KAPA Library Quantification Kit from Roche (cat. no. 07960255001, kit code KK4844, Hoffmann-La Roche, Basel, Switzerland), normalized to a concentration of 7 nM, and pooled with 5 µl of each normalized library. Biological replicates from DNA extraction controls, TE buffer, and mock bacterial communities (ZymoBIOMICS™ Microbial Community DNA Standard, Zymo Research, Irvine, CA, USA) were simultaneously processed alongside the tick samples. The purified library was sequenced in a single run of an Illumina MiSeq sequencing instrument using reagent kit v2 (500 cycles) with 2 X 250 bp output at the University of Mississippi Medical Centre (UMMC) Genomics Core Facility.

Data processing and analysis

Unless otherwise stated, all preprocessing was performed following the video tutorial of the Quantitative Insights into Microbial Ecology (QIIME2) pipeline (Bolyen et al., 2019). Briefly, demultiplexed fastq files were unzipped using the “unzip -q casava-18paired-end-demultiplexed.zip” command, followed by merging individual forward and reverse fastq files into a single fastq file. The Atacama soil microbiome pipeline was incorporated to control demultiplexed paired-end reads using the DADA2 plugin as previously described (Callahan et al., 2016). Low-quality sequences were trimmed and filtered out, and subsequent merging of paired-end reads ensured 20 nucleotide overhangs between forward and reverse reads. Chimeric sequences were removed from the sequence table. The qiime2 feature-classifier plugin v. 7.0 was used for taxonomic assignment against the pre-trained Silva_132_99% classifier (Quast et al., 2013). The “qiime diversity core-metrics-phylogenetic” command was used to compute diversity metrics, after which the “QIIME diversity alpha-group-significance” command was used to explore microbial diversity using Faith_pd and evenness metrics as measures of community richness and evenness, respectively. To determine the taxa that significantly differ between groups, we used differential abundance analysis (DAA) with the OTU table (at genera-level, taxon = Genus) as input data. Post-hoc test (the false discovery rate, FDR) was used to determine significantly different (FDR p-value < 0.05) taxa between groups. Raw sequences were submitted to DRYAD repository (https://doi.org/10.5061/dryad.r2280gbgg).

Functional characterization of microbial communities

To determine the functionality of the microbial genome, we predicted the functional pathways using the Phylogenetic Investigation of Communities by Reconstruction of Unobserved States 2 (PICRUSt2) pipeline (Ye and Doak, 2009; Mirarab et al., 2012; Louca and Doebeli, 2018; Czech et al., 2020; Douglas et al., 2020). We assigned functional predictions by normalizing the bacterial sequences and OTUs to the KEGG Ortholog (KO) databases and Clusters of Orthologous Groups (COGs) to predict metabolic pathways.

Visualization

Microbiome Analyst, a web-based interface, was used for data visualization using taxonomy and metadata tables generated from data processing as input files (Dhariwal et al., 2017; Chong et al., 2020). Low count and prevalence data were filtered from the OTU table by setting values of 10 and 20, respectively. A filtered abundance table was exported and used in generating histograms of bacterial abundance in Microsoft Excel 2016 (Microsoft, 2018). Boxplots of alpha diversity were generated with GraphPad Prism version 8 for Windows (GraphPad Software, La Jolla CA, USA, www.graphpad.com) using the pairwise table generated from the “qiime diversity alpha-group-significance” command. Network correlation maps were constructed based on the sparse correlations for compositional data (SparCC) approach (Friedman and Alm, 2012). This approach uses the log-transformed values to carry out multiple iterations, which subsequently identifies taxa outliers to the correlation parameters (Chong et al., 2020). To compare the differences in the microbiome between tick groups based on measures of distance or dissimilarity, a matrix was generated from log-transformed sequence data, and ordination of the plots was visualized using both Principal Coordinates Analysis (PCoA) and non-metric multidimensional scaling (NMDS). The Bray–Curtis distance matrix was used to visualize compositional differences in the microbiome across all groups.

Statistical analysis

The Wilcoxon rank-sum test estimated significant correlations in diversity indexes across developmental stages and tissue data sets. A Kruskal–Wallis test followed by Dunn’s multiple comparison test was used in comparing the effect of blood meal and R. parkeri on microbial richness between two different groups (unfed vs. fed or uninfected vs. infected). Statistically significant data were represented as P <0.05.

ResultsSequencing results

A total of 91 samples, representing at least three biological replicates across developmental stage (egg, larvae, nymph, and adult), tissue (salivary gland, midgut, and ovary), and body fluids (saliva) from uninfected and R. parkeri-infected ticks, generated a total of 4,090,811 forward and reverse reads, with a mean of 45,453.5 and a maximum of 112,242 reads per sample (Additional File, Table S1). Rarefaction analysis of the individual samples from a sequencing depth of 500 to 4000 confirmed that there was adequate sequence coverage relative to the number of operational taxonomic units (OTUs) and a plateau for individual tick samples (Additional File, Figures S1A, B). Overall, 711 unique OTUs were identified, with 254 identified from each developmental stage and tissue.

Microbiome dynamics across developmental stages

Proteobacteria represented the most abundant phylum throughout ontogeny, with a relative abundance of 84%, followed by Actinobacteria, Firmicutes, and an uncharacterized phylum (Additional File, Figure S2A). An extended breakdown of each identified phylum and the number of member bacteria at the genus level is represented using pie charts for each phylum, as shown in Additional File, Figures S2B–D. The breakdown of the top 10 most frequent bacterial genera shows that Francisella and Rickettsia were the predominant members of the Am. maculatum microbiome across all developmental stages (Additional File, Figure S3A). We observed a progressive decrease in the relative abundance of the genus Cutibacterium throughout each developmental stage, from egg to adult male and female ticks (Additional File, Figure S3A). Other bacterial genera, such as Rhodococcus, Pseudomonas, and Candidatus_Midichloria, were also detected at relatively lower abundance. The abundance of representative bacterial genera was also compared between blood-fed and unfed life stages. Uptake of a blood meal was observed to impact the relative abundance of Francisella and Cutibacterium (Additional File, Figure S3B). Fed larvae showed a decrease in Francisella abundance and an increase in the abundance of Cutibacterium compared to unfed larvae.

Interestingly, opposing trends were observed when blood-fed nymphs were compared with unfed nymphal ticks. The relative abundance of Cutibacterium was reduced, whereas Francisella abundance was observed to increase (Additional File, Figure S3B). To examine the modulation of microbial assemblages in the presence of pathogenic microbes, we compared the relative abundances of bacteria between uninfected and R. parkeri-infected samples from each developmental stage. There is a clear trend of decreasing Francisella abundance, and the corresponding increase in Rickettsia and Candidatus_Midichloria across all R. parkeri-infected stages compared with uninfected samples (Additional File, Figure S3C). Notably, the Cutibacterium load was negatively impacted by the presence of R. parkeri when infected and uninfected nymphs were compared (Additional File, Figure S3C).

Eggs exhibited higher bacterial richness and evenness, regardless of blood meal or infection (Wilcoxon rank-sum test p < 0.05). At the same time, nymphs showed the least microbial richness and evenness (Wilcoxon rank-sum test p < 0.001) (Figures 1A, B). No significant microbial richness in the blood-fed samples was observed; however, blood meal was seen to exert a significant impact on species proportions in unfed nymphs when compared with unfed larvae (Kruskal–Wallis test followed by Dunn’s multiple comparison test, p < 0.05) and in the fed nymphs compared with unfed larvae (Kruskal–Wallis test followed by Dunn’s multiple comparison test, p < 0.001) when the Pielou_e evenness index was estimated (Figures 1C, D) Overall, blood-feeding had no significant effect on microbial richness compared with unfed ticks across all developmental stages.

Figure 1 Alpha diversity of microbial richness across all developmental stages and isolated tissues. Diversity analysis using Faith phylogenetic distance (species richness) and Pielou_e metrics (evenness) in (A, B) all developmental stages, (C, D) unfed and fully fed developmental stages, (E, F) R. parkeri-infected and uninfected developmental stages, (G, H) all tissues, (I, J) unfed and fully fed tissues, and (K, L) R. parkeri-infected and uninfected tissues. Each bar on the abundance plot represents average data from 3–5 individual replicates. Each point of the box plots represents average data from 3–5 individual replicates. Unfed (UF), fed (F), larvae (LV), nymph (N), male (M), female (F), clean (CL), infected (IN). Boxplot height indicates the upper and lower quartiles. Whiskers indicate the maximum and minimum values. The black horizontal line indicates the median. Asterisks denote significance (*P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001).

The presence of R. parkeri infection was not found to affect microbial richness and evenness between the developmental stages. However, a reduction in microbial richness in R. parkeri-infected nymphs compared with uninfected eggs and infected eggs was noted (Figure 1E). Significant reduction in specie evenness was observed when uninfected and infected eggs was compared to uninfected and infected nymph respectively (Figure 1F). NMDS analysis of community structure using the Bray–Curtis (BC) distance matrix (PERMANOVA, F-value = 22.818; R2 = 0.75015; p-value < 0.001 [NMDS] Stress = 0.12544) revealed two distinct clustering patterns in which the earlier developmental stage (egg and larvae), were uniquely clustered from the adult tick (male, female) (Additional File, Figure S4). The BC distance matrix measured dissimilarity in microbial communities based on compositional differences among our samples. The observed clustering pattern further supports the differences observed in the presence of Cutibacterium in the eggs, larvae, and nymphs but not in adult male and female ticks. Overall, both Francisella and Rickettsia predominated in the microbiome across all the developmental stages. Both blood meal and R. parkeri infection negatively affected the abundance of Francisella, while blood meal reduced microbial richness in earlier developmental stages.

Microbiome composition and diversity across tissue niches

As with developmental stages, the phylum Proteobacteria had the greatest abundance of bacteria, with a relative abundance of 86%, followed by Firmicutes, Actinobacteria, and an unknown phylum, as illustrated in the Additional File, Figures S5A–D. Representative sequences generated from salivary gland, midgut and ovarian tissues showed a similar trend to what was found for developmental stages. Francisella and Rickettsia predominated in the microbiomes of isolated tissues. In addition, the genus Candidatus_Midichloria was detected at relatively high abundance across the isolated tissues, except for the saliva microbiome, which contained the genus Bacillus. By contrast, the genus Cutibacterium was present in all tissues, except for ovarian samples (Additional File, Figure S6A). Following blood uptake, the relative abundances of certain bacteria identified in the partially fed salivary gland was higher than in the unfed salivary gland. Interestingly, blood meal reduced the abundance of Francisella, whereas the prevalence of Rickettsia, Cutibacterium, and Candidatus_Midichloria increased (Additional File, Figure S6B). The effect of blood meal on the microbial abundance in the midgut tissues was different from salivary gland; in the reverse of what was observed in the salivary gland, blood meal increased the relative abundance of Francisella and Candidatus_Midichloria, while reducing that of Rickettsia. Due to the unique developmental trajectory of the female tick, in which ovarian development did not commence until 4–5 days post-blood-feeding, ovarian tissues from unfed female ticks were not included in this analysis (Additional File, Figure S6B). Bacterial communities were significantly different between unfed and partially blood-fed salivary gland (PERMANOVA F-value = 22.766; R2 = 0.7779; p-value < 0.001) and midgut tissues (PERMANOVA F-value = 16.167; R2 = 0.69785; p-value < 0.001), according to PCoA analysis of the Bray–Curtis’s distance matrix (Additional File, Figures S7A, B).

The presence of R. parkeri was associated with the corresponding increase in the relative abundance of Candidatus_Midichloria in all infected compared with uninfected tissues, with the corresponding decrease in the prevalence of Francisella (Additional File, Figure S6, Figure 2C). There was an increase in the relative abundance of Rickettsia detected in R. parkeri-infected saliva compared with saliva from uninfected ticks. Bacterial communities were significantly different between uninfected and infected salivary gland (PERMANOVA F-value = 20.316; R2 = 0.75761; p-value < 0.001), midgut (PERMANOVA) F-value = 36.225; R2 = 0.83806; p-value < 0.001) and ovarian tissues (PERMANOVA F-value = 29.509; R2 = 0.89397; p-value < 0.008), according to PCoA analysis of the Bray–Curtis distance matrix (Additional File, Figures S8A–C).

Figure 2 DeSeq2 differential abundance analysis of OTUs (genus-level) between (A) all life stages (B) R. parkeri infected and uninfected, and (C) unfed and partially fed Am. maculatum. For multiple comparisons an estimate of the FDR was calculated to determine differentially abundant families with a significance threshold of p-value < 0.05. IN (R. parkeri infected), UN (uninfected), F (blood fed), UF (unfed).

Overall, microbial richness was highest in midgut tissues (Wilcoxon rank-sum test p < 0.001) and lowest in the ovary (Wilcoxon rank-sum test p < 0.05), based on the Faith_phylogenetic distance index. At the same time, no significant differences were found in microbial proportions across different tissues, including saliva (Figures 1G, H). The impact of blood meal on alpha diversity between unfed and partially fed tissues was assessed. Unfed midgut tissues had significantly higher specie richness compared to saliva and ovarian tissues, while partially fed midgut tissues displayed a higher species richness compared to the saliva (Figure 1I). We observed no significant effect of blood meal on species proportions (Figure 1J). Similarly, there was no significant difference in either species richness or evenness between R. parkeri-infected and uninfected tissues (Figures 1K, L). Here we found that blood meal and R. parkeri infection shaped the relative abundance of Am. maculatum microbiome in a tissue-dependent manner.

Distinct microbial genera are affected with blood meal or R. parkeri infection in developmental stages and tissues

Differential abundance analysis (DAA) was performed on OTU at the genera level to determine which members of the bacterial community were associated with the significant changes blood meal and R. parkeri infection. To get the useful information, we selected bacterial genera with overall abundance more than 0.1% to avoid sequencing errors, and any genera for which taxonomy was not available.

Differential abundance analysis (DAA) analysis on the OTU table at genus level indicated higher relative abundance of Francisella in female and nymphs, while Stenotrophomonas was highest in female and male ticks (Figure 2A). DAA of ticks following R. parkeri infection indicated Francisella was uniquely abundant in uninfected ticks while Candidatus_Midichloria was uniquely abundant in R. parkeri infected ticks (Figure 2B). In addition, Francisella was least abundant in eggs while Escherichia, Stenotrophomonas and Brachybacterium were all relatively abundant in unfed female, male and nymph (Figure 2C). DAA of tissue-derived OTUs indicated a higher relative abundance of Francisella in ovary and salivary gland while Rickettsia was higher in midgut and saliva (Figure 3A). DAA indicated higher relative abundance of Rickettsia in unfed and partially fed midgut (Figure 3B). The presence of R. parkeri was associated with a higher relative abundance of Candidatus_Midichloria and Rickettsia in infected tissues and a higher relative abundance of Francisella in uninfected tissues following DAA analysis of tissue-derived OTUs (Figure 3C). On the other hand, results of DAA on OTU table at the genus level were in accordance with the observations made about bacterial community structure at the tissue level and indicated the synergism between R. parkeri and Candidatus_Midichloria and antagonism between R. parkeri and Francisella (Additional File, Figures S6A-C).

Figure 3 DeSeq2 differential abundance analysis of OTUs (genus-level) between (A) whole tissues (B) R. parkeri infected and uninfected, and (C) unfed and partially fed Am. maculatum. For multiple comparisons an estimate of the FDR was calculated to determine differentially abundant families with a significance threshold of p-value < 0.05. PF (partially blood fed), UF (unfed), IN (R. parkeri infected), UN (uninfected).

Microbiota interactions across developmental stages

We performed network analysis to visualize significant positive and negative correlations between the bacterial genera identified. We observed a total of 286 significant partial interactions (measure of the degree of relationships between genera) between 43 bacterial genera, 52% (n=150/286) of which showed positive correlations, while the remaining 48% (n=136/286) exhibited negative correlations (Additional File, Table S2). All identified bacteria were from the phyla Actinobacteria (10 genera), Firmicutes (12 genera), Spirochaetota (1 genus), Proteobacteria (12 genera), Cyanobacteria (1 genus), Bacteroidota (4 genera), and Fusobacteria (1 genus) as well as an unknown phylum (1 genus) (Figure 4). The genera Francisella, Rickettsia, and Candidatus_Midichloria were found within the same network cluster where they share significant positive correlations. Individually, they both interacted with other genera as is the case with Francisella and Rickettsia OTUs sharing negative correlations with Listeria, Coxiella, Bacillus, and Enterococcus, and positive correlations with Cutibacterium. The genus Candidatus_Midichloria on the other hand was positively correlated with OTUs belonging to the Mitochondria genus. It is worth noting that Candidatus_Midichloria had fewer taxa with which it interacted than Francisella or Rickettsia (Additional File, Table S2). Another interesting observation made from the network interactions was the significant negative correlations in the interactions between Coxiella and Francisella and between Coxiella and Rickettsia. A more surprising correlation was identified between Bacillus, Listeria, Enterococcus, and Staphylococcus, all of which are members of the phylum Firmicutes. They were all negatively correlated with Rickettsia, Francisella, and Candidatus_Midichloria from Proteobacteria (Figure 4).

Figure 4 Correlation network analysis across the life stages of Am. maculatum ticks. Correlation network generated using the SparCC algorithm. Correlation network with nodes representing taxa at the family level and edges representing correlations between taxa pairs. Node size correlates with the number of interactions in which a taxon is involved, and the line color represents either positive (red) or negative (blue) correlations. The color-coded legend shows the bacterial family to which each taxon belongs.

To identify the significantly represented microbial genera between different life stages, we used linear discriminant analysis effect size (LEfSe) to determine the respective LDA scores for the identified bacteria observed in the samples. LDA analysis showed that, out of 46 bacteria taxa, 13 genera were significantly correlated at different life stages. Several of the bacteria either increased or decreased in more than one life stage. By contrast, the bacterial taxon Enhydrobacter is the only genus represented exclusively at the egg stage (Additional File, Table S3). Likewise, the taxa Enterococcus, Lactobacillus, and Coxiella were not detected at any developmental stages. Another interesting observation was the increased presence of the genera Corynebacterium and Cutibacterium exclusively at early developmental stages (Additional File, Figure S9).

Tissue-driven microbial interactions

Network analysis in tissues showed 220 significant correlations between 46 bacterial genera, each belonging to the phyla Proteobacteria (20 genera), Firmicutes (13), Actinobacteriota (6 genera), Cyanobacteria (2 genera), Bacteroidota (3 genera), and Fusobacteriota (1 genus) as well as an unknown phylum (1 genus). Seventy-two (33%) were negative correlations, while 148 (67%) were positive correlations (Additional File, Table S4). We observed that OTUs belonging to Candidatus_Midichloria, Francisella and Rickettsia were positively correlated and represent a network cluster that shares both positive and negative correlations with OTUs from other bacteria genera. Interestingly, bacteria belonging to the phylum Firmicutes, such as Listeria, Staphylococcus, Bacillus, and Enterococcus were negatively correlated with Francisella, Rickettsia, and Candidatus_Midichloria which are bacteria from the Proteobacteria phylum (Figure 5). LDA analysis revealed that Francisella was significantly enriched in ovaries and salivary glands, whereas Rickettsia was increased significantly in the midgut tissues and saliva (Additional File, Table S5; Figure S10). The saliva microbiome was also the only sample exclusively enriched with Alcaligenes, Moraxella, Pseudomonas, and Spiroplasma. The first three belong to the phylum Proteobacteria and the latter from the phylum Firmicutes.

Figure 5 Correlation network analysis across tissue stages of Am. maculatum ticks. Correlation network generated using the SparCC algorithm. Correlation network with nodes representing taxa at the family level and edges representing correlations between taxa pairs. Node size correlates with the number of interactions in which a taxon is involved, and the line color represents either positive (red) or negative (blue) correlations. The color-coded legend shows the bacterial family to which each taxon belongs.

Microbial interactions in the salivary gland

A total of 608 significant partial interactions between 48 bacteria genera were identified in network analysis of the salivary gland, with 176 (29%) negative correlations and 432 (71%) positive correlations (Additional File, Table S6). Most interactions were seen with Francisella, Rickettsia, and Candidatus_Midichloria (Figure 6). The most surprising correlation was the negative correlation observed between Rickettsia and Candidatus_Midichloria, and the most interesting was that Rickettsia shares a correlation network only with Candidatus_Midichloria (Figure 6). No direct interaction was observed between Francisella and Candidatus_Midichloria or between Francisella and Rickettsia. LDA analysis of significantly enriched bacteria revealed a unique pattern between R. parkeri-infected and uninfected salivary glands. Uninfected salivary glands were exclusively enriched with Francisella, Escherichia_Shigella, Cutibacterium, and Pseudomonas. In contrast, Candidatus_Midichloria, Enterococcus, and Cutibacterium were the only bacteria significantly enriched in R. parkeri-infected salivary glands (Additional File, Figure S11A). Interestingly, Francisella was the only bacterium to be significantly enriched in unfed salivary gland. Francisella was also enriched in the partially fed salivary gland, albeit at a much lower expression compared to the unfed salivary gland. Several bacteria were also significantly enriched in partially fed salivary gland which include Escherichia_Shigella, Bacillus, Staphylococcus, Enterococcus, Cutibacterium, and Pseudomonas (Additional File, Figure S11B).

Figure 6 Correlation network analysis in the salivary glands of Am. maculatum ticks. Correlation network generated using the SparCC algorithm. Correlation network with nodes representing taxa at the family level and edges representing correlations between taxa pairs. Node size correlates with the number of interactions in which a taxon is involved, and the line color represents either positive (red) or negative (blue) correlations. The color-coded legend shows the bacterial family to which each taxon belongs.

Microbial interaction in the midgut

We observed many fewer (56) significant partial interactions (measure of the degree of relationships between genera) in the midgut tissues between 29 bacteria genera, all belonging to the phyla Proteobacteria (15 genera), Firmicutes (3 genera), Actinobacteriota (4 genera), Cyanobacteria (1 genus), and Bacteroidota (2 genera) as well as an unknown phylum (1 genus), 16 (29%) of which were negative correlations and the remaining 40 (71%) positive correlations (Additional File, Table S7). The network map of bacterial interactions identified five distinct networks, with Rickettsia and Candidatus_Midichloria observed to be part of the same network with positive correlations. Francisella was found to be in a different network (Figure 7). A striking observation was the ubiquitous presence of bacteria in the phylum Proteobacteria across all five identified networks (Figure 7). The LDA score from LEfSe identified the enrichment of Rickettsia, Escherichia_Shigella, and Cutibacterium in both uninfected and R. parkeri-infected midgut tissues. Compared with the R. parkeri-infected midgut, only Francisella was significantly enriched in the uninfected midgut, while Bacillus was enriched in the R. parkeri-infected midgut (Additional File, Figure S11A). Assessing the effect of blood meal on bacterial enrichment identified a shared enrichment of Rickettsia, Escherichia_Shigella, and Cutibacterium in both unfed and partially fed midgut. We observed Pseudomonas and Brevibacterium in the unfed midgut, while Francisella and Staphylococcus were enriched in partially blood-fed midgut tissues (Additional File, Figure S11B). Microbiota–microbiota interactions across developmental stages were driven mainly by three of the most abundant bacteria: Francisella, Rickettsia, and Candidatus_Midichloria. Other interactions revealed that certain taxonomic groups were inversely correlated, as seen between the phyla Firmicutes and Proteobacteria. Many specific interactions were observed in the tissues, suggesting tissue-dependent influences on how microbes interact in different tissue niches.

Figure 7 Correlation network analysis in the midgut of Am. maculatum ticks. Correlation network generated using the SparCC algorithm. Correlation network with nodes representing taxa at the family level and edges representing correlations between taxa pairs. Node size correlates with the number of interactions in which a taxon is involved, and the line color represents either positive (red) or negative (blue) correlations. The color-coded legend shows the bacterial family to which each taxon belongs.

Functional characterization of bacterial communities

Functional prediction of microbiota genes compared with COG functional categories identified genes involved in ribosomal structure and biogenesis, nucleotide transport, and metabolism as the most abundant across developmental stages (Additional file, Table S8). Several genes involved in metabolisms, such as lipid transport, coenzyme transport, amino acid transport, and carbohydrate transport, were also identified (Additional files, Figure S12). Prediction of metabolic genes and pathways using the KEGG Orthology (KO) categories assigned 19% of the metabolic functions to genes in carbohydrate metabolism pathways and 18% of the functions were encoded by genes in amino acid metabolism pathways. By comparison, genes encoding cofactors in the vitamin metabolism pathways and energy metabolism pathways represented 15% and 13% of the metabolic pathways, respectively (Additional file, Figure S12).

Differences in functional community profiles within the metabolic pathways were further examined between the different feeding stages and infection status across dissected tissues based on the KEEG metabolism categories. The partially fed salivary gland and midgut had a higher proportion of sequences assigned to energy metabolism and carbohydrate metabolism than the unfed salivary gland and midgut (Figures 8A, B). The presence of R. parkeri led to an increase in the proportion of glycan biosynthesis and metabolism genes, carbohydrate metabolism genes, and genes involved in the metabolism of cofactors and vitamins in the salivary glands. A slight increase was also seen in genes involved in the metabolism and biodegradation of xenobiotic in the R. parkeri-infected salivary gland (Figure 9A). A slight increase in the proportion of energy metabolism genes was associated with R. parkeri in infected midgut compared with the uninfected midgut (Figure 9B). Interestingly, genes involved in amino acid metabolism were seen to be reduced in R. parkeri-infected salivary glands and midgut compared with uninfected tissues (Figures 9A, B). No differences between metabolic predictions were detected between R. parkeri-infected and uninfected ovarian tissues (Additional file, Figure S13).

Figure 8 Histogram indicating the functional differences between the fed and partially fed Am. maculatum microbiota in (A) salivary gland and (B) midgut tissues. KEGG metabolic categories were obtained from 16S rRNA gene sequences using the PICRUSt2 pipeline.

Figure 9 Histogram indicating the functional differences between the uninfected and R. parkeri-infected Am. maculatum microbiota in (A) salivary gland and (B) midgut tissues. KEGG metabolic categories were obtained from 16S rRNA gene sequences using the PICRUSt2 pipeline.

Discussion

This study is the first to characterize microbiome changes throughout the ontogeny of the Gulf Coast tick, a competent vector of the spotted fever group rickettsiae, R. parkeri. Using a high-throughput metagenomic approach that allowed us to barcode DNA extracted from individual ticks and process them in a single run, the present study set out with the aim of assessing the core microbial assemblages present in Am. maculatum ticks and determine the changes under certain physiological conditions vis-à-vis development, hematophagy, and pathogen infection. The results show that changes to the core microbiome occur in response to a blood meal and pathogen association but not as a function of the developmental stage. Similarly, as expected, microbe-microbe interactions depended on the specific tick tissues in which those interactions took place.

The Amblyomma maculatum microbiome remains stable across ontogeny

The typical life cycle of Am. maculatum ticks encompass critical physiological and environmental events conferred by intrinsic and extrinsic factors, such that the persistence of the microbiome and every similar micro-community that utilizes the tick as a microenvironment across these dynamic events is continuously tested. Two very well-studied factors that exert significant pressure on the tick microbiome are blood meal and pathogen replication (Abraham et al., 2017; Swei and Kwan, 2017; Adegoke et al., 2020; Wei et al., 2021), albeit with varying outcomes reported. In the current study, a shift favoring the establishment of a single bacterial genus along developmental lines (eggs to adult male and female) is exhibited (Additional File, Figures S3A–C). Even so, the uptake and breakdown of nutritious blood meals and the intricacies of tick–pathogen interactions involve several detailed mechanistic processes that have been shown to significantly impact tick physiology (de la Fuente et al., 2016; Hoxmeier et al., 2017; Budachetri et al., 2018), the current study shows that members of the Am. maculatum microbial communities exhibit reduced plasticity thus maintaining a stable microbial core. This observation raises the question of whether the tick microbiome is self-regulated in ticks and tick tissues.

Notably, while the core microbiome persists across developmental stages, certain bacterial genera were found to be associated with both blood meal and pathogen infection, as seen with Cutibacterium and Candidatus_Midichloria, and whether the introduction of these bacteria in fewer numbers significantly impacts overall tick biology is a question beyond the scope of this study (Abraham et al., 2017; Budachetri et al., 2018; Adegoke et al., 2020). To further understand how changes across tick developmental stages shape microbial assemblages, it is imperative to assess microbial assemblages and community measures of diversity, such as richness and evenness (Brown et al., 2020). Community diversity is important because the smallest differences at the taxon level induced by blood meal or pathogen interaction could interfere with community structure, which would otherwise not be reflected in the abundance profile. We observed a continuous decrease in community diversity towards later developmental stages (Figures 2A-F). For both diversity metrics employed, a reduction in community richness likely indicated the corresponding decline in all or certain specific bacterial taxa. At the same time, a decrease in evenness would imply that certain bacterial taxa are allowed to increase in abundance at the expense of other bacteria. This finding is consistent with observations made in the kissing bugs, where there is a decrease in diversity associated with ontogenic development (Mann et al., 2020), which is in disagreement with Oliveira et al. (2018), who reported an increase in diversity with ontogenic development. A possible explanation for this discrepancy may be attributed to internal and external changes that occur during metamorphosis. Shedding of the midgut lining and the external cuticles are two of the most important phenomena during development. These may come with significant loss of resident and non-obligate bacteria whose absence poses no detrimental effect on tick biology (Engel and Moran, 2013). There are, however, other possible explanations for the reduction in diversity seen in this study. One such explanation is that bacteria that are not vertically maintained from the female to her eggs are replaced by vertically maintained bacteria. The preservation of vertically maintained species would explain the loss of Cutibacterium in adult male and female ticks. By extension, there would be a reduction in overall diversity. Previous findings have highlighted the importance of initial exposure to bacteria in eliciting a robust immune response upon subsequent exposure to the same bacteria or bacteria from a similar group, a phenomenon referred to as immune priming (Trappeniers et al., 2019; Powers et al., 2020). This could explain the fact that immature developmental stages of ticks are easily infected with bacteria from the environment because of their naïve immune status. These bacteria will not be carried over to the adult stage due to a more developed immune status. The persistence of Francisella, Rickettsia, and Candidatus_Midichloria to the adult male and female further emphasizes their position as candidates for vertical maintenance, as previously reported for bacterial species from the same genus (Baldridge et al., 2009; Wright et al., 2015; Azagi et al., 2017; Mukhacheva and Kovalev, 2017; Budachetri et al., 2018).

Microbial assemblages in tissues are driven by pathogen and host blood

To gain in-depth insight into how different arthropod tissues contribute to tick biology and, most importantly, their role in vector competence, it is imperative to know how resident microbial communities differ from one tissue to the other. Several factors have been identified that significantly affect how the microbiome is shaped in different arthropod tissues, albeit with varying outcomes and findings. In this study, we asked whether differences in microbial composition occur within tissues and whether changes influence these differences in the state of blood-feeding and infection status. Our results clearly show a difference in the abundance and composition of the bacteria identified across tick tissues. A conclusion from our findings is the higher bacterial proportion in saliva, which does not correspond to the bacterial balance in the salivary glands. Tick saliva simply contains a snapshot of the microbes secreted by the tick and is possibly not a good representation of all microbes delivered during prolonged feeding. These differences can be explained in part by the fact that some of the microbial populations identified in the saliva, such as Bacillus and Staphylococcus, are known environmental bacteria, have also been reported to be found on the skin and fur of the animal hosts (Dastgheyb and Otto, 2015; Onyango and Alreshidi, 2018; Ravine, 2019), and may have been acquired during the feeding process. Some of these bacteria could also have been mechanically transferred during saliva collection from the ticks. The presence of the microbiome within tick saliva, while preliminary, could implicate the tick microbiome in modulating the host immune response to tick bite at the host–tick interface. The observation of a higher diversity seen in midgut tissues compared with salivary gland, ovary, and saliva contrasted with the findings of Zolnik et al. (2016), who reported a higher diversity in salivary gland compared with midgut tissue. The tick midgut and associated peritrophic membrane play an important role in blood breakdown and its components by several catalytic proteins. The production of toxic reactive oxygen species (ROS) and reactive nitrogen species (RNS), which follow blood digestion, is also dampened by these midgut proteins, thus enabling midgut-resident microbial replication. The heme component of the blood meal derived by catalytic proteins also counters the potential tissue damage associated with a heme-induced oxidative burst. Another plausible explanation for the higher microbial diversity observed in the midgut reported from this study is that the midgut is the only tissue that directly encounters host blood, and it provides the first line of contact with any bacterium coming in with the blood meal (Houk, 1977; Tarnowski and Coons, 1989; Gaio et al., 2011; Lejal et al., 2019). However, prolonged blood feeding and, by extension, the initiation of blood meal digestion might provide an enabling environment for bacterial replication and growth, hence impacting the dynamic and increasing overall diversity of the midgut microbiome as seen from this study. Due to the importance of the ovarian tissues in developmental success, it is expected that several physiological barriers will be put in place to confer protection and reduce the tick’s interaction with circulating microbes to the barest minimum. The reduced microbial diversity observed in the ovary is consistent with this, as seen in Figure 1G. Several studies have outlined the contribution of invertebrate hemocytes to immune functions while also identifying distinct populations of hemocytes as either circulating or sessile (Kadota et al., 2003; King and Hillyer, 2013; Sigle and Hillyer, 2018). The sessile population of hemocytes is usually found attached to cavity walls, such as the abdomen, thorax, and internal organs. While extensive studies have described how arthropod hemocytes play a significant role in immunomodulation (Pereira et al., 2001; Hillyer, 2009; Liu et al., 2011; Hillyer and Strand, 2014; Fiorotti et al., 2019; Kwon and Smith, 2019), few studies have defined the functional role of sessile populations in the biology of the organism (Sigle and Hillyer, 2016). Hence, it could be hypothesized that circulating or sessile hemocytes would form immunological barriers as part of an extensive physiological barrier to protect a vital organ, such as the ovary, as seen with several human and mammalian blood–tissue barriers (Fröhlich, 2002)

It has already been established that blood uptake, which usually lasts several days in ticks, can introduce a plethora of pathogenic and non-pathogenic microbes into ticks (Bhatia et al., 2018; Zolnik et al., 2018; Landesman et al., 2019). While ticks acquire host blood for nutritional requirements, several host factors present in the blood meal could also potentially contribute to compositional changes that are likely to occur to the microbiome throughout the extended feeding process. For instance, the tick midgut responds to the presence of the blood meal via intricate protein expression patterns (Rachinsky et al., 2008; Schwarz et al., 2014; Perner et al., 2016; Oleaga et al., 2017; Liu et al., 2018) while also undergoing certain morphological and structural adjustments (Tarnowski and Coons, 1989; Franta et al., 2010), all of which could have an impact on the microbial environment in the midgut, as previously reported (Zhang et al., 2014; Zolnik et al., 2018; Landesman et al., 2019). Results from our study show that microbial diversity is impacted following uptake of the blood meal in both salivary gland and midgut tissues (Figure 1K). We also observed a combinatorial effect of blood meal and R. parkeri infection in reducing the abundance of Francisella while increasing Rickettsia abundance in the partially fed and infected salivary gland and, vice-versa, in the partially fed and infected midgut. This suggests that there is the trafficking of Rickettsia from the midgut into the salivary gland during feeding. While studies that assessed the impact of blood meal on the tick microbiome are limited to developmental stages, it is interesting to note that our studies corroborate a great number of these studies, albeit in life stages. Zolnik and colleagues (2016) observed a significant reduction in alpha-diversity throughout the feeding period of nymphal ticks compared with the unfed nymphs of black-legged ticks, and Zhang et al. (2014) also reported an altered microbial composition following a prolonged blood meal in Ixodes persulcatus ticks. A similar observation was reported in the malaria mosquito, in which Muturi et a