Sampling and tick species identification

A total of 1855 ticks were collected, including questing ticks (n = 752) and engorged ticks (n = 1083) from livestock cattle and goats. 20 ticks (WQrt) were collected from the long-tailed ground squirrels (Table 1). EDTA-anticoagulated blood specimens were collected from parasitized livestock cattle (ALSKcb, n = 25) and goats (WQgb, n = 79). We also collected blood samples from great gerbils (ALSKrb n = 10) and long-tailed ground squirrels (WQrb n = 29), and 368 (WQhb) local herdsman blood samples were collected.

Among the 1855 ticks examined, we identified three genera and six species. Notably, 1132 ticks, comprising both questing ticks (n = 600) and engorged ticks from livestock goats (n = 512) and the long-tailed ground squirrels (n = 20), were collected in the WQ area. These ticks were conclusively identified as Dermacentor nuttalli. In our study, we observed the presence of two or three different tick species on each livestock animal. Among the 571 engorged ticks examined, we found that 49.91% were Hyalomma asiaticum (n = 285), 40.28% were Dermacentor marginatus (n = 230), and 8.07% were Hyalomma detritum (n = 46), while only a small number of ticks belonged to Dermacentor nuttalli (n = 5) and Rhipicephalus sanguineus (n = 2). Three ticks could not be identified. Of the 152 questing ticks, 147 were identified as Hyalomma asiaticum, 4 were identified as Dermacentor nuttalli, and 1 was identified as Hyalomma marginatum (n = 1) (for details see Table 1).

16S partial region sequencing

Out of the 111 tick pools comprising a total of 1855 ticks, we obtained a total of 4,770,842 quality-filtered reads. On average, each sample yielded approximately 42,981 reads, with a standard deviation of 22,548 and a range spanning from 19,265 to 65,352. The quality of these reads, indicated by an average quality score of 36, was deemed satisfactory. It is noteworthy that all libraries generated from tick samples exhibited ample sequencing depth for subsequent analysis. This is evident from the rarefaction curves plotting the number of observed OTUs, which reached a plateau when considering a sequencing depth ranging from 5000 to 10,000 sequences (Additional file 2). This plateau suggests that our samples had achieved sufficient coverage, validating their suitability for further analysis.

Microbiome profile

At the phylum level, tick pools in Alataw City carried more diverse bacteria. that those in Wenquan County The shared bacteria between ticks from the two locations belong to 9 phyla: Acidobacteriota, Actinobacteriota, Bacteroidota, Deinococcota, Firmicutes, Gemmatimonadota, Proteobacteria, Fusobacteriota and Verrucomicrobiota, and there were 11 exclusive phyla in ALSK ticks but only one exclusive phylum in WQ ticks (Additional file 3). Regarding abundance differences among groups, ticks from the two environments showed bacterial differences in 14 phyla out of all 21 phyla (Additional file 4). Specifically, questing ticks from ALSK carried more Acidobacteriota (W = 325, P = 0.0117), Actinobacteriota (W = 447, P < 0.0001), Bacteroidota (W = 426.5, P = 0.0005), Cyanobacteria (W = 297.5, P = 0.0315), Gemmatimonadota (W = 297.5, P = 0.0315) and Patescibacteria (W = 420, P < 0.0001) than WQ ticks. In addition, engorged ones from ALSK carried more Actinobacteriota (W = 807.5, P < 0.0001), Bacteroidota (W = 775, P < 0.0001), Fusobacteriota (W = 516, P = 0.0367) and Patescibacteria (W = 577.5, P = 0.0016) but less Deinococcota (W = 188, P = 0.0004) and Proteobacteria (W = 269, P = 0.0215). Regarding questing ticks and engorged ticks, engorged tick pools carried more diverse bacteria than questing tick pools. There were 13 shared bacterial phyla and no exclusive phyla in questing ticks but 8 exclusive phyla in engorged ticks (Additional file 5). Furthermore, questing ticks from ALSK carried more Acidobacteriota (W = 363, P = 0.0291) and Deinococcota (W = 493, P = 0.0003) but less Fusobacteriota (W = 217.5, P = 0.0265) and Verrucomicrobiota (W = 232.5, P = 0.0494) than engorged ticks. The questing ticks from WQ carried more Deinococcota (W = 525, P = 0.0071) but less Bacteroidota (W = 232.5, P = 0.0059) and Firmicutes (W = 209, P = 0.0067) than engorged tick pools.

