Introduction

During the last decades, biological invasion studies have strongly benefited from the use of neutral molecular markers to disentangle routes of invasion1,2,3. This has helped to identify source populations, frequency of invasion events, and the geographical patterns and demographic consequences of invasive species spread4,5. Reconstruction of past events of invasion can help in the recognition of the mechanisms of dispersal to prevent the prevalence of invasion or further spread to additional areas1,6,7,8. However, understanding local dispersal remains a central challenge to prevent and control the economic and biodiversity costs of biological invasions.

Initially, molecular markers were used to identify sources and frequency of dispersal events of invasive species to non-native regions1,2,3,9. Whereas in some cases, these patterns explain the lower genetic variation of invasive species relative to native areas1,10,11,12, combination of migration events from different source populations have sometimes increased genetic variation within invaded regions, potentially increasing the risk of fast adaptation to novel conditions13,14,15,16. Once invasive species are established within non-native areas, understanding further local dispersal is a major challenge to identify potential environmental barriers and recommend control management programs7,17. It is possible to use highly variable molecular markers for understanding local scale dispersal and the entangled nature of the species invasion phenomena.

The cactus moth, Cactoblastis cactorum, is an oligophagous herbivore during its larval stage that consumes the inner tissue of cladodes of plants within the genus Opuntia, negatively affecting the plant survival18. This moth, native to South America, was used in 1924 for biological control purpose against exotic Opuntia species with successful results in Australia. Later, was intentionally introduced to South Africa (1933), New Caledonia (1933) and Hawaii (1950). In 1957 moth larvae, from South Africa and Australia, were introduced in the island of Nevis, seeking to control the exponential growth of native Opuntia populations caused by deforestation and cattle ranging19,20. According to a recent study, this has been the single introduction event of C. cactorum to the Caribbean area21. Since its introduction in Nevis and latter in Saint Kitts in 1957 and until 1987-1989, the moth was detected in other areas as far as Florida22, increasing the risk of extinction of native species23. The entire Florida peninsula is now invaded and the moth has spread through the Gulf of Mexico coasts20. According to the Texas Invasive Species Institute the moth has been detected in Brasoria County (Texas, US) in 2018 (around 800 km from the Northern Mexican border). Following the Atlantic coast of North America, it has dispersed up to Charleston County in South Carolina24,25,26. The rate of dispersal of the moth according to females flight distances was 3-6 km and 16-24 km in 2.5 years in South Africa and Australia respectively27,28. Hence the rapid spread throughout the Caribbean and North America probably involved additional factors. Previous reports suggest that spread among Hawaiian Islands probably occurred through island-hopping29 and tropical storms30. In turn, human-mediated dispersal and climatic events like hurricanes have been proposed as sources of dispersal in the Caribbean21,30,31,32. The threat for North American deserts is that they have the highest diversity of Opuntia cacti species, which may suffer from this oligophagous moth18. Furthermore, domesticated Opuntia in Mexico have a high cultural and economic value as a food resource33. Thus, given the risk of further spread to continental areas and the ecological, social and economic associated costs, the present study examined possible routes of local dispersal in the Caribbean and Florida to better understand the factors that may favor additional introductions of this invasive species to the continent.

A previous study using a mitochondrial gene (COI) reported that commercial transportation of ornamental cacti from Dominican Republic and Puerto Rico was the most likely route of invasion of C. cactorum to Florida31. However, a more recent study with the same marker and a more extended population survey suggested that hurricanes may have also contributed to the dispersal of the moth within the region32. Whereas commercial transportation of ornamental cacti has been controlled ever since the detection of this invasive species18, further introductions within the continent may occur via climatic events. Accordingly, the goal of this study was to further advance in the comprehension of the patterns of dispersal of C. cactorum by using nuclear microsatellites and computational Bayesian approaches. This approach will allow us to add a new piece of evidence to determine whether dispersal through trade, climatic events (hurricanes) or both, better account for the current genetic structure of C. cactorum in North America.

Material and Methods

Sample collection

From 2011 to 2012, a total of 229 larvae of the cactus moth, C. cactorum, were collected from 10 invasive sites within the Caribbean region and Florida (USA) (N = 185) (Fig. 1), and two native ones in Argentina (N = 44) from where the moth was initially collected and transported for biological control of ruderal Opuntia species in Australia in 1924 (Table 1). Within each collecting site, only one larva was selected per individual host of Opuntia spp. to avoid over-representation of genotypes belonging to the same genetic family. All samples were stored in 95% ethanol until DNA extraction.

