Climate change influences on the geographic distributional potential of the spotted fever vectors Amblyomma maculatum and Dermacentor andersoni

Abstract

Amblyomma maculatum (Gulf Coast tick), and Dermacentor andersoni (Rocky Mountain wood tick) are two North American ticks that transmit spotted fevers associated Rickettsia. Amblyomma maculatum transmits Rickettsia parkeri and Francisella tularensis, while D. andersoni transmits R. rickettsii, Anaplasma marginale, Coltivirus (Colorado tick fever virus), and F. tularensis. Increases in temperature causes mild winters and more extreme dry periods during summers, which will affect tick populations in unknown ways. Here, we used ecological niche modeling (ENM) to assess the potential geographic distributions of these two medically important vector species in North America under current condition and then transfer those models to the future under different future climate scenarios with special interest in highlighting new potential expansion areas. Current model predictions for A. maculatum showed suitable areas across the southern and Midwest United States, and east coast, western and southern Mexico. For D. andersoni, our models showed broad suitable areas across northwestern United States. New potential for range expansions was anticipated for both tick species northward in response to climate change, extending across the Midwest and New England for A. maculatum, and still farther north into Canada for D. andersoni.

Introduction

Beside the tick Dermacentor variabilis, Amblyomma maculatum (Gulf Coast tick), and D. andersoni (Rocky Mountain wood tick) are three North American ticks that transmit spotted fever (Boorgula et al., 2020; CDC, 2018; CDC, 2019). Spotted fever rickettsioses (spotted fevers) are a group of bacterial pathogens that cause disease to humans by exposure to infected ticks or mites (CDC, 2019). In the United States, there are several spotted fevers: Rocky Mountain spotted fever (RMSF), which is the most documented spotted fever, caused by Rickettsia rickettsii; R. parkeri rickettsiosis caused by Rickettsia parkeri; rickettsialpox caused by R. akari; and Pacific Coast tick fever caused by R. philippi (CDC, 2019). Amblyomma maculatum transmits R. parkeri, and Francisella tularensis which cause diseases in humans, and Hepatozoon americanum, which causes health problems in dogs (Sonenshine, 2018). Dermacentor andersoni transmits R. rickettsii, Anaplasma marginale, Coltivirus (Colorado tick fever virus), and Francisella tularensis (Alkishe, Raghavan & Peterson, 2021; Dantas-Torres, Chomel & Otranto, 2012).

Amblyomma maculatum and D. andersoni have different geographic distributions: A. maculatum occurs throughout the southern states of the Gulf Coast and Mid-Atlantic states (Cumbie et al., 2020), whereas D. andersoni occurs throughout the Rocky Mountain region, Nevada, California, and southwestern Canada (British Columbia, Alberta, and Saskatchewan; Animal Diversity Web (ADW), 2021). Those different geographic ranges are associated with different climate conditions for each tick species (Fig. 1). Amblyomma maculatum is found in different months in different states with fall and winter considered as low-activity seasons for this species (Nadolny & Gaff, 2018), whereas D. andersoni is found in hot and dry areas in summer (Wilkinson, 1967).

Figure 1: Occurrence points and calibration areas for Amblyomma maculatum (blue dots and blue buffer) and Dermacentor andersoni (red dots and red buffer) in geographic and environmental space. DOI: 10.7717/peerj.13279/fig-1

Climate warming is warming North America dramatically. Mean global temperature has increased more than 1°C owing to anthropogenic greenhouse gas emissions (Djalante, 2019). This increasing temperature has caused mild winters with increasing rain more than snow during winter, and more extreme drier periods during summers (Wuebbles et al., 2017). Increasing temperature can also affect vector disease survival, abundance, and activity as well as transmission dynamics, re-emergence of vector-borne diseases, and geographic expansions (Rocklöv & Dubrow, 2020).

Here, we used ecological niche modeling (ENM) to assess the geographic potential of these two medically important vectors of diseases in North America under current conditions and then transfer those models to the future under different scenarios, with special interest in highlighting potential range new expansion areas. We also assess the model uncertainty for projected future models to highlight areas with high versus low confidence of geographic expansions.

Methods Data preparation