At the genus level, the abundance table consisted of 495 genera. Similar to the phylum level, tick pools in Alataw City carried more diverse bacteria than those in WQ County, and engorged tick pools carried more diverse bacteria than questing tick pools (Additional file 6–7). For the 4 groups, there were 203.5 ± 147.63 (mean ± standard deviation) genera for each. Thirty-nine (7.88%) were shared in all groups (Fig. 1A). The shared genera occupied 86.21% of the total abundance. The engorged tick pools in Alataw City had the most exclusive genera (238, 78.81% in all exclusive genera), and 17 among these (7.15%) were from exclusive phyla of this group. While the engorged tick pools in WQ County had the least exclusive genera (3, 0.99% in all exclusive genera), one of these belongs to the only exclusive phylum of this group. In general, the genera with the top 10 relative abundances were Coxiella (22.99%), Rickettsia (17.90%), Acinetobacter (12.23%), Francisella (9.65%), Anaplasma (5.94%), Psychrobacter (5.25%), Thermus (3.19%), Pseudomonas (2.92%), Anoxybacillus (2.39%) and Rickettsiales bacterium Ac37b (1.74%) in all 111 samples. These high abundance genera are included in 39 shared genera except the genus Anaplasma, and 8 in 10 belong to the phylum Proteobacteria, which showed that shared genera formed a major component of the tick-carrying bacterial microbiome. The diversity of the bacteria within individual samples (alpha diversity) shows richer in engorged tick pools in ALSK (P = 5.4e-6) (Fig. 1B). However, the alpha diversity of bacteria in WQ ticks displayed no difference between questing and engorged tick pools (P = 0.37). The beta diversity between groups was assessed using the Bray–Curtis distance metrics and visualized using NMDS (Fig. 1C). In beta diversity observations, the genera of the 4 groups were significantly distinct from each other (ANOSIM, p < 0.01). LEfSe determined genera that had different abundances between groups (Fig. 1D-E). It showed 7 differential genera between questing tick-carrying bacteria from two places and 22 between engorged tick-carrying bacteria (Fig. 1D). Among these, the top differential genera according to LDA score were intracellular bacteria. However, the most differential genera were environmental bacteria. Two common differential genera that occurred in both places were Coxiella and Thermus, which are more abundant in questing ticks. For different engorged statuses, there were 16 and 6 differential genera between questing ticks and engorged ticks from ALSK and WQ (Fig. 1E). Similarly, intracellular bacteria were the top differential genera, and the most differential genera were environmental bacteria.

Fig. 1

Profiles of the bacterial microbiome from ticks at the genus level. (A) Venn diagram of bacterial genera for 4 groups of ticks: questing ticks from ALSK (ALSKqt), engorged ticks from ALSK (ALSKet), questing ticks from WQ (WQqt) and engorged ticks from WQ (WQet); (B) Alpha diversity of bacterial genera for the 4 groups of ticks above; (C) Beta diversity of bacterial genera for the 4 groups of ticks above; (D) Bacterial genera with different abundances between environments; (E) Bacterial genera with different abundances between engorgement statuses

The predicted function of microbial communities

Four groups were used to detect differences in the functions of microbial communities to inspect the relation of environment and parasite statuses on observed microbial communities. First, 40 samples of the ALSKet group were divided into two subsets, both of which contain 20 samples but have different tick compositions. In the two subsets, one consisted of 0.69% D. nuttalli, 61.51% Hy. asiaticum, 27.84% D. marginatum, 9.62% Hy. detritum and 0.34% R. sanguineus (n = 291), and the other subset consisted of 1.08% D. nuttalli, 38.27% Hy. asiaticum, 53.79% D. marginatum, 6.5% Hy. detritum and 0.36% R. sanguineus (n = 277). Statistical analysis showed no significant change in the functional abundance of microbial communities (Additional file 8), suggesting that tick species have a limited influence on the formation of microbial functions in ticks.

For ALSK and WQ tick pools, there were 222 differential function pathways in bacteria (Additional file 9). When engorged statuses were considered, there were 97 detected in questing ticks and 237 detected in engorged ticks (Additional file 10–11). This means that there are many metabolic differences among tick-carrying bacteria in different environments. By inspecting microbial function abundance between ALSKqt and ALSKet, WQqt and WQet, 106 and 83 differential function pathways were discovered (corrected Ps < 0.05, Additional file 12–13), respectively. Among these pathways, there were 9 shared differential function pathways (Fig. 2), including chorismate biosynthesis I, chorismate biosynthesis from 3-dehydroquinate, the superpathway of β-D-glucuronide and D-glucuronate degradation, glycogen degradation I (bacterial), GDP-mannose biosynthesis, D-fructuronate degradation, the superpathway of aromatic amino acid biosynthesis and tRNA processing, which were all less abundant in engorged ticks, and one pathway of biotin biosynthesis, which was more abundant in engorged ticks.