Figure 1
figure 1

Geographic location of the invading Caribbean populations of Cactoblastis cactorum. The location is on the conductance map for the Caribbean and Atlantic sea obtained with the CIRCUITSCAPE program based on hurricane incidence values, the darkest areas show those with the highest probability of flow. This conductance map was built using the isopletes shown at Lugo et al.48 and available in NOAA (https://www.nhc.noaa.gov/).

Full size image
Table 1 Description of studied invasive and native populations of Cactoblastis cactorum sampled in the Caribbean, Florida and in the native region. Average number of alleles (NA), observed heterozygosity (Ho), expected heterozygosity (He), allelic richness (RA), and Fixation Index (FIS) are provided. Standard deviations are indicated in parenthesis. *significant P < 0.05 FIS values. Statistics calculated with 10 microsatellites.
Full size table

DNA extraction and microsatellite analyses

We performed total DNA extraction with DNEasy blood and tissue kits (Qiagen, MD, USA) following the manufacturer protocol. Population genetic variation was determined using nuclear microsatellites developed by Genetic Marker Services (http://www.geneticmarkerservices.com). We removed loci that did not amplify among the 29 primers pairs tested. We chose 14 potential polymorphic loci (Table 2) and labeled with fluorophores (Applied Biosystems) for fragment analysis. Each multiplex PCR mixture (5 µL) contained 2.5 µL RadyTaq (Qiagen cat. 206143) and 1 µL DNA template (20 ng) 5pmol for each fluorescent labeled primers. PCR were performed through touchdown reaction starting with initial heat activation at 95 °C for 10 min followed by 6 denaturation cycles of 94 °C for 1 min, annealing for 1 min and 1 min of extension at 72 °C; annealing cycle temperature began at 60 °C or 57 °C and decreased 1 °C every cycle. The PCR reaction ended with two stages of 12 cycles each (57° and 56°, respectively) and final elongation of 72 °C for 5 min. The PCR product was diluted and run on an ABI 3730xl automated capillary sequencer. The allele size was manually scored using a Liz 600 size standard (Applied Biosystems) in GeneMarker (V2.2.0) (Soft Genetics LLC, State College, Pennsylvania, USA).

Table 2 Fourteen nuclear microsatellites amplified on Cactoblastis cactorum (subscripts correspond to the groups for multiplex PCR reactions).
Full size table

Basic population estimators of genetic variation

Genetic variation within sampled populations was characterized using the mean number of alleles per locus (NA), the allelic richness (AR), the mean expected- (He)34 and observed heterozygositiy (Ho), the Fixation Index (Fis) and estimation of FST35 between each pair of sampling site using FSTAT 2.9.3.36. Mean allelic richness (AR) was calculated using the rarefaction method of Leberg37. Exact test for deviation from Hardy-Weinberg equilibrium were performed with GENEPOP38. We used FreeNA to determine the frequency of null alleles using EM algorithm39. Pairwise population differentiation was tested using only those loci that were in H-W equilibrium in more than 50% of the populations.

Genetic clustering of populations was examined using STRUCTURE40 (v 2.3.3). We chose the admixture model with correlated allele frequencies. Two analyses were performed, one including all native and invaded population samples, and one including only populations from the invaded area. The first analysis examined whether genetic differentiation has occurred after the human-mediated journey of C. cactorum from its native Argentina to the Caribbean (1957). The second analysis was performed to explore genetic structuring within the invaded Caribbean only. Because groups of larvae were collected from discrete sites, sampling locations were used as prior information41. Each run consisted of a burn-in period of 105 Markov Chain Monte Carlo (MCMC) iterations, followed by 106 iterations. 20 replicated runs were carried out for each value of the potential number of clusters (K) set, between 1 and 9. STRUCTURE HARVESTER42 (available at http://taylor0.biology.ucla.edu/structureHarvester/) was used to collating output results from STRUCTURE and to determine the uppermost level of structure using the method of Evanno et al.43. Because STRUCTURE assumes the absence of null alleles and H-W equilibrium for all loci, analyses were performed using 10 out of 14 microsatellites that fulfilled the required assumptions.

Isolation genetic and connectivity in the invaded region

The hypothesis that genetic differentiation between a pair of populations is a linear function of the geographic distance between them was examined through a Mantel test that correlates genetic and geographic distance matrices. Paired genetic distances were estimated as [FST/ (1-FST)] and geographic distances were estimated as the logarithm of the Euclidean distances (blog) using GenAlEx44 (v.6.4). Isolation was also tested considering that hurricanes may also shaped the genetic structure of populations in this area. To explore the relationship between genetic distances and the frequency of hurricanes incidence across pairs of populations, a resistance matrix was constructed using the CIRCUITSCAPE program45 (Fig. 1). Based on a previous subdivision of the geographic area using a grid with the probability of incidence of hurricanes (categories 1 to 5) following NOAA46,47,48. Pairs of populations connected by cells with high incidence of hurricanes were assigned a high conductance value (i.e., low resistance to migration). This rationale was applied to all pairs of populations to construct a resistance matrix (indicated by light areas in Fig. 1). To correlate the resistance matrix against the genetic distance matrix controlling for the geographic distance among populations, a partial Mantel test with 10,000 permutations was implemented.

Approximate Bayesian Computational analysis

To provide a quantitative evaluation of the dispersal routes of the cactus moth into Florida (USA), an Approximate Bayesian Computational Analysis was performed with DIYABC49. Analyses were performed with ten microsatellites that were at Hardy-Weinberg equilibrium. This approach uses genetic information provided by microsatellites data and historical data given by dates of first observation within invaded regions. Three hypothetical invasion scenarios were examined (Fig. 2), and parameter values were drawn from prior distributions (Table 3). In all three cases we used the native Argentinian population of Ayuí as the initial source of invasion. This population was genetically close in average with the invasive populations (Paired FST = 0.225). Since moths were taken from Argentina to Australia, and then to South Africa, the analysis considered a “ghost” population representing a non-sampled population that may represent this step. Given the location and date of introduction to the Caribbean, the sampled population of Antigua was considered the source population that initiated the invasion in the region. Antigua was preferred over Nevis because of the larger number of individual larvae in our collection, the high genetic similarity between both islands (Paired FST = 0.012, N.S.) and the short difference in time of introduction between them (1 year). The following step in the invasion process was characterized by the population La Romana (Dominican Republic), given that the moth was detected in the island of Desecheo, located between Puerto Rico and Dominican Republic in 1963. After this step, hypotheses differ in how the cactus moth reached Florida. In the first hypothesis, Cuba, represented by the sampled population of Santiago, is the source of the Floridian invasive population (Fig. 2, Scenario 1). Due to US embargo to Cuba, any dispersal of the moth from Cuba to Florida is more likely related to environmental agents (such as hurricanes) rather than commercial transportation. The second hypothesis considers that moths belonging to La Romana (Dominican Republic) directly invaded Florida (Highlands) (Fig. 2, Scenario 2). This route constitutes the recorded trajectory of commercial traffic to Florida of ornamental cacti infested with C. cactorum from Dominican Republic and Puerto Rico detected in 198950. However, we cannot rule out the possibility that hurricanes also contributed to dispersal among these areas. The third hypothesis considers that the population of Highlands (Florida) was the result of the admixture of genotypes belonging to Cuba and Dominican Republic (Fig. 2, Scenario3).

Figure 2
figure 2

Scheme of the three competing scenarios evaluated with the ABC model. The prior distribution values of parameters are described in Table 3 (0 = year 2012, assuming two generations per year back in time). The scenarios was built and evaluated with DIYABC (V. 2). The edition was made in the Inskape program.

Full size image
Table 3 Parameters used for data simulation in the three competing scenarios using ABC analyses.
Full size table

Given that C. cactorum usually has two generations per year27,28,51, this generation time was used to scale coalescent time in all scenarios. Dates of introduction to Australia (Tg) and Antigua (T2) were given as fixed values based on historical data of the intentional introduction events18. In the other areas (T1 (Argentina), T3-5 (Dominican Republic, Cuba and USA); Table 3), a minimum and maximum value was assigned to each time parameter, with the dates of first observation of the insect as the lower boundary (Table 3). We simulated three million microsatellite data sets (one million for each scenario). Subsequently, the probability for each scenario was inferred by polychotomous logistic regression on the 1% of the simulated data sets closest to the observed data set49. The selected scenario was chosen as that with the highest posterior probability value. To evaluate the robustness of the analysis, we used pseudo-observed simulated data sets to quantify the type I error rate (risk to exclude the focal scenario when it is the true one) and the type II error rate (risk to select the focal scenario when it is false)49. All simulations and ABC analyses were carried out in DIYABC (v.2) sofware52.

Results

No evidence of linkage disequilibrium was detected for all 14 microsatellite loci. Four loci were excluded from the analysis because they had more than 20% of null alleles and were not at Hardy-Weinberg equilibrium (i.e., cc7b, cc16, cc65, cc6b). With the remaining ten microsatellites loci we found that the native population from Yuquerí (Argentina) had on average the highest number of alleles (3.8), while the lowest average number of alleles (1.6) was observed in one of the most recent invasive population (Lee Co., USA) (Table 1). Allelic richness (AR) was higher in native range (2.016–2.185) than invaded populations (1.453–1.992; Table 1). A deficit of heterozygotes was detected in one of the native populations (Ayuí), and in three of the Caribbean invaded populations (Antigua, Nevis and Santiago; Table 1). On average, lower levels of heterozygosity were recorded for invaded than native populations (Table 1).

Pairwise comparisons indicated significant differentiation between the native Argentinian populations and all invaded populations in the Caribbean and Florida (mean FST = 0.229, range = 0.128–0.49). Within the invaded region, genetic differentiation ranged from almost no differences to rather high levels of genetic differentiation (Table 4). The average level of differentiation within the invaded region was FST = 0.185. In general, early invaded populations (Nevis and Antigua) presented lower levels of differentiation with the rest of the populations in the region than recently invaded populations. The Jamaican population (Palisadores) presented the greatest level of genetic differentiation within the Caribbean (mean FST = 0.310; Table 4), whereas Florida populations had a low level of differentiation with population La Romana (Dominican Republic) and low to moderate differentiation with two of the closest Cuban populations (Pinar del Río and Trinidad; Table 4).

Table 4 Pairwise FST of Cactoblastis cactorum populations. Values in bold were not statistically significant. The adjust for multiple comparisons was α = 0.0006.
Full size table

Clustering genetic analyses using STRUCTURE on the whole dataset revealed the presence of two main groups (K = 2, with the delta K method) distinguishing native from all invasive populations (Fig. 3B), and tend to confirm a single introduction event. When analyzing only the group of invasive populations, K = 2 is also selected (Fig. 3A). While individuals from Puerto Rico, Palisadores and Pinar del Rio belong to a well-defined cluster, the rest of the populations showed an admixed pattern. Analyses of isolation by distance and connectivity including only the invading populations indicated a low but almost significant level of isolation by distance (r = 0.165; P = 0.05). On the other hand, the isolation by resistance hypothesis based on hurricane and wind currents was significantly supported (r = 0.435; P = 0.008).

Figure 3
figure 3

Ancestry estimation based on the Bayesian clustering method STRUCTURE in the Cactoblastis cactorum samples. (A) Genetic clustering of the 10 introduced populations, assuming two population clusters (K = 2). (B) Genetic clustering of the 10 introduced populations and the 2 native populations, assuming two population clusters (K = 2). Note: each vertical line corresponds to an individual and the shades of gray represent the probability of belonging to a group. Individuals are grouped by population sample.

Full size image

Results from Approximate Bayesian Analyses (ABC) support quiet well the hypothesis described in the third scenario (Fig. 2), due to (1) a higher consistency between simulated and sampled data (i.e. the highest posterior probability; P = 0.7798), (2) non-overlapping confidence intervals, and (3) low Type I and mean Type II errors (0.038 and 0.000; respectively). This dispersal scenario considers that the invasion to Florida (represented by the Highlands population) was funded by insects belonging to at least two sources of dispersal, one from Dominican Republic and the other from Cuba (represented in the analysis by the Santiago population; Fig. 2). The posterior probability estimation of parameters of model three assigned a median admixture rate of R = 0.4.

Discussion

During the 1981–1991 decade, North American customs authorities intercepted Opuntia plants infested with C. cactorum. Seizures occurred in plant shipments sent from the Dominican Republic to Miami and in luggage to Dallas Texas International Airport53. Commercial trade of ornamental cacti from the Dominican Republic to the United States of America coupled with the constant influx of American tourists to the Caribbean islands lead to the proposal that human-mediated dispersal has been the major agent of migration to the continent. Similarly, the comparison of the mitochondrial haplotypes from eastern populations of C. cactorum of Florida and those of Dominican Republic supports this idea21,31,32. Despite other potential mechanisms of dispersal as island-hopping and tropical storms or hurricanes30, the present study support the genetic similarity and likelihood of migration between La Romana (Dominican Republic) and Highlands (Florida) populations. However, the distribution of genetic variation and ABC analyses of invasion routes within the Caribbean region also suggest the presence of other nonexclusive routes of dispersal.

In the Caribbean and North Atlantic region, hurricanes are a ubiquitous temporary phenomenon with consistent wind current patterns. Studies of genetic variation have shown that hurricanes and marine currents influence the patterns of dispersal and constitute a source of variation in ecological and demographic processes54. Wind has been identified as one of the most important factors promoting the spread and long-distance dispersal of invasive insect species7,55 and to enhance distance and/or speed of dispersion as in the case of the monarch butterfly56, the cattle egret57 or the mealy bug, Oracella aculata58. Also, the grasshopper Eumetopina flavipes, vector of the sugarcane virus Ramu, entered Australia from New Zealand favored by wind currents59,60. Wind is also responsible for long-distance dispersal of the wasp Megastigmus schimitscheki in France7. Although previously suggested, the role of hurricanes on local dispersal of C. cactorum within the Caribbean was not examined until recently32. Our previous results using mitochondrial COI32 coupled with the findings of this study using nuclear microsatellites, are consistent and provide new pieces of evidence supporting the role of hurricanes on the spread of C. cactorum toward mainland North America. In comparison with the isolation by distance model, the isolation by resistance migration approach was a much better predictor of the current pattern of genetic variation of C. cactorum in the Caribbean and Florida. This finding suggests that hurricanes were one of the main drivers of dispersal from the Caribbean to mainland North America as previously shown using a more conserved mitochondrial marker32. In addition, ABC analyses of the invasion routes supported our hypothesis that moths entered Florida both through commercial transportation from the Dominican Republic and from Cuba as a consequence of climatic events (hurricanes); whereas trade with Cuba was suspended as a result of the embargo. Moreover, this finding suggests that hurricanes may have also dispersed C. cactorum from its initial introduction in Nevis, Antigua and Montserrat to other islands onto the Caribbean region. Historical records of hurricanes in the region (https://www.nhc.noaa.gov/) further support this statement. After the initial human-mediated introduction of C. cactorum, three hurricanes connected the Lesser Antilles and the Dominican Republic between 1963 and 1967, four impacted the south of Cuba and the Dominican Republic between 1975 and 1980, and two passed from Cuba and impacted the USA in 1985 and 1987. Nevertheless, the presence of the moth in the lower Bahamas in 198361, and the possibility of island hopping, as recorded in the Hawaiian archipelago, can also represent a possible route of dispersal towards Florida. During the last decades at least three hurricanes impacted the Bahamas before reaching the Atlantic coast of Florida (Mitch in 1988; Andrew in 1992, and Katrina in 2005). Hence a more extended survey of populations in this area of the Caribbean will provide further insight on the dispersion routes across the region.

The process of invasion is generally associated with the loss of genetic variation due to successive demographic bottlenecks as the invasion of new locations proceeds62. This process can be counteracted if repeated invasion events occur in a given location, thus increasing the amount of genetic variation to equivalent, or larger levels of variation, than that observed on the native region17,63. Our analyses revealed that populations of C. cactorum from the Caribbean and Florida, sustained less genetic variation than in the native region, likely as a consequence of genetic drift. Given that present populations of C. cactorum outside the Caribbean region are located more that 8,000 km away in other biogeographical areas and continents, the genetic structure results support previous analyses using mitochondrial markers suggesting that the insect entered the Caribbean probably once in 195721. Similar to the Lepidopteran tomato pest, Tuta absoluta, in Africa and the Mediterranean Basin62, the occurrence of just one introduction of the cactus moth was enough to allow its spread throughout the Caribbean despite the lower genetic variability than the source population from South America. Unlike the tomato pest and other invasive insect species1,9 significant genetic differentiation was detected in the invaded region suggesting again a possible effect of drift. A more complete survey of invaded populations in Australia, South Africa and in other continents will help to decipher whether the magnitude of reduction in genetic variation will affect the adaptive potential of the moth into other continental areas with different climatic conditions as those of the Caribbean and Florida.

The combined effect of hurricanes and genetic drift reducing genetic variation through space will likely affect the evolutionary potential of the moth during future introductions to the continent. Hence, monitoring programs of climatic events like hurricanes in the region and environmental genomics surveys will add new insight to better understand possible barriers to the expansion of the moth into continental areas.