We obtained totals of 255 and 586 occurrence points for A. maculatum and D. andersoni, respectively. Those data were obtained from various sources: Global Biodiversity Information Facility (GBIF; http://www.gbif.org), VectorMap (http://vectormap.si.edu/), and BISON (https://bison.usgs.gov) (sources summarized in Supplemental Information 1). We followed Cobos et al. (2018) in cleaning the data to remove errors that clearly fall outside of the known geographic distribution of the species, duplicate records, and localities with missing or meaningless coordinates such as zero degrees latitude and zero degrees longitude (0°N, 0°E), or georeferencing errors (records in ocean and far from coast). We used the spTthin R package to reduce the data spatially based on a 50 km distance filter for several reasons: based on precision of the occurrence points in the area, environmental heterogeneity that present in the area, and to avoid problems with autocorrelation (Aiello-Lammens et al., 2015). In the end, we had 93 and 82 occurrence points for A. maculatum and D. andersoni, respectively. We divided the final occurrence data randomly into two sets: 50% for model calibration and 50% for evaluation steps involved in model calibration. For producing final models, we used the entire cleaned occurrence points.

Delineate the calibration area

The accessible area (M) is the set of places to which the species has had access over relevant time periods, and depends on the dispersal of the species from populations (Barve et al., 2011). Since the movement of tick species is associated with the movement of host species, we assumed ample dispersal abilities for the ticks (Nadolny & Gaff, 2018; Sonenshine, 2018). As such, we created 200 km buffer areas around the known occurrence points for each species (Fig. 1).

Environmental variables

For the current time, bioclimatic variables were downloaded from WorldClim version 1.4, at 10′ spatial resolution (Hijmans et al., 2005) (available at http://www.worldclim.org/). We removed variables 8, 9, 18, and 19 because of known spatial artefacts (Escobar, 2020). The 15 remaining variables were masked to the calibration area (M) for each species. We then used principal component analysis (PCA) to reduce dimensionality and multicollinearity among those variables. After having PCA results, we created 11 sets of environmental variables that represent all possible combinations of the first four principal components to test them with other parameter settings to choose best models during model calibration, following Cobos et al. (2019) (see below).

For future climatic conditions, we used five general circulation models (GCMs) under two representative concentration pathway scenarios (RCP 4.5, and RCP 8.5). Future climate data layers were downloaded from the Climate Change, Agriculture and Food Security (CCAFS) database at 10′ resolution (available at: http://www.ccafs-climate.org/data_spatial_downscaling). GCMs used were (1) National Center for Atmospheric Research (NCAR_CCSM4); (2) Met Office Hadley Centre (HadGEM2); (3) Model for Interdisciplinary Research on Climate (MIROC5); (4) Institut Pierre Simon Laplace (IPSL_CM5A); and (5) Russian Institute for Numerical Mathematics Climate Model Version 4 (INM_CM4). GCM choice was based on frequency of use in other such research applications, and on full availability of scenarios for both RCP scenarios.

Ecological niche modeling and model transfers

The combination of 11 sets of environmental variables, 15 feature classes (all combinations of linear = l, quadratic = q, product = p, hinge = h), and 17 regularization multiplier values (0.1 to 1 at intervals of 0.1, and 2 to 10 at intervals of 1) resulted in 2,805 candidate models for each species. We evaluated candidate models based on statistical significance (partial ROC, P ≤ 0.05; Peterson, Cobos & Jiménez-García, 2018), predictive performance (omission rates, <5%; Anderson, Lew & Peterson, 2003), and a criterion of minimum complexity (Akaike Information Criterion corrected for small sample sizes, AICc; Warren & Seifert, 2011). Specifically, we used differences between particular AICc values and the minimum values (ΔAICc < 2) to select best model parameter settings with which to produce final models.

Final models

For creating final models, we used the complete set of occurrences and the parameterizations selected during model calibration. We created 10 bootstrap replicates, and transferred the models across North America (Mexico, United States and Canada) under current and future scenarios. We calculated medians of all replicate medians from final predictions for each calibration area in which final models were produced to summarize model results. Then, we binarized models using a threshold of allowable omission error rate (E) of 5%, assuming that as a percentage of data may have included errors that misrepresented environments used by the species.

We calculated differences in suitable areas between current and the two future scenarios RCP (4.5, and 8.5). For representing changes of suitable areas, we used the agreement of changes (stable, gain, loss) among the five GCMs per RCP scenario. Simply, for each RCP scenario, we took all projections to future conditions based on distinct GCMs and compared against the current projection, and quantified the agreement of gain and loss of suitable areas, as well as the stability of suitable and unsuitable conditions. All modeling analysis steps were done in R 3.5.1 (R Core Team, 2018) using Maxent 3.4.1 (Phillips et al., 2017), implemented in the kuenm package (Cobos et al., 2019).

Uncertainty in model projections

We used the mobility-oriented parity metric (MOP, considering the nearest 5% of reference cloud) (Owens et al., 2013) to assess strict extrapolation risk. We also calculated variance arising from distinct sources (replicates, parameter settings, GCMs, and RCPs) in our model projections (Peterson, Papeş & Soberón , 2018). Both model variability and strict extrapolation were represented geographically following Owens et al. (2013) and Cobos et al. (2019), respectively.

Results Model calibration results

From among 2,805 candidate models for each of A. maculatum and D. andersoni, 2,728 and 2,554 were significantly better than random expectations, respectively (pROC test, p ≤ 0.05). Of these models, 2,129 and 761 met the omission rate criteria, (i.e., OR ≤ 0.05) respectively. Based on AICc, 55 and two models were selected as best models for A. maculatum and D. andersoni, respectively. For A. maculatum, models performed better with the variables in set 1 (PC1, PC2, PC3, PC4), set 2 (PC1, PC2, PC3), set 3 (PC1, PC2, PC4), and set 6 (PC1, PC2), whereas for D. andersoni variables in sets 4 (PC1, PC3, PC4) and 7 (PC1, PC3).

Current and future potential distribution Amblyomma maculatum

Current model predictions for A. maculatum showed suitable areas across the southern United States (Florida, Georgia, South and North Carolina, Virginia, West Virginia, Maryland, Delaware, Kentucky, Tennessee, Arkansas, Alabama, Mississippi, Louisiana, Oklahoma, and Texas), and in the Midwest (Missouri; eastern Kansas; southern Illinois, Indiana, and Ohio), and restricted areas of northeastern states (New Jersey and Pennsylvania). Suitable areas extend to include areas in western states (Arizona, California, Oregon, and Washington), although those areas are not likely accessible to the species (Fig. 2). Our models also showed suitable areas for the species across parts of eastern, western, and southern Mexico (Quintana Roo) (Fig. 2).

Figure 2: Left panel: potential suitable areas of Amblyomma maculatum based on binarized (5% threshold) models under current conditions (in blue and gray), and future (blue = no longer suitable, red = newly suitable) conditions. Right panel: agreement in strict extrapolation areas among the five general circulation models. Results are presented for RCP 4.5 (top) and RCP 8.5 (bottom). DOI: 10.7717/peerj.13279/fig-2

Future model transfers showed stable suitable areas (i.e., suitable in current time and in the future time) across the South, Midwest, and the Northeast, in the United States (Fig. 2). Areas of range reduction (loss) were in restricted areas in Kansas, Oklahoma, and Texas. Range expansion (gain) was anticipated in the northeastern (Pennsylvania, New York, Connecticut, Rhode Island, Massachusetts, Vermont, New Hampshire, Maine) and midwestern states (Kansas, Missouri, Nebraska, Iowa, Illinois, Indiana, Ohio, Michigan, Wisconsin) (Fig. 2). In general, we noted greater agreement among models in terms of losses and gains in the RCP 8.5 scenario compared to RCP 4.5.

Dermacentor andersoni

Current-time range predictions for D. andersoni showed broad suitable areas across Washington, Idaho, Oregon, California, Montana, Nevada, Utah, Wyoming, and Colorado, in cases where this species is known to occur in the United States. Climatically suitable areas extended across the Midwest, and Northeast and in some southeastern states (Fig. 3), although these areas are not likely accessible to the species. Currently suitable areas were also observed in parts of central and western Canada (British Columbia, Alberta, Saskatchewan, and restricted areas in Manitoba) (Fig. 3).

Figure 3: Left panel: potential suitable areas of Dermacentor andersoni based on binarized (5% threshold) models under current conditions (in blue and gray), and future (blue = no longer suitable, red = newly suitable) conditions. Right panel: agreement in strict extrapolation areas among the five general circulation models. Results are presented for RCP 4.5 (top) and RCP 8.5 (bottom). DOI: 10.7717/peerj.13279/fig-3

Future model transfers showed stable suitable areas across the states listed above, with some degree of reduction in suitable areas in the western states including much of Washington, Oregon, California, Nevada, Arizona, New Mexico, and Utah, and restricted areas in Colorado, Idaho, and Montana (Fig. 3). Predictions for the two RCP scenarios showed closely similar patterns of range stability, expansion, and loss, with more agreement among models in the RCP 8.5 scenario (Fig. 3).

Model uncertainty

MOP results for A. maculatum showed that strict extrapolative areas among future scenarios were concentrated in northern parts in North America, particularly in Canada, and in some restricted areas of the United States and southern Mexico (Fig. 2). Model variability results showed almost no variation coming from replicates and RCPs, but high contribution to variation from GCMs and parameter choice (Fig. S3).

In D. andersoni, we noted high agreement of strictly extrapolative areas in both southern and northern North America, and in lesser degree in the eastern United States and Canada (Fig. 3). High model variability came mainly from parameter choice in the eastern United States; we noted low variation deriving from GCMs, RCPs, and replicates (Fig. S4).

Discussion

The geographic distributions of A. maculatum and D. andersoni are much wider today than they were in the recent past. For example, A. maculatum has expanded its geographic range from the southeastern United States to become well established in the Northeast in Connecticut (Molaei et al., 2021), and in the Midwest in southern Illinois (Jolley, 2020). Beside the movement of tick adults for long distances via their hosts to new areas, immature A. maculatum can also access new areas with the help of migratory birds; larvae and nymphs can move thousands of miles during bird migratory seasons from the southern United States north to southern Canada (Florin et al., 2014; Teel et al., 2010). Cuervo et al. (2021) showed similar suitable ranges using current time predictors, and demonstrated that levels of niche conservatism differed among different members of A. maculatum group (A. tigrinum and A. triste).

This study is the first to assess the geographic distributions of the spotted fever vectors A. maculatum and D. andersoni in North America under current and future climate conditions. We included uncertainty analyses (MOP analysis and model variability) in our future model projections to detect areas with strict extrapolation, and to assess variation coming from multiple sources, such as different GCMs and RCPs. We considered only abiotic climatic variables such as temperature and precipitation as predictors that may influence the geographic distributions of those tick species.

Our models predicted that suitable areas for A. maculatum will remain stable in most southern and Midwestern states, whereas few reductions in suitable areas were anticipated only in western parts of Texas, Oklahoma, and Kansas (Fig. 2). Most importantly, our models predicted newly suitable areas northward in the United States successfully, to cover areas that were recently discovered to hold new populations in Connecticut and Illinois (Jolley, 2020; Molaei et al., 2021) (Fig. 4). For D. andersoni, our models showed broader suitable areas beyond its known range (from Washington state to Colorado). Midwestern and eastern states; however, most of the anticipated reduction in ranges were in areas not known to hold this tick species. Most of the anticipated expansions in range were in northward in Canada (Fig. 3).

Figure 4: Detail of Fig. 2, showing the most recent confirmed established populations of Amblyomma maculatum in counties of Illinois and Connecticut in the United States (light blue boundaries) (Jolley, 2020; Molaei et al., 2021). Gray indicates suitable areas under current and future conditions. Red color indicates newly suitable areas with climate change. DOI: 10.7717/peerj.13279/fig-4

Our projections suggested higher potential of A. maculatum to invade new areas outside its native range mainly in the southeastern United States. For D. andersoni, suitable areas were mostly in northern North America (United States and Canada). We also noted more extensive strict extrapolative areas for D. andersoni than A. maculatum, especially in the eastern United States, which suggested caution about interpreting those areas as suitable for D. andersoni (Fig. 3).

Several significant limitations and caveats regarding predictions emerging from ecological niche modeling that should be considered. First, a species faces dispersal limitations and biotic interactions that may prevent it from occupying the full suitable area that corresponds to its fundamental ecological niche. Second, the variation in spatial precision associated with different occurrence data records, which can cause problems for model results. Third, data availability in which biases in sampling in regions more than others can cause biases in model output (Peterson, 2014). All these points have been considered in the design of our methodology to achieve the most robust model possible.

In the United States, numbers of documented spotted fever cases have increased in recent years, especially in 2017, with 6248 new cases (CDC, 2021). Previous analyses have noted overlap between reported cases in some states and suitable areas for spotted fever vectors including Dermacentor variabilis (Alkishe, Raghavan & Peterson, 2021; Boorgula et al., 2020). Spotted fever case data collected by the Centers for Disease Control and Prevention, and used by Alkishe, Raghavan & Peterson (2021) however, were lacking in full detail on the type of pathogen and associated tick species, which made it difficult to interpret the source of the infection.

In summary, using ecological niche modeling allowed us to highlight suitable areas of two medically important tick species in North America. We showed the potential for expansion of those tick vectors into new areas that were not suitable in the past with emphasis on the newly discovered dispersal of A. maculatum to those newly suitable areas in Illinois and Connecticut. We also showed the uncertainty and variability that can come from projection models to different times and places. These results help to recognize the uncertainty and source of variability in predicting suitability.

Additional Information and Declarations

Competing Interests

The authors declare there are no competing interests.

Author Contributions

Abdelghafar Alkishe conceived and designed the experiments, performed the experiments, analyzed the data, prepared figures and/or tables, authored or reviewed drafts of the paper, and approved the final draft.

A. Townsend Peterson conceived and designed the experiments, analyzed the data, authored or reviewed drafts of the paper, and approved the final draft.

Data Availability

The following information was supplied regarding data availability:

The raw data is available in the Supplemental File.

Funding

This research was supported by a grant from the National Science Foundation (OIA-1920946). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

References

Aiello-Lammens ME, Boria RA, Radosavljevic A, Vilela B, Anderson RP. 2015. spThin: an R package for spatial thinning of species occurrence records for use in ecological niche models. Ecography 38:541-545

Alkishe A, Raghavan RK, Peterson AT. 2021. Likely geographic distributional shifts among medically important tick species and tick-associated diseases under climate change in north America: a review. Insects 12:225

Anderson RP, Lew D, Peterson AT. 2003. Evaluating predictive models of species’ distributions: criteria for selecting optimal models. Ecological Modelling 162:211-232

Animal Diversity Web (ADW). 2021. Dermacentor andersoni Rocky Mountain wood tick. Animal Diversity Web.

Barve N, Barve V, Jiménez-Valverde A, Lira-Noriega A, Maher SP, Peterson AT, Soberón J, Villalobos F. 2011. The crucial role of the accessible area in ecological niche modeling and species distribution modeling. Ecological Modelling 222:1810-1819

Boorgula GD, Peterson AT, Foley DH, Ganta RR, Raghavan RK. 2020. Assessing the current and future potential geographic distribution of the American dog tick, Dermacentor variabilis (Say) (Acari: Ixodidae) in North America. PLOS ONE 15:e0237191

CDC. 2018. Rocky Mountain spotted fever (RMSF). Centers for Disease Control and Prevention.

CDC. 2019. Other spotted fever group rickettsioses. Centers for Disease Control and Prevention.

CDC. 2021. Tickborne disease surveillance data summary. Centers for Disease Control and Prevention.

Cobos ME, Peterson AT, Barve N, Osorio-Olvera L. 2019. KUENM: an R package for detailed development of ecological niche models using Maxent. PeerJ 7:e6281

Cuervo PF, Flores FS, Venzal JM, Nava S. 2021. Niche divergence among closely related taxa provides insight on evolutionary patterns of ticks. Journal of Biogeography 48(11):2865-2876

Cumbie AN, Espada CD, Nadolny RM, Rose RK, Dueser RD, Hynes WL, Gaff HD. 2020. Survey of Rickettsia parkeri and Amblyomma maculatum associated with small mammals in southeastern Virginia. Ticks and Tick-Borne Diseases 11:101550

Dantas-Torres F, Chomel BB, Otranto D. 2012. Ticks and tick-borne diseases: A One Health perspective. Trends in Parasitology 28:437-446

Djalante R. 2019. Key assessments from the IPCC special report on global warming of 1.5 °C and the implications for the Sendai Framework for disaster risk reduction. Progress in Disaster Science 1:100001

Escobar LE. 2020. Ecological niche modeling: an introduction for veterinarians and epidemiologists. Frontiers in Veterinary Science 7:713

Florin DA, Brinkerhoff RJ, Gaff H, Jiang J, Robbins RG, Eickmeyer W, Butler J, Nielsen D, Wright C, White A. 2014. Additional US collections of the Gulf Coast tick, Amblyomma maculatum (Acari: Ixodidae), from the state of Delaware, the first reported field collections of adult specimens from the state of Maryland, and data regarding this tick from surveillance of migratory songbirds in Maryland. Systematic and Applied Acarology 19:257-262

Hijmans RJ, Cameron SE, Parra JL, Jones PG, Jarvis A. 2005. Very high resolution interpolated climate surfaces for global land areas. International Journal of Climatology 25:1965-1978

Jolley T. 2020. Disease-carrying coastal tick established in Illinois.

Molaei G, Little EA, Khalil N, Ayres BN, Nicholson WL, Paddock CD. 2021. Established population of the Gulf Coast tick, Amblyomma maculatum (Acari: Ixodidae), infected with Rickettsia parkeri (Rickettsiales: Rickettsiaceae), in Connecticut. Journal of Medical Entomology 58:1459-1462

Nadolny RM, Gaff HD. 2018. Natural history of Amblyomma maculatum in Virginia. Ticks and Tick-Borne Diseases 9:188-195

Owens HL, Campbell LP, Dornak LL, Saupe EE, Barve N, Soberón J, Ingenloff K, Lira-Noriega A, Hensz CM, Myers CE. 2013. Constraints on interpretation of ecological niche models by limited environmental ranges on calibration areas. Ecological Modelling 263:10-18

Peterson AT. 2014. Mapping disease transmission risk: enriching models using biogeography and ecology. Baltimore: Johns Hopkins University Press.

Peterson AT, Cobos ME, Jiménez-García D. 2018. Major challenges for correlational ecological niche model projections to future climate conditions. Annals of the New York Academy of Sciences 1429:66-77

Peterson AT, Papeş M, Soberón J. 2008. Rethinking receiver operating characteristic analysis applications in ecological niche modeling. Ecological Modelling 213:63-72

Phillips SJ, Anderson RP, Dudík M, Schapire RE, Blair ME. 2017. Opening the black box: an open-source release of Maxent. Ecography 40:887-893

R Core Team. 2018. A language and environment for statistical computing. Vienna: R Foundation for Statistical Computing. software

Rocklöv J, Dubrow R. 2020. Climate change: an enduring challenge for vector-borne disease prevention and control. Nature Immunology 21:479-483

Sonenshine DE. 2018. Range expansion of tick disease vectors in North America: implications for spread of tick-borne disease. International Journal of Environmental Research and Public Health 15:478

Teel P, Ketchum H, Mock D, Wright R, Strey O. 2010. The Gulf Coast tick: a review of the life history, ecology, distribution, and emergence as an arthropod of medical and veterinary importance. Journal of Medical Entomology 47:707-722

Warren DL, Seifert SN. 2011. Ecological niche modeling in Maxent: the importance of model complexity and the performance of model selection criteria. Ecological Applications 21:335-342

Wilkinson P. 1967. The distribution of Dermacentor ticks in Canada in relation to bioclimatic zones. Canadian Journal of Zoology 45:517-537

Wuebbles DJ, Fahey DW, Hibbard KA, Arnold JR, De Angelo B, Doherty S, Easterling DR, Edmonds J, Edmonds T, Hall T. 2017. Climate science special report: fourth national climate assessment (NCA4), Volume I. Washington, D.C.: U.S. Global Change Research 470 pp