Fig. 2

The 9 shared differential metabolic pathways between questing ticks and engorged ticks from two environments. Bar plot with extended error bar. Red and blue objects represent groups of questing ticks and engorged ticks, respectively

Difference between microbial co-occurrence networks (CoNets) of questing and engorged ticks

Microbial co-occurrence networks representing different groups were constructed, and the degree parameters of CoNets are displayed in Fig. 3. The vertices and average degrees of the CoNets of tick-carrying bacteria on questing and engorged ticks of ALSK are less than those of WQ. Similarly, the vertices and average degrees of the CoNets of tick-carrying bacteria on engorged ticks are less than those on questing ticks. These results suggest that the CoNets of tick-carrying bacteria were affected by the environments where ticks grow and in the engorgement status of the ticks.

Fig. 3

Comparison between microbial co-occurrence networks for tick-carrying bacteria at different levels. (A) microbial co-occurrence networks for ticks from ALSK and WQ, left for group questing ticks, right for engorged ticks; (B) microbial co-occurrence networks for questing and engorged ticks, left for group ALSK, right for group WQ; (C) dominated genera which account for 50% total degree of net of each group, orange lines indicated shared key genera between two groups. Each ball represents a genus, and the size of the ball indicates the degree of this genus. Red and green balls indicate dominant genera of different groups, which account for 50% of the total degree of net. Orange balls indicate shared dominant genera for both groups. Links represent co-occurrence relationships between genera, and their colours show which group displayed this relationship. Values in the dashed box indicate the average degree of the color corresponding group in the cooccurrence network

All groups were composed of different taxonomic profiles and were dominated by the phylum Proteobacteria (Additional file 14). For the abundance of bacteria in CoNets, ticks from two places showed bacterial differences in 5 phyla out of all 10 phyla. Specifically, questing ticks from ALSK carried more Actinobacteriota (W = 417, P < 0.0001), Bacteroidota (W = 393, P = 0.0004), Firmicutes (W = 329, P = 0.0391), and Patescibacteria (W = 391, P < 0.0001) and less Proteobacteria (W = 134, P = 0.0176) than WQ ticks. In addition, engorged ones from ALSK carried more Actinobacteriota (W = 749.5, P < 0.0001), Bacteroidota (W = 709, P < 0.0001), and Patescibacteria (W = 540, P = 0.0018) and less Deinococcota (W = 172, P = 0.0005) and Proteobacteria (W = 112, P < 0.0001). Additionally, questing ticks from ALSK carried more Deinococcota (W = 458, P = 0.0002) and Proteobacteria (W = 418, P = 0.0028) than engorged ticks. Questing ticks from WQ carried more Deinococcota (W = 487, P = 0.0079) but less Bacteroidota (W = 463, P = 0.0224) and Firmicutes (W = 188, P = 0.0059) than engorged ticks. For dominant genera that account for 50% of the total degree of net, few of the same genera were detected between different places or different engorgement statuses. These results showed that the tick-carrying bacteria included in CoNets between the two environments had significant differences in taxonomy and abundance. Meanwhile, in same environment, no pattern for bacterial differences between tick-questing and engorged ticks was detected.

Pathogen-related genera correlated to symbionts Coxiella and Francisella

After extracting local correlations for pathogen-related genera and 2 high abundance symbionts (Coxiella and Francisella) from the SparCC correlation matrix, a correlation network was summarized, as shown in Fig. 4. There were 63 connected to Coxiella in a total of 85 genera and 28 connected to Francisella. Among these, 44 environmental genera, 9 pathogen-related genera, 6 gut-related genera, 2 symbiont genera and 2 undetermined genera were connected to Coxiella. In addition, 23 environmental genera, 2 pathogen-related genera, 1 gut-related genus and 2 symbiont genera were connected to Francisella. The remaining 21 genera connected to Coxiella were mediated by one genus, but Francisella still had 10 genera that were not connected to Coxiella. Pathogen-related genera that were positively correlated with Coxiella were Anaplasma, Mycoplasma, Rickettsia, Roseomanos, Spiroplasma and Ehrlichia. Arsenophonus, Moraxella, Samonella, and Trueperella were negatively correlated with Coxiella. In addition, Ehrlichia and Mycoplasma were directly correlated with Francisella, both positively, and 5 indirectly correlated with Coxiella. These results revealed that the main correlations of symbionts were environmental genera, and some important pathogen-related genera correlated to symbionts where these genera showed different correlations with each symbiont. Additionally, no stable pattern of relationships between the two symbionts exists in different environments or engorged statuses.

