Genetic structure and Rickettsia infection rates in Ixodes ovatus and Haemaphysalis flava ticks across different altitudes

Ixodid ticks, such as Ixodes ovatus and Haemaphysalis flava, are important vectors of tick-borne diseases in Japan, such as Japanese spotted fever caused by Rickettsia japonica. This study describes the Rickettsia infection rates influenced by the population genetic structure of I.ovatus and H. flava along an altitudinal gradient. A total of 346 adult I. ovatus and 243 H. flava were analyzed for the presence of Rickettsia by nested PCR targeting the 17kDA, gltA, rOmpA, and rOmpB genes. The population genetic structure was analyzed utilizing the mitochondrial cytochrome oxidase 1 (cox1) marker. The Rickettsia infection rates were 13.26% in I. ovatus and 6.17% in H. flava. For I. ovatus, the global FST value revealed significant genetic differentiation among the different populations, whereas H. flava showed non-significant genetic differentiation. The cox1 I. ovatus cluster dendrogram showed two cluster groups, while the haplotype network and phylogenetic tree showed three genetic groups. A significant difference was observed in Rickettsia infection rates and mean altitude per group between the two cluster groups and the three genetic groups identified within I. ovatus. No significant differences were found in the mean altitude or Rickettsia infection rates of H. flava. Our results suggest a potential correlation between the low gene flow in I. ovatus populations and the spatially heterogeneous Rickettsia infection rates observed along the altitudinal gradient. This information can be used in understanding the relationship between the tick vector, its pathogen, and environmental factors, such as altitude, and for the control of tick-borne diseases in Japan.

articles from PLOS ONE for specific examples.
This statement is required for submission and will appear in the published article if the submission is accepted.Please make sure it is accurate.

Unfunded studies
Enter: The author(s) received no specific funding for this work.Tick here if the URLs/accession numbers/DOIs will be available only after acceptance of the manuscript for publication so that we can ensure their inclusion before publication.
Powered by Editorial Manager® and ProduXion Manager® from Aries Systems Corporation

Introduction
Tick-borne diseases are a significant public health concern in Japan and are transmitted by a diverse range of tick species, such as Ixodes ovatus [1] that potentially transmit Borrelia sp.
causing Lyme disease [2] and Haemaphysalis flava [3], which transmits Rickettsia japonica and is a suspected vector of severe fever with thrombocytopenia syndrome virus [4][5][6].Since ticks are small, their dispersal is closely related to the movement of their hosts [7][8].The dynamics of tickborne pathogens are influenced by the habitat distribution and dispersal behaviors of vectors and hosts along environmental gradients [9].Therefore, understanding the complex interaction between these factors is important in understanding the spread of tick-borne diseases in Japan.
Tick population genetic analysis provides data that help identify the dispersal pattern of ticks based on gene flow between local populations [10].Tick dispersal of many three-host Ixodid ticks depends on the movements of their vertebrate hosts, which influences each tick's potential to spread its pathogen [7][8].For example, contrasting patterns in the population genetic structures of I. ovatus and H. flava in the Niigata Prefecture of Japan suggest that host mobility during the immature stages of tick development may influence the genetic structure of adult ticks by affecting survivability into their adult stages [11][12].I. ovatus populations had greater genetic divergence possibly due to the limited dispersal of their small mammalian hosts during the immature development stage; H. flava populations showed a more homogenized structure possibly due to the larger mobility of their large mammalian hosts and avian-mediated dispersal [11].Other studies have also revealed low gene flow in ticks with low-mobility hosts (e.g., small mammals) and higher gene flow in ticks with highly mobile hosts (e.g., large mammals and birds) [10; 13-16].
The spatial distribution and movement of the vector (i.e., ticks) may determine the spatial distribution of the pathogen's (i.e., Rickettsia) infection rate [17].Previous studies have shown that the pathogen infection rate can be influenced by many factors, such as the vector's genetic diversity, gene flow, and spatial structure [5-6; 11; 18-25].For example, previous studies have shown that strong gene flow between local vector populations tends to reduce the spatial heterogeneity of pathogen infection rates between populations [26][27].Thus, the spatial distribution and movement of the vector may affect the spatial distribution of the pathogen.To our knowledge, no previous studies have examined the relationship between the spatial heterogeneity of Rickettsia infection rates and population genetic structure in ticks.
Environmental factors may relate to the population genetic structure of ticks [28][29], with limited gene flow increasing genetic variation between populations along an altitudinal gradient, as reported in several studies on other species [30][31][32].Thus, in this study, we expected to see a highly divergent population genetic structure in I. ovatus along an altitudinal gradient due to the limited movement of their mammalian hosts along that gradient; whereas H. flava should show a less divergent structure along the altitudinal gradient due to the higher mobility of its hosts.To our knowledge, no studies have focused on tick gene flow along an altitudinal gradient.
In the study by [11], no significant influence of environmental factors, including altitude,