Fig. 4

Genera correlated with endosymbionts Coxiella and Francisella. Each ball represents a genus, and the size of the ball indicates the degree of this genus. Orange, red, purple, blue and green balls represent endosymbionts, pathogen-related genera, environmental bacteria, gut bacteria and unknown genera, respectively. Names of endosymbionts and pathogen-related genera are displayed

Multiple linear regression model (MLRM) revealed a quantitative relation between Coxiella and its correlated genera

For the two symbionts, only the dataset of Coxiella passed the evaluation of MLRM in gvlma, which suggested that the abundance change between Francisella and its correlated genera did not satisfy a linear relation. The MLRM formula for the abundance of Coxiella comprised 7 intracellular parasites, 6 gut-related genera and 28 environmental genera (Fig. 5). The abundance of these 41 genera interpreted a 99.97% change in Coxiella abundance. For the test dataset, the R squared between the prediction of the model and observation was 0.7280 after the negative predicted value was corrected to 0. The top 3 coefficients with the most weights, which belong to the genera Phenylobacterium, Moraxellla and Mycoplasma, were all negative. The top 3 positive coefficients with the most weights were possessed by the genera Herminiimonas, UCG-009 and Cloacibacterium. The weight proportions of the 6 genera in the formula were 21.96%, 17.92%, 7.43%, 7.18%, 6.34% and 3.46%. These results revealed a quantitative relation between Coxiella and its correlated genera, and environmental genera were the most important factors influencing the abundance of symbionts of Coxiella.

Fig. 5

The 41 genera and their coefficients that have abundance correlations with Coxiella revealed by a multiple linear regression model. Negative and positive coefficients indicate negative and positive correlations with Coxiella, respectively. The absolute value of the coefficient indicates the weight of the correlation with Coxiella

Pathogens and symbionts in ticks

A total of 8 microorganisms with different carrier rates were detected in samples from 10 sources, including pathogens and tick symbionts: Ehrilichia spp. (0.42%, 10/2336), Bartonella spp. (0.77%, 18/2336), Babesia spp. (5.77%, 135/2336), Borrelia spp. (7.49%, 175/2336), Francisella spp. (19.34%, 452/2336), Anaplasma spp. (21.48%, 502/2336), Rickettsia spp. (30.52%, 713/2336), and Coxiella spp. (52.48%, 1226/2336). In addition, 86.36% of the tick samples carried one or more of the detected bacterial microorganisms, and the carrying rates from low to high were 20% (4/20) for WQrt, 75.66% (432/571) for ALSKet from livestock cattle, 87.33% (524/600) for WQqt, 95.70% (490/512) for WQet from goats and 100% (152/152) for ALSKqt. Moreover, 43.77% (812/1855) of tick samples carried bacterial pathogens, and nearly 18% of blood samples carried one or two detected pathogens. More specifically, the bacterial pathogen carrying rates of the blood samples from low to high were 0% (0/29) for the WQrb, 0% (0/10) for the ALSKrb, 1.6% (6/368) for the WQhb, 46% (13/28) for the ALSKcb and 94% (74/79) for the WQgb (for details see Table 2).

Table 2 Detection the statistics of individual microbes within each sample type Ehrilichia spp.

Sanger sequencing analysis showed that these Ehrilichia sequences were the same and had 100% homology with 16S rRNA sequence of Ehrilichia (JX402605.1) found in Hy. asiaticum ticks in Xinjiang Province (northwestern China), the Ehrilichia 16S rRNA sequence KX987325.1 carried by Boophilus microplus ticks in Hubei Province (central China) and the Ehrilichia 16S rRNA sequence KY046298.1 carried by Rhipicephalus microplus ticks in Malaysian [44]. By aligning the sequences with the Basic Local Alignment Search Tool (BLAST), it was found that this Ehrilichia sequence had 100% homology with MT875368.1 (Candidatus Ehrlichia Hainanensis), KJ513197.1 (Ehrlichia canis), MZ733621.1 (Candidatus Ehrlichia pampeana) and MT738235.1 (Ehrlichia ruminantium). These results suggest that this species may be a pathogen that is prevalent in humans, livestock, and rodents.

Bartonella spp.

Bartonella spp. were detected only in tick samples collected from Wenquan County. The results of Sanger sequencing showed that Bartonella spp. had 100% homology with Bartonella henselae. Bartonella henselae is a zoonotic pathogen causing neurological disease [45]. However, this pathogen was not detected in the human blood and goat blood samples collected from Wenquan County. Since the felines and canines in Wenquan County were not included in the scope of this investigation, the situation was still unknown.

Babesia spp.

Babesia spp. was detected in a total of 136 cases in ticks, livestock and human blood samples from Wenquan County and Alataw City. By sequencing and sequence alignment, it was found that the sequence of the Babesia spp. in questing ticks and parasitic ticks had 100% homology with LC553515.1, which is a Babesia strain closely related to the zoonotic Babesia gibsoni [46].

Borrelia spp.

Borrelia spp. were detected in a total of 175 ticks. The results of Sanger sequencing showed that all Borrelia spp. had the same sequence, which was 100% homology with the NR_170496.1 sequence of Borrelia maritima. B. maritima was a new species of the Borrelia burgdorferi sensu lato complex [47, 48].

Francisella spp.

Francisella spp. was detected in a total of 452 samples (5 cases in human blood samples and others in tick samples). The results of Sanger sequencing showed that 4 cases of Alataw City engorged ticks from livestock cattle carried Francisella tularensis subsp. (CP009653.1), and the others carried Francisella-like symbionts (KX852466), which is a tick symbiont that is transmitted vertically through the maternal line. These results suggest that Francisella-like symbionts are stable in the population of ticks in Alataw City.

Anaplasma spp.

Anaplasma spp. was detected in a total of 502 samples. Two species were identified by Sanger sequencing and sequence alignment, including Anaplasma ovis (NZ_CP015994.1), which was found in 416 cases, and Anaplasma phagocytophilum (NC_021880.1), which was found in 86 cases. Anaplasma ovis is the Anaplasma (Rickettsiales: Anaplasmataceae) of susceptible goat. The positive rate of goat blood samples collected from Wenquan County was 93.67% (74/79), and the positive rate of engorged ticks from goat, was 63.08% (323/512), suggesting that the goats that are raised in Wenquan County should be further tested to prevent the spread of pathogens. Moreover, Anaplasma phagocytophilum, an intracellular parasite that can cause human granulocytic anaplasmosis, was detected in tick samples from both places [49], suggesting that the prevention and control of A. phagocytophilum infection in both Wenquan County and Alataw City need to be reinforced.

Rickettsia spp.

Rickettsia spp. was detected in a total of 713 samples. Among those cases, only 1 case was from cattle blood samples in Alataw City, and the rest were from tick samples. Six species were identified by next-generation sequencing, and the sequence alignment results indicated that 5 of the six species were zoonotic pathogens, including Rickettsia raoultii (229/713) [50], Rickettsia sibirica (53/713) [50], Rickettsia sp. Hme_HirooL009 (4/713) (LC544134.1), uncultured Rickettsia sp. (KM587631.1; KM587632.1; KM587633.1; MT434895.1; KX591658.1) (111/713) [51] and Rickettsia sp. BJ-90 (2/713) [52]. Rickettsia sp. BJ-90 is a pathogen identified as belonging to the Rickettsia sibirica clade by evolutionary analysis. Taken together, the carrying rate of the Rickettsia pathogen in ticks was 55.96%. The Rickettsia species of the remaining 314 samples could not be determined, accounting for 44% of the total sequences.

Coxiella spp.

A total of 1226 Coxiella spp. infection cases were detected, and all these cases were ticks collected from the two counties. All Coxiella spp. had the same sequence, which was homologous to the sequence of uncultured Coxiella sp. clone Dx-56 and that of JX432012.1, indicating that they were tick Coxiella-like symbionts. The carrying rates of Coxiella spp. in ticks collected were 80.83% (485/600) for WQqt, 76.75% (393/512) for WQet from goats, 10% (2/20) for WQrt, 42.91% (245/570) for ALSKet from livestock cattle, and 66.45% (101/152) for ALSKqt. These results suggest that the Coxiella-like symbionts are widely distributed and stable in the tick populations in the two regions. According to the number of detected microorganisms, most of the pathogen-carrying blood samples only carried one detected microorganism, and only one case was detected in the Wenquan goat blood sample carrying both Ehrilichia spp. and Anaplasma ovis. However, tick samples often carry up to 6 pathogens (Fig. 6).

Fig. 6

The proportion of the different coexist number of microorganisms in different sample types. Different color column represents the different coexist number of microorganisms in that sample type. WQqt: questing ticks from WQ; WQet: engorged ticks from WQ; WQgb: goat blood from WQ; WQst: squirrel’s ticks from WQ; WQsb: squirrel blood from WQ; WQhb: human blood from WQ; ALSKqt: questing ticks from ALSK; ALSKet: engorged ticks from ALSK; ALSKcb: cattle blood from ALSK; ALSKgb: gerbil blood from ALSK