Materials and Method
Published data of [6;11] In this study, we used cox1 sequence data from [11] for I. ovatus (n = 307) and H. flava Prefecture.The 38 sites surveyed in the previous study by [6] include the 30 sites that were used in this present study (S1 Table ).Please refer to [6;11] for more information about the study sites, collection, sampling identification, DNA extraction, PCR amplification, and sequencing methods used in each respective study.
To strengthen our analysis, we also added data collected from April to October 2018 from additional individuals (n=62; ) sampled at 30 sites across the Niigata Prefecture, including two sites not previously sampled by [11].At these sites, ticks were collected 2-14 times from six core sites among the 30 sites, while ticks were collected once at the remaining sites.The altitude at each site ranged from 8 to 1402 meters above sea level (m.a.s.l.), with a mean altitude of 348 m.a.s.l.

Unpublished data from 2018 collection
Ticks collected were stored at 4 C in microcentrifuge tubes with 70% ethanol.Each collected tick was morphologically identified using a stereo microscope following the identification keys of [36].Genomic DNA was extracted using Isogenome DNA extraction kits (Nippon Gene Co. Ltd. Tokyo, Japan) following the manufacturer's recommended protocol.
In this study, we combined previously published data from [6] with our newly collected data to calculate the Rickettsia infection rate, which is the percentage of Rickettsia-infected ticks from each obtained population.We analyzed the obtained tick DNA for spotted fever group Rickettsia (SFGR) detection and host identification and amplified the mitochondrial gene cox1 for population genetic analysis.We performed nested PCR targeting the following genes for the detection of Rickettsia sp.: 17-kDA antigen gene (17-kDA); citrate synthase gene (gltA); spotted fever group (SFG)-specific outer membrane protein A gene (rOmpA); and outer membrane protein B gene (rOmpB) as described and analyzed in [6].Briefly, we first amplified the 17-kDa protein.If the results were positive, then PCR was performed to target gltA.Samples that were positive with both 17-kDA and gltA were regarded as positive for SFGR and a nested PCR was performed to target the rOmpA and rOmpB gene, samples that are positive for 17-Kda, gltA, rOmpA, and rOmpB genes were sequenced to identify the Rickettsia species.The amplified PCR products were purified using AMPure XP (Beckman Coulter Co., Japan) and sequenced using the Big Dye Terminator Cycle Sequence Kit (Thermo Fisher Scientific).
The cox1 mitochondrial gene was amplified by PCR for cox1 (658 base pairs) using the primer pairs LCO-1490 (5′-AAACTTCAGGGTGACCAAAAAATCA-3') and HCO1-2198 (5′-AAACTTCAGGGTGACCAAAAAATCA -3) for phylogenetic analysis and tick species identification [37].The PCR amplification profile included an initial denaturation of 94 C for 2 min, followed by denaturation at 94 °C for 30 s, then annealing at 38 °C for 30 s, followed by an extension of 72 °C for 1 min for 30 cycles, and a final extension of 72 °C for 10 min.The obtained PCR products were purified using the QIAquick 96 PCR Purification Kit (Qiagen, Germany) following the manufacturer's instructions and were sequenced by Eurofin Genomics, Inc. (Tokyo, Japan).
Each forward and reverse read was assembled using CodonCode Aligner version 1.2.4 software (https://www.codoncode.com/aligner/).Low-quality bases were removed in the aligned sequences, and no ambiguous bases were detected.We used the MAFFT alignment online program (https://mafft.cbrc.jp/alignment/server/) to perform multiple alignments using the default settings.The sequences were checked for similarities with the deposited reference sequences from GenBank for sequence quality and tick species confirmation using BLAST (https://blast.ncbi.nlm.nih.gov/Blast.cgi).All final aligned sequences were checked in Mesquite version 3.5 [38].The protein-coding genes were translated to amino acids to confirm the absence of stop codons.

Population genetic analysis
Multiple sites that are in proximity to each other were combined for population genetic analysis if less than eight individuals were obtained per site, which resulted in 8 populations labeled A to H (S1 Table ).Three sites were excluded from the population genetic analysis because of the limited number of obtained individuals (<8) and the lack of a nearby site to combine into a single population.
The final cox1 sequences of the tick species: I. ovatus and H. flava were individually analyzed using DNASp version 6.12.03 to determine the haplotype diversity per species [39].The level of genetic divergence between each population was quantified per species using global FST.
The genetic relationship between the I. ovatus populations was visualized using the unweighted pair group with the arithmetic mean (UPGMA) cluster method using the APE package [41] for the RStudio software (R Development Core Team, 2016).A cluster dendrogram was created using pairwise FST values genetic distance matrix from GenAlEx.

Haplotype network and phylogenetic analyses
We constructed a haplotype network analysis using PopART program version 1.7 (http://popart.otago.ac.nz/index.shtml)on cox1 I. ovatus and H. flava sequences to assess haplotype relationships using the median-joining network algorithm [42].Briefly, we constructed a Bayesian phylogenetic tree of cox1 haplotypes for I. ovatus and H. flava, respectively, using BEAST version 1.10.14[43].We used the Hasegawa-Kishino-Yano substitution model with estimated base frequencies.We employed a strict clock model and used the coalescent prior as the tree prior.A maximum clade credibility tree was acquired using TreeAnnotator version 1.10.14 using trees from BEAUti version1.10.14 with 90% of the trees as the burn-in.We viewed the constructed maximum clade credibility tree using FigTree version 1.4.4.

Statistical analysis
To determine whether there was a significant difference in the Rickettsia infection rate between haplotype groups for I. ovatus and H. flava, we performed a z-score test at p < 0.05.To determine whether there were differences in the mean altitude between the haplotype groups, we used the Welch t-test at p < 0.05.

Results
The total number of positive (pos) and negative (neg   Rickettsia asiatica in the 17kDA, gltA, and rOmpB markers [6], while an additional 19 adult I. ovatus from the 2018 collection were also positive with R. asiatica.Two haplotypes were found in the rOmpB and 17kDA markers, respectively.One haplotype was found in only one individual (17369).
Two out of the 15 Rickettsia-infected H. flava ticks were found to have the same haplotypes in the 17Kda, gltA,and rOmpA markers, and were identified as Rickettsia sp.(LC461063).The remaining 13 Rickettsia-infected H. flava ticks were identified as Rickettsia sp.Four out of the six Rickettsiainfected I. monospinus were identified as Rickettsia helvetica with the same haplotype found in the 17kDA, gltA, and rOmpB markers.
Based on the population genetic analysis of the cox1 sequences, we found a significant global FST of 0.4154 at p < 0.05 for I. ovatus (Table 1 1).There were 59 and 66 cox1 haplotypes found among the 346 I. ovatus and 243 H. flava individuals, respectively.We found the following number of cox1 haplotypes in the remaining species: I. monospinus (n = 9), I. asanumai (n = 4), I. nipponensis (n = 4), I. persulcatus (n = 3), and The cox1 haplotype network of I. ovatus (Figure 1) revealed four genetic groups, wherein three genetic groups (1, 2, and 3) were distributed along different altitudinal gradients, as shown in ) were distributed at lower altitudes (Figure 2).The Welch t-test revealed a significant difference between the mean altitudes of genetic groups 1 and 3 at p < 0.05 (Table 2); however, no significant difference was observed between genetic groups 1 and 2 or groups 2 and 3. We found a significant difference in the Rickettsia infection rates between I. ovatus genetic groups 1 and 2 based on the z-score test, but no significant difference between groups 1 and 3 or groups 2 and 3 (Table 2).The mean altitude between the Rickettsia-infected (=273.72 m.a.s.l.) and non-infected I. ovatus (=369.61m.a.s.l.) revealed a significant difference based on the Welch t-test at p < 0.05 (Figure 3).
The UPGMA dendrogram of I. ovatus revealed two genetic clusters, 1 and 2 using the genetic distance among the seven populations excluding one population due to the limited number of samples (Figure 4).The z-score test showed a significant difference at p < 0.05 between the Rickettsia detection rates in haplotype groups 1 and 2 (indicated by ab) and in cluster dendrogram groups 1 and 2. The Welch t-test at p < 0.05 revealed a significant difference in the mean altitude of haplotype groups 1 and 3  ).No significant difference between the Rickettsia infection rates of H. flava genetic groups 1 and 2 was observed using the z-score test at p < 0.05 (S2 Table ).

Discussion
Our findings support our hypothesis that a genetically structured tick population, such as I.
ovatus, can cause the Rickettsia infection rate to be spatially heterogenous due to limited gene flow along an altitudinal gradient.The significant global FST estimate of 0.4154 among the I. ovatus populations revealed genetic differentiation between the populations as supported by the occurrence of two genetic clusters in the cluster dendrogram.Our results were consistent with our previous study [11] which suggested that the low mobility of the host species for immature I. ovatus contributed to low gene flow in the tick populations.The low I. ovatus gene flow along the altitudinal gradient might have caused the spatial heterogeneity of Rickettsia infection rates among these populations, which is supported by the significant difference found in Rickettsia infection rates between genetic clusters 1 and 2. A similar pattern was observed in the studies of [44][45][46] which found that the infection rate of Borrelia burgdorferi, the causative agent of Lyme disease, decreased in ticks along the altitudinal gradient.Low gene flow can cause infected and uninfected ticks to have limited opportunities to traverse a wider spatial area thus causing a heterogeneous Rickettsia infection rate [47][48].
We found two genetic groups in the H. flava haplotype network, but no significant difference in the Rickettsia infection rates between the two groups.These results might be due to the high gene flow observed in the H. flava populations, which enable Rickettsia-infected and uninfected H. flava individuals to traverse between the study sites.The high mobility of the large mammalian hosts used by adult H. flava and avian-mediated dispersal during their immature stage probably contributed to their homogenized population genetic structure [11], and the resulting homogenized Rickettsia infection rates.Large mammalian hosts and birds have a wide dispersal range that enables the broader movement of Rickettsia-infected ticks, as observed in previous studies of Amblyomma Americanum [49][50][51], H. flava [11], and I. ricinus [52].Birds are especially good at dispersing over large areas since they can easily traverse landscape barriers such as mountains, fences, glaciers, and oceans that would be difficult for mammals to cross [53].
In addition to host mobility, environmental factors can influence the population genetic structure of ticks [54].Each tick species has preferred environmental conditions that are conducive to completing the tick life cycle, thus influencing the geographical distribution of ticks and the risk areas for tick-borne diseases [55][56].Altitude may influence the population genetic structure of ticks through the effect of altitude on the distribution and abundance of ticks and/or their hosts [57][58][59].Altitudinal differences between populations can affect genetic divergence [60], through, for example, ecological isolation, which causes natural selection against maladapted immigrants and limits gene flow [61][62].For example, organisms adapted to low altitudinal sites that cannot tolerate the lower temperatures at higher altitudes would not survive if they dispersed to those higher altitude sites, thus restricting gene flow across the altitudinal gradient [60].
The difference in mean altitude between I. ovatus cluster groups 1 and 2 might be due to adaptive divergence along the altitudinal gradient.Individuals from cluster group 1 were distributed in higher altitude areas of the northern mountainous area of Niigata Prefecture, while cluster group 2 individuals were found in the lower altitude regions of northern Niigata Prefecture.Populations along an altitudinal gradient are prone to differentiating selection pressures, which result in local adaptation [63].Altitudinal gradients may also cause the varying ambient temperature, precipitation, and humidity levels essential to ticks' development and survival [64].These varying environmental factors may cause ticks to have difficulty dispersing over a wide habitat range.Thus, ticks may need to adapt to extreme habitats, such as extreme altitude or precipitation levels, for survival.For example, Rhipicephalus compositus were found in altitudes of 1000-2500 m, but optimal conditions were between 1200-2500 m [64][65][66].
The different Rickettsia infection rates and altitudinal ranges between the I. ovatus phylogenetic groups may be explained by adaptive evolutionary theory, which states that organisms adjust to new or severe changes in their environment to become better suited to their habitat [67][68].Based on the relationship between the I. ovatus phylogenetic groups and their mean attitudes, I. ovatus might be undergoing local adaptation along the altitudinal gradient due to the higher genetic differentiation between populations as supported by the significant global FST (0.4154) found in I. ovatus.Based on isolation by environment (IBE), genetic differentiation will increase with increased environmental differences independent of geographic distances [33-34;69].When environmental conditions differ, the success of immigration in a new habitat is reduced, which may increase the genetic fixation rate due to a lower chance of outcrossing; thereby enhancing genetic isolation [70].Thus, the lower gene flow along the altitudinal gradient reduced the spatial homogeneity of Rickettsia infection rates among the I. ovatus tick populations, thus causing the different Rickettsia infection rates obtained.
In this study, only a few (n < 20) I. monospinus, I. asanumai, I. nipponensis and H. japonica were collected; however, a high Rickettsia infection rate (31.58%-100%) was found despite these limited numbers.In previous studies, high Rickettsia infection rates despite a few individuals were also observed in I. monospinus in several prefectures in Japan [71][72][73].The extent of tick habitat distribution may vary between species, which may have influenced the widely different Rickettsia infection rates observed among the species in this study.The high Rickettsia infection rate in these ticks is probably due to effective transovarial transmission [74].Furthermore, Rickettsial endosymbionts are inclined to have a high infection rate in some tick species populations because they contribute to the tick microbiome which might imply nutritional symbiosis in the tick species with high Rickettsia infection despite the low number of samples tested [75][76][77].Additionally, the few collected individuals of these species may be habitat specialists that thrive in a narrow habitat range in contrast to I. ovatus and H. flava, which may be habitat generalists capable of thriving in a wide habitat range [78].The narrow habitat range of habitat specialist ticks might have enabled increased interactions between infected and uninfected Rickettsia ticks within the small spatial scale compared to the habitat generalist ticks, thus causing their high Rickettsia infection rates.This conclusion should be considered with caution, since the limited number of individual ticks collected across the Niigata Prefecture may have caused a bias in the estimated Rickettsia infection rates determined in this study.We suggest conducting future research on Rickettsia infection rates in a wider sampling range with an increased number of I. monospinus, I. asanumai, I. nipponensis, and H. japonica individuals to further increase our understanding of the association between Rickettsiae and their tick vectors.
The occurrence of local adaptation in tick populations could affect the future of the tickborne disease landscape.Environmental factors, such as precipitation, temperature, and altitude, have been shown to drive population differentiation in insects, such as Anopheles mosquitoes and Drosophila flies [79][80][81], but studies on the environmental adaptation of ixodid ticks, such as I.
ovatus, and its Rickettsia infection susceptibility, have not yet been performed.[7] suggested that environmental conditions that affect bird hosts can also affect the local adaptation of ticks.Few studies have assessed such local adaptation in multiple organisms with varying dispersal abilities [82][83][84][85][86] and are an area in need of future research.
One of the limitations of this study is the use of one mitochondrial gene cox1, despite this, we were able to determine the relationship between the tick population genetic structure and Rickettsia infection rates as influenced by the altitudinal gradient.The mitochondrial cox1 gene has been widely used for population genetic analysis of many tick species and was proven to be informative in determining the relationship from the subfamily to the population level [87][88][89][90][91].
Mitochondrial genes have a mutation rate that is useful in species-level phylogenetics and can be used for wide geographic ranges however its resolution is not fine enough to study species selection [10].In future studies, we suggest including additional mitochondrial genes and or nuclear genes.
Since ticks are blood-sucking ectoparasites, they directly influence their mammalian hosts and the pathogens they transmit [94][95][96].The interaction between the vector (tick), host, and pathogen (Rickettsia) is essential in understanding and predicting the risk and transmission of tickborne diseases [97].Understanding the genetic structure of ticks can serve as an alternative indicator to infer the potential spread of its pathogen [98].Our study found relationships between (1) the population genetic structure of ticks and the corresponding Rickettsia infection rates, (2) altitude and the population genetic structure of ticks, and (3) altitude and Rickettsia infection rates.We found that host mobility may influence the genetic structure of ixodid ticks.This information can be used to design more effective tick-borne disease control programs that focus on screening and detecting pathogens found in ticks and their mammalian hosts.For example, patterns of disease transmission from ticks with a high genetic divergence and less mobile hosts, such as I. ovatus, are likely due to the movement of infected hosts rather than infected ticks.Thus, screening prospective tick hosts for Rickettsia infection would be more suitable in this example.
We suggest screening the hosts of immature I. ovatus, such as small rodents, instead of screening ticks.The Rickettsia infection rate in tick genetic groups can predict the spread of tick-borne diseases caused by Rickettsia, such as the Japanese spotted fever.We also found that altitude may influence the Rickettsia infection rate of I. ovatus genetic groups.This information can be used to determine the high-risk areas (e.g., lowland, mountains, etc.) of tick-borne diseases along an altitudinal gradient.Genetically structured arthropod vectors, such as ticks, can have different vector competencies, and environmental factors, such as altitudinal gradients, that can influence the vector's ability to acquire, transmit, and maintain the pathogen infection [99].
Funded studies Enter a statement with the following details: Initials of the authors who received each award • Grant numbers awarded to each author • The full name of each funder • URL of each funder website • Did the sponsors or funders play any role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript?• NO -Include this sentence at the end of your statement: The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.• YES -Specify the role(s) played.• * typeset Competing Interests Use the instructions below to enter a competing interest statement for this submission.On behalf of all authors, disclose any competing interests that could be perceived to bias this work-acknowledging all financial support and any other relevant financial or nonfinancial competing interests.This statement is required for submission and will appear in the published article if the submission is accepted.Please make sure it is accurate and that any funding sources listed in your Funding Information later in the submission form are also declared in your Financial Disclosure statement.Format for specific study types Human Subject Research (involving human participants and/or tissue) Give the name of the institutional review board or ethics committee that approved the study • Include the approval number and/or a statement indicating approval of this research • Indicate the form of consent obtained (written/oral) or the reason that consent was not obtained (e.g. the data were analyzed anonymously) • Animal Research (involving vertebrate animals, embryos or tissues) Provide the name of the Institutional Animal Care and Use Committee (IACUC) or other relevant ethics board that reviewed the study protocol, and indicate whether they approved this research or granted a formal waiver of ethical approval • Include an approval number if one was obtained • If the study involved non-human primates, add additional details about animal welfare and steps taken to ameliorate suffering • If anesthesia, euthanasia, or any kind of animal sacrifice is part of the study, include briefly which substances and/or methods were applied • Field Research Include the following details if this study involves the collection of plant, animal, or other materials from a natural setting: Field permit number • Name of the institution or relevant body that granted permission • Data Availability Authors are required to make all data underlying the findings described fully available, without restriction, and from the time of publication.PLOS allows rare exceptions to address legal and ethical concerns.See the PLOS Data Policy and FAQ for detailed information.Yes -all data are fully available without restriction Powered by Editorial Manager® and ProduXion Manager® from Aries Systems Corporation A Data Availability Statement describing where the data can be found is required at submission.Your answers to this question constitute the Data Availability Statement and will be published in the article, if accepted.Important: Stating 'data available on request from the author' is not sufficient.If your data are only available upon request, select 'No' for the first question and explain your exceptional situation in the text box.Do the authors confirm that all data underlying the findings described in their manuscript are fully available without restriction?Describe where the data may be found in full sentences.If you are copying our sample text, replace any instances of XXX with the appropriate details.If the data are held or will be held in a public repository, include URLs, accession numbers or DOIs.If this information will only be available after acceptance, indicate this by ticking the box below.For example: All XXX files are available from the XXX database (accession number(s) XXX, XXX.).• If the data are all contained within the manuscript and/or Supporting Information files, enter the following: All relevant data are within the manuscript and its Supporting Information files.• If neither of these applies but you are able to provide details of access elsewhere, with or without limitations, please do so.For example: Data cannot be shared publicly because of [XXX].Data are available from the XXX Institutional Data Access / Ethics Committee (contact via XXX) for researchers who meet the criteria for access to confidential data.The data underlying the results presented in the study are available from (include the name of the third party • All data generated or analyzed during the study are included in this published article and its supplementary additional files.All the newly generated sequences are available in the Genbank database under the accession numbers _______.Under submission Powered by Editorial Manager® and ProduXion Manager® from Aries Systems Corporation and contact information or URL).This text is appropriate if the data are owned by a third party and authors do not have permission to share the data.• * typeset Additional data availability information: was observed in the genetic structures of I. ovatus and H. flava based on the mantel test, but the study did not use any other robust analytical methods to thoroughly examine the influence of altitude on tick genetic structure.In another study, major spotted fever group Rickettsia (SFGR) prevalence was analyzed in a total of 3,336 immature and adult ticks across the Niigata Prefecture, Japan in the following tick species: Dermacentor taiwanensis, Haemaphysalis flava, Haemaphysalis hystricis, Haemaphysalis longicornis, Haemaphysalis megaspinosa, Ixodes columnae, Ixodes monospinosus, Ixodes nipponensis, Ixodes ovatus, and Ixodes persulcatus[6].Three SFGR species namely Rickettsia asiatica, R. helvetica and R. monacensis were detected in H. flava, Haemaphysalis longicornis, Ixodes monospinus, Ixodes nipponensis, and Ixodes ovatus, no spatial distribution of Rickettsia infection rates was found among the local populations.To our knowledge, no previous studies have considered the influence of environmental factors on tick population genetic structure or the spatial distribution of Rickettsia infection rates along an altitudinal gradient in local tick populations.Here, we determine the relationship between the vector tick population genetic structure and Rickettsia infection rates as influenced by an altitudinal gradient.Based on the isolation by environment (IBE) theory, genetic differentiation increases with environmental variation, regardless of geographic distance[33][34][35]. Thus, we hypothesized that limited gene flow among tick populations along an altitudinal gradient might cause spatially heterogenous Rickettsia infection rates in ticks.The present study used cytochrome oxidase 1 (cox1) mitochondrial gene sequences and Rickettsia-infected and uninfected data for I. ovatus (n = 307) and H. flava (n = 220) individuals collected in 2016-2017 as reported in[6; 11].Additionally, we added new cox1 sequences and Rickettsia-infected/uninfected ticks from I. ovatus (n= 39)and H. flava (n=23) individuals and five other species collected in 2018.

Figure 2 .
Figure 2.These four genetic groups were concordant with the four clusters found in the I. ovatus

Figure 1 .
Figure 1.Median-joining network of the 59 cox1 haplotype sequences of Rickettsia positive and

Figure 2 .
Figure 2. The influence of altitude on the habitat distribution of I. ovatus.The points indicate the

Figure 3 .
Figure 3.The relationship between the mean altitude of Rickettsia positive (n = 46) and negative

Figure 4 .Table 2 .
Figure 4.An unweighted pair group method with the arithmetic mean (UPGMA) dendrogram of I.

Table 1 .
). Summary of the haplotype and Rickettsia infection rates among the 7 ixodid tick species Abbreviations: ns no. of sampling sites; n sample size; nh no of haplotypes; r Rickettsia infection rate per species; *p < 0.05** not enough samples for analysis***tick species identification is based on molecular identification using the cox1 marker and BLAST resultsWe detected SFGR in 78 (12.44%) out of 627 ixodid ticks, with the highest detected in I. ovatus(46/346; 13.29%), H. flava (15/243; 6.17%), and I. monospinus (6/19; 31.58%) as summarized in Table 1.Out of the 46 Rickettsia-infected I. ovatus ticks, 25 displayed a 100% identity match with and cluster dendrogram groups 1 and 2 (indicated by ab ).Haplotype group 4 was not included in the analysis due to its low sample size.
The cox1 haplotype network of H. flava displayed two genetic groups (S2 Figure) consistent with the H. flava phylogenetic tree (S3 Figure).Based on the Welch t-test, we found no significant differences in mean altitude (205 and 165 m.a.s.l.) at p < 0.05 between the two H. flava genetic groups (S3 Table