What Limits the Distribution of Liriomyza huidobrensis and Its Congener Liriomyza sativae in Their Native Niche: When Temperature and Competition Affect Species' Distribution Range in Guatemala.

Factors limiting distribution range for most species are generally unknown regardless of whether they are native or invasive. We studied factors that could enable or restrict the distribution of two cosmopolitan invasive leafminer fly species, Liriomyza huidobrensis (Blanchard) and Liriomyza sativae (Blanchard) in their native niche. In order to test which ecological and environmental factors affect leafminer distribution we conducted thermal tolerance assays, sampled along elevation gradients and modeled species distribution. Findings from the field and rearing chambers showed a physiological restriction due to high temperatures for L. huidobrensis at 28-29 °C, above which adult emergence is compromised. We also found that maximum temperatures below 22 °C, typical of tropical highlands, favored L. huidobrensis. L. sativae was found across a wider temperature range (i.e., from 21 to 36 °C) in Guatemala. Our finding of a physiological threshold in temperature for L. huidobrensis may enable us to predict its invasive risk when combined with the environmental conditions at horticultural ports of entry and the global agricultural landscape. Further, it strengthens our predictions on shifts in distribution of the leafminer fly under future climate. We also found a temperature mediated competitive exclusion interaction between the two herbivore species, where L. sativae occurred at temperatures < 22 °C only in the absence of L. huidobrensis. We show that parasitoids had a negative effect on the leafminer flies, which varied with host plant. Finally, we show the importance of taking a multiaspect approach when investigating what limits distribution and invasiveness of a species.

The threat of invasive species has increased with the globalization of trade. Invasions have consequences for biodiversity, and in the case of herbivores they have consequences for plant productivity and local economies (Vila et al. 2010 but see, Webber et al. 2012). As an example of the severity, 40% of major North American insect pests are exotic (Kim and McPheron 1993).

As species retain physiological and ecological requirements (i.e., niche conservatism), we may be able to predict new distributions of the species and estimate the risk of establishment of an invasive species (Wiens and Graham 2005). It is not uncommon to predict invasion vulnerable areas based on the physiological and ecological requirements from the invasive species native niche (e.g., Peterson and Viegalis 2001, Bradley et al. 2010, Zhu et al. 2012). In fact, using the thermal tolerance of eggs, Kearney et al. (2008) were able to predict the risk of expansion of the cane toad to nonanalog climates in Australia. In addition, studies of top-down effects and the conditions under which native natural enemies are effective controlling herbivore pest populations have rendered successful biological control of pests such as the introduction of the beetle Rodolia cardinalis (Mulsant) (Coleoptera: Coccinellidae) for the control of the cottony cushion scale Icerya purchasi (Maskell) (Hemiptera: Margarodidae) (reviewed by Caltagirone and Doutt 1989). Successful management of biological invasions benefits from baseline studies of the invasive species in their native niche. Physiological tolerance to temperature is useful in making global predictive models based on physiological constraints of invasive species and for an understanding of the environmental conditions under which biocontrol can be successful (Baker et al. 2000, Peterson 2003, Jímenez-Valverde et al. 2011).

In the case of invasive species, shifts in tolerances from the native niche may occur through two mechanisms: 1) Species acquires increased physiological tolerance to new environmental conditions (i.e., changes in its fundamental niche, primarily climate; Holt et al. 2005), for example laboratory experiments with the fruit fly Drosophila have found more plasticity in the evolution of tolerance to lower temperatures and desiccation than to higher temperatures (Hoffmann 2010); and 2) Species are restricted by ecological interactions at the native niche which no longer apply in the introduced niche; for example, strong competition from neighboring plants reduces the growth of Centaurea stoebe (Asteraceae) populations in its native range, but competition by neighboring plants has an insignificant effect on invasive populations of C. stoebe in Montana (Callaway et al. 2011). This mechanism has been found to be a norm rather than an exception on introduced niches (Early and Sax 2014). Ecological and biogeographic studies of invasive species in their native niche, in combination with the introduction history, enable us to understand the constraints to the distribution and population growth that the species has evolved under. This information has the potential to improve our prediction of the areas that are at risk of invasion (Strubbe et al. 2014). Further, understanding the thermal physiological limits of an invasive species in combination with the trophic interactions allows us to make more accurate predictions of its future distribution under climate change. With detailed studies on physiology, distribution records, and thermal tolerance assays, herbivore pests present a great opportunity to study factors limiting populations locally, regionally, and globally. In this study, we aim to understand the distribution of an invasive leafminer fly, Liriomyza huidobrensis (Blanchard) (Diptera: Agromyzidae), as shaped by climate and temperature tolerance, as well as by the presence of a closely related species, Liriomyza sativae (Blanchard), in Guatemala, Central America.

Study System

There are currently 456 Liriomyza species described (ITIS 2016) with over 20 of these recorded as pests worldwide (Spencer 1973). However, three Liriomyza species cause the vast majority of all damage caused to horticultural field crops: L. sativae, L. huidobrensis, and L. trifollii (Burgess; Spencer 1973). All species are polyphagous and native to the Americas, and have spread globally in recent decades (i.e., reports from Europe and Asia date from 1989 to 1995), acquiring a global distribution status over the past 20 years (Spencer 1973, Weintraub and Horowitz 1996, Shepard and Braun 1998, Kang et al. 2009). The current distribution of these species results from a combination of the species’ polyphagy, the change in land use toward more agricultural areas and the globalization of vegetable trade (Scheffer and Lewis 2005).

Host plant preference has been documented for both L. huidobrensis and L. sativae even though they are known as polyphagous species (Wei et al. 2000, Zhao and Kang 2003). Furthermore, different adult size and pupal weight have been related to the nutritional value of host plants (Musundire et al. 2012). We focus our study on two of these polyphagous species because we found that in Guatemala, L. sativae was the only species sharing host plants and habitats with L. huidobrensis (i.e., both species emerged from beans, Phaseolus vulgaris, broccoli, Brassica oleracea, lettuce, Lactuca sativa, and tomato, Solanum lycopersicum, G.R.C. unpublished data). This same pattern of host plant sharing by these two leafminers was also observed in Indonesia (Rauf et al. 2000). Moreover, in Guatemala, L. huidobrensis and L. sativae are the only agromyzid leafminers commonly found in horticultural crops throughout the highlands, although Liriomyza trifolii has only been found in a few lowland greenhouse collections (MacVean and Pérez 1996,1997).

L. huidobrensis has a worldwide distribution but appears to be primarily restricted to cool climates (Lanzoni et al. 2002, Scheffer and Lewis 2001, Weintraub and Horowitz 1996) and has not (yet) established in the USA (Weintraub et al. 2017). This may be attributed to regulatory detentions of vegetables at the port of entry to Miami, FL, where most fresh produce enters the USA from Central and South America and where ∼85% of all Liriomyza spp. interceptions are recorded (Borchert 2006) (the other ports of entry are San Juan, Puerto Rico with 11% of interceptions and Fort Lauderdale, FL with 4%). Since 1980s there have been quarantines of fresh produce due to species of leafminers in goods or the residues of insecticides used to control them (MacVean and Pérez 1996, 1997; Pérez et al. 1997).

The cooler highlands of South and Central America have been described as the possible origins of L. huidobrensis as well as the source of its spread to other parts of the world, although the native range of L. sativae is recorded for both North and South America (Spencer 1973, 1981; Parrella 1987; Weintraub et al. 2017). Phylogeographic studies reveal that Guatemala may be the origin of incursion for Asian populations of L. huidobrensis and L. sativae as DNA from the Guatemalan specimens was found to be the closest to the invasive populations in Asia, Israel and Africa (Scheffer 2000, Scheffer and Lewis 2005). At locations where L. sativae is introduced it has not been observed at high elevations or temperate latitudes whenever L. huidobrensis is present. Furthermore, studies of the thermal tolerant genes of these two species have revealed that only L. huidobrensis thrives at cooler temperatures (Huang et al. 2007). Spencer (1981) reported that L. sativae replaces L. huidobrensis at lower altitudes in South America, which matches collection records in Guatemala (MacVean and Pérez 1996, 1997) and provided additional justification for focusing on these two species.

We ask first which factors predict leafminer species’ occurrence based on environmental variables and the co-occurrence of a close species. We then look within the region of leafminer occurrence under a trophic cascade frame and ask which factors affect their abundance in crops. Finally, we use physiological constraints to predict the risk of invasion and establishment of L. huidobrensis in southern USA.

Materials and Methods

We first took a mechanistic approach and tested for physiological constraints of L. huidobrensis to changes in daily temperature regimes. Then we studied which factors affect presence/absence of the two leafminer species through environmental niche models (ENM) using data from extensive field sampling across temperature gradients. We use the term species distribution in the study to encompass both the environmental niche of the species (determined by abiotic factors, namely climate variables, and ENMs) and the biotic factors that affect species distribution (i.e., preferred host plants, pressure from parasitoids, and competitive interference). Finally, we used the trophic cascade framework and analyzed the emergence of each fly species and ask which factors affect the abundance of both leafminer flies, on crops.

Factors Limiting the Geographical Distribution (i.e., Presence/Absence) of the Two Leafminer Species

Physiological Constraints of L. huidobrensis to Daily Maximum Temperatures

We reared caged populations of L. huidobrensis in controlled rearing chambers at the Universidad del Valle’s campus in Guatemala from 1999 to 2000. The parental generation for all experiments consisted of 475 pupae originated from a single collection of infested snowpeas (Pisum sativum L.) in highland Guatemala (1,840 m). Mined leaves containing larvae were transferred to the laboratory and into sealed containers with a mesh top for ventilation. Leaves inside the containers were placed in layers of paper toweling to avoid drowning of larvae due to excess humidity. Ambient temperature conditions were similar to the field site of initial collection, with a daily high of ca. 24°C and daily low of 13°C. The leaves were checked every other day for pupal emergence. When pupae were found they were removed with a fine brush and placed on open Petri dish bottoms with filter paper and then placed in rearing cages. All rearing was carried out on potted snowpeas which were placed in rearing cages 25–30 days after planting, and which served as the substrate for oviposition and larval development. Response variables included percentage of adult emergence and overall reproductive rates based on the number of pupae produced in each generation (R, net per capita reproduction per generation; Price 1997, Wilson and Bossert 1971).

Caged populations (cohort) of flies were reared with a range of experimental regimes of temperatures, humidity, photoperiod that span Guatemalan climate conditions as well as those of potential introduction areas in southern Florida (Table 1). The Guatemala regimes correspond to the central highland province of Chimaltenango where permanent populations of L. huidobrensis are found, most of them thriving in snowpea plantations. All rearing was conducted in programmable rearing chambers with smooth transitions between night and day temperatures together with photoperiod and humidity settings. For each temperature regime (treatment) two or three cages of flies were utilized with 50–100 pupae in each as the initial cohorts, emerged and reared through the adult stage of the next generation. Treatments were replicated over time with the total number of cohorts (n) shown in Figure 1a and b. At each temperature regime, cohorts were reared for at least six generations or as many as possible until reduced adult emergence diminished the population to extinction.

Table 1.

Environmental conditions used for rearing L. huidobrensis in controlled growth chambers

Maximum temperature (°C, ±0.5°)	Minimum temperature (°C, ±0.5°)	Photoperiod (L:D) h	% relative humidity day	% relative humidity night	Regime mimics climatic conditions at: (locality, month)
32	25	14:10	80	90	Miami, FL. July
30	23	14:10	80	90	Miami, FL. May–June
29	21	14:10	70	90	Miami, FL. April–May
28	20	13:11	80	90	Miami, FL. April
28	13	13:11	65	95	Chimaltenango, Guatemala. April
27	13	13:11	60	90	Chimaltenango, Guatemala. April
26	13	13:11	65	95	Chimaltenango, Guatemala. May
25	16	11:13	70	90	Miami, FL. February
24	11	12:12	65	95	Chimaltenango, Guatemala. February

Correspondence is shown between a particular regime and climatic conditions in either the native breeding area of the Guatemalan highlands (province of Chimaltenango) or in Miami, Florida. Historical weather records for Miami were obtained from the Southeastern Regional Climate Center (SERCC 1999), and from the National Institute for seismology, volcanology, meteorology and hydrology (INSIVUMEH 1999) for Guatemala. For each temperature regime cages of 50–100 pupae were used as the initial cohorts. See Figure 1 for n values (number of cohorts) at each temperature.

Fig. 1.

Performance of Liriomyza huidobrensis when multiple generations were reared at different temperature regimes. (a) Mean ± standard error (SE) percent emergence per intervals of temperature 23–32 °C. A log 10 transformation was effective to normalize percent emergence data prior to analysis of variance (ANOVA). F4, 56 = 11.6 **. Tukey HSD test was used to separate means with P < 0.05. (b) Mean ± standard error (SE) population growth rate (R) measured at different temperature regimes ranging from 23 to 32 °C. N values (number of cohorts) are shown above temperatures. A square root transformation was effective to normalize R values prior to ANOVA. F7, 36 = 7.09**. Tukey HSD test was used to separate means with P < 0.05.

We found that when six or more generations with R at or above 1.0 (increasing reproductive rate) were obtained, the population could be reared in stable fashion for 15 or more generations. Conversely, when a temperature regime was detrimental, R values below 1.0 and population crashes occurred before six generations were completed. Preliminary experiments showed that subjecting flies to the highest temperature regime (daily maximum of 32°C) directly from control conditions yielded the same results as a gradual transition across the range of temperatures over many generations (MacVean 1999); therefore the data for percent emergence and reproductive rate (R) were pooled for a given high-temperature regime.

We did not conduct this experimentation for L. sativae as its widespread distribution in the USA suggested that it favored lowland hot and dry conditions (Spencer 1983). These conditions may overlap with the niche of L. huidobrensis which, in Guatemala, we have found on highland snowpeas (MacVean and Pérez 1997). Haghani et al. (2007) found that the development of L. sativae on cucumber, Cucumis sativus (a common host plant) increased with temperature over the range typically found in vegetable growing areas of Guatemala (15–35 °C).

Field Collection Methods

We focused on the elevational distribution of L. huidobrensis and L. sativae in vegetable plots across broad elevation gradients in the Guatemala highlands, in order to include a wide array of environmental temperature regimes. The correlation factor between elevation and maximum temperature of hottest month (r2 = 0.85, P > 0.0001) was sufficiently high to use elevation measures as a surrogate of temperature. We sampled and reared these flies during a period of approximately three years (2001–03). Our sampling points spanned the central highlands and volcanic chain of Guatemala (altiplano) from which a large percentage of forested area has been transformed to community-based, commercial horticultural plantations. The elevation gradients included a range from 200 m near the Pacific coastline up to 2,600 m in the highlands; all samples were from a range of plantation sizes, from 10 m2 up to 1,000 m2. No greenhouse collections were included in the analysis (for more information on sampling locations and elevation ranges, see Supp Table 1 [online only]).

We sampled leaves with larval fly mines from all the different species of vegetable host plants represented at each farm, recorded the GPS location, name of the site or farm, the host plant species, and collection date (Supp Table 1 [online only]). Samples from each individual plant were placed in separate bags marked with sample number, date and location and handled as described in Physiological constraints section in the laboratory at the Universidad del Valle in Guatemala City. When pupae were found, they were removed with a fine brush, counted and placed on labeled Petri dishes with filter paper. When adults emerged they were collected, counted and placed in 75% ethanol. All samples were then keyed to species; a subset of these samples was mounted or prepared as genetic material for species confirmation. Fly samples were sent to Sonia Scheffer at the USDA Systematic Entomology Lab in Maryland for molecular identification as described in Scheffer (2000) and Scheffer and Lewis (2001, 2005). Parasitoids which emerged from fly pupae, that is, endoparasitoids, were keyed and confirmed with the collaboration of R. Cave, then entomologist at the Escuela Agricola Panamericana in Zamorano, Honduras, expert in parasitoid taxonomy. Exoparasitoids, which develop externally on fly larvae, are not reported since percent parasitism cannot be readily calculated (>1 parasitoid may have come from a given fly pupa). Reference specimens were retained both at the U. del Valle and Zamorano entomological collections.

From the field collection data, we extracted the modal temperature at which the leafminer host plants were grown at in Guatemala. To obtain a range, we also extracted the maximum and minimum temperatures at which the vegetables were grown.

Environmental (Climatic) Niche Modeling of the Geographical Distribution of the Leafminer Species

Presence data for L. huidobrensis (n = 39) and L. sativae (n = 65) from field sampling were used as input for niche modeling to get a better understanding of the environmental limitation to the distribution of both leafminers. Locality data were georeferenced with an average precision of the nearest degree and transformed to a shapefile using ArcGIS9. In total nineteen variables related to temperature and precipitation (PPT) were obtained at a scale of ∼1 km2 from Worldclim (http://worldclim.org; accessed 21 June 2017; Hijmans et al. 2005). To avoid over fitting, we only included six environmental variables that had the lowest correlation indices. We also favored variables that would match the temperature regime conducted in laboratory studies to evaluate how well laboratory studies predict presence/absence in the field. We did not include elevation in the model as a predictor as many of these bioclimatic variables strongly correlated with elevation (Supp Table 2 [online only]). We created two models for both species using the same environmental variables: mean temperature of the warmest quarter (MTWQ), maximum temperature of the warmest month (MAXTEMP), temperature seasonality, PPT of the driest month, PPT of the driest quarter, and the mean PPT during the coldest quarter (for correlation matrix see Supp Table 2 [online only]). We first created an environmental-only model which modeled the climatic niche, for both species individually using the maximum entropy algorithm (MaxEnt, Phillips et al. 2004).

Model outcome can vary based on the geographical extent used to train the model (Anderson and Raza 2010). It is also based upon the order with which variables are entered in the models, an automated process carried out by MaxEnt (Elith et al. 2011). We modeled the distribution of both species in Guatemala, however we extended the geographical extent of the study region ∼500 km to the east and west of Guatemala to include the entirety of neighboring countries Belize, El Salvador, and Honduras and nearly all of Nicaragua and large parts of Mexico up to the state of Oaxaca. This was done to make sure that the breadth of the climatic conditions used to generate the predictions captured the climatic niche of the species in Guatemala. Furthermore, this extent was deemed neither too restrictive nor too large to negatively affect the model performance, as established by Barbet-Massin et al. (2012). The maximum number of background points, which MaxEnt generated from this extent and treats as absence data, was set at 10,000, which is recommended (Phillips and Dudik 2008, Barbet-Massin et al. 2012). We used the default convergence threshold (10−6) and maximum number of iterations (500) values. Hinge features were applied.

The area under the curve (AUC) of a receiver operating characteristic (ROC) plot was used to test the accuracy of the models. Furthermore, we manually selected 30% of the occurrence data to set aside as calibration data by means of randomized partition to compare the AUC of the calibration models with those of the training models that consisted of the remaining 70% of the occurrence data. The models were replicated 23 times, using cross validation. MaxEnt could not be forced to run more than these 23 models because of the number of occurrences included in the model. We however believed 23 models were enough to capture model variations. We used the average of these runs for further use and analyses to capture the model variability, based on recommendations by Marmion et al. (2009). We then created a competitive-exclusion model which resembles the realized niche of the species to a better extent than the environmental-only model. In order to incorporate how the occurrence of one species affects the distribution of the other, we added the average probability distribution range of the 23 environmental-only models of the competitor species (1), as a predictor variable in addition to the bioclimatic data to the competitive-exclusion model of species (2). Similarly, we ran a competitive-exclusion model of species (1) by adding the average probability distribution range of the 23 environmental-only models of species (2) as a predictor variable to the competitive-exclusion model, following Hof et al. (2012). These competitive-exclusion models were also run 23 times for each species. MaxEnt uses a continuous scale to predict the suitability of a region which we transformed into binary suitable/unsuitable area by applying a cut-off threshold that was based upon the smallest difference between sensitivity and specificity (Liu et al. 2005). From the resulting output, we calculated the differences in area between the environmental-only and the competitive-exclusion models for both species and calculated the allopatric and sympatric regions at which the flies occurred. We used R version 3.2.2 (R Development Core Team 2014) for the analyses of the results generated by MaxEnt (version 3.1; http://www.cs.princeton.edu/∼schapire/maxent; Phillips et al. 2004).

Environmental and Ecological Effects to the Abundance of the Leafminer Species, Under the Trophic Cascade Framework

We used results from our field sampling to test for factors limiting the abundance of the two leafminer species. Given that not all host plants are grown at all locations, an inherently incomplete data set emerges. However, despite some loss of statistical power, we could conduct the analyses as outlined below.

Host Plant

We asked if host plant and ambient temperature affected the abundance of the leafminer flies and conducted a multivariate analysis of variance (MANOVA) that included host plant and MAXTEMP. We followed with a profile analysis, to test for possible interactions between host plants, MAXTEMP and the emergence of the two species of leafminers.

Co-occurrence

We separated the data into areas and host plants where each leafminer species occurred in isolation from the areas and host plants where they co-occur. We then conducted an analysis of covariance (ANCOVA) to test the effect of MAXTEMP on the emergence at areas of isolation for each species versus areas of co-occurrence (analyses were conducted with SAS v.9.2) (SAS institute, Inc., Cary, NC).

Top-down Effect From Endoparasitoids

To calculate the pressure on leafminer abundance by parasitoids we exclusively looked at endoparasitoids which emerged from the pupal stage of the leafminers, so that we knew that each parasitoid emerged from a single leafminer individual. We calculated endoparasitism by dividing the number of endoparasitoids emerged by the total number of pupae. To separate the effects of parasitism for a single leafminer species, we excluded samples in which >1 species of leafminer emerged.

Path Analysis of Trophic Cascades, for Different Host Plant Species

We constructed possible interactions among host plant, herbivore and natural enemies and analyzed the potential effect of maximum temperatures on the interactions to generate new hypotheses and lead future studies of trophic interactions and food web dynamics with leafminers. With our natural history knowledge, we used path analysis to test for possible relationships among these variables (Shipley 2000). We only modeled the host plants for which we had >20 samples with emergence from the two species of leafminer flies. This included host plants that represented a broad range of temperatures: lettuce, onion (Allium cepa), tomato, and beans (Fig. 2b). We used maximum temperature of the hottest month, percent emergence from the two co-occurring leafminer species and parasitism rates. We then tested for possible causal relationships between environmental and ecological variables and the emergence of both leafminer species by using proc CALIS in SAS v.9.2 (SAS institute, Inc., Cary, NC).

Fig. 2.

Temperature intervals at which herbivores and host plants thrived. (a) Above maximum, minimum and modal temperatures of the locations we found the two leafminer species in Guatemala are shown. Below the temperature limits found for Liriomyza huidobrensis (this study) and L. sativae (Haghani et al. 2007). (b) Maximum, minimum and modal temperatures of the locations we found: Lettuce (Lactuca sativa), onions (Allium cepa), tomato (Solanum lycopersicum) and beans (Phaseolus vulgaris).

Results

Biogeographically, L. huidobrensis was restricted to plants grown in the altiplano (i.e., Guatemalan highlands; shown in green; Fig. 3a), where the maximum temperatures of the warmest month were on average 23.6 ± 0.23°C (Fig. 3b). This species was limited by the high temperatures typical of mountain valleys or the coasts at lower elevations (i.e., areas with orange-yellow; Fig. 3b). L. sativae had a wider distribution range across the Guatemalan vegetable plots (Fig. 3a; shown in gray). Its presence was limited by low temperatures, typical of the Guatemalan altiplano (i.e., areas in blue; Fig. 3b). During the three year sampling we found three samples in which three other species of leafminers emerged: Liriomyza sabaziae (Spencer), Liriomyza commelinae (Frost), and Liriomyza irazui (Spencer) feeding on weeds. We also found L. trifolii, but only on 13 locations from lowland crops including bell peppers, Capsicum annuum, and cucumbers. We also found that host plants were grown in different parts of the temperature gradients (Fig. 3b), making it ambiguous to discern if the differences in leafminer abundance per species is due to host plant preference or temperature.

Fig. 3.

Modeled environmental niche of the two species and the influence of maximum temperature of warmest month. (a) Environmental niche of Liriomyza huidobrensis (in green/black) and Liriomyza sativae (in gray) in Guatemala. (b) Map of the maximum temperature of warmest month gradient that goes from red (34.5–37.5 °C) to blue (10.6–18.7 °C).

Factors Limiting the Geographical Distribution of the Two Leafminer Species

Our controlled rearing experiments yielded results consistent with an upper developmental threshold at 28–29 °C (Fig. 1a). Above this limit, adult emergence dropped sharply and populations survived only a single generation or less, with lethal effects apparent during adult emergence. Although oviposition and larval development occurred, adult emergence rarely progressed beyond the initial opening of the puparium, and most adults died before emerging. No adults emerged successfully at temperatures of 32 °C. Within the range of threshold temperatures (28–29 °C) populations survived for variable numbers of generations, but percent emergence was low and growth rate (R) values were always < 1, insuring extinction (Fig. 1b). Below the threshold, with daily maximum temperatures of between 24 and 27 °C, populations were viable and maintained R values > 1. Our experiments provided no evidence for either genotypic selection or acclimation as mechanisms for increasing the upper temperature threshold in L. huidobrensis. In the field, we found no L. huidobrensis on sites with MAXTEMP above 27 °C (Fig. 2a). Both sets of experiments yield results consistent with a narrow threshold of L. huidobrensis at 27–29 °C (Fig. 2a) which prevents the species from occurring in Guatemala at elevations below 1,000 m.

For L. huidobrensis there was a strong influence of environmental factors setting a limit to its distribution across elevation in Guatemala (Fig. 3a). Mild maximum temperatures of the hottest month (MAXTEMP, Fig. 3b) were a strong predictor of the presence of this species.

L. sativae’s distribution model (shown in gray, Fig. 3a), was limited by PPT of the coldest quarter, MTWQ, and the presence of L. huidobrensis.

According to the environmental models (Supp Table 3 [online only]; AUC = 0.98) the best predictor of L. huidrobensis presence was MAXTEMP (71.45% of the AUC), followed by the PPT of the driest month (26.30% of the AUC). The ENM for L. sativae (Supp Table 3 [online only]; mean AUC = 0.96) generated by MaxEnt using environmental variables found that the PPT of the coldest quarter (46.98% of the AUC), the MTWQ (23.32% of the AUC), and the PPT during the driest quarter (20.96% of the AUC) were the strongest predictors.

When accounting for geography and elevation, we found that while 95.4% of the environmental niche (ENM) for L. huidobrensis is suitable for L. sativae, only 17.8% of L. sativae’s environmental niche is suitable for L. huidobrensis (Fig. 4). L. huidobrensis dominates the highlands while L. sativae dominates the lowlands of Guatemala (Fig. 1b); L. huidobrensis therefore had a narrower environmental niche than L. sativae in Guatemala. At mid to high elevations (1,000–2,300 m), where MAXTEMP ranges from 29 to 35°C, both species often co-occur, and have been found to share host plants (Fig. 4a and b). The results of the “competitive-exclusion model” for predicting the presence of L. huidobrensis (Supp Table 3 [online only]; AUC = 0.98) showed that MAXTEMP (71.3% of AUC) and the PPT of the driest month (25.6% of the AUC) were the most important variables in predicting presence. Emergence of L. sativae did not have an effect on whether L. huidobrensis was present or absent. Interestingly, the “competitive-exclusion model” for L. sativae showed that the PPT of the coldest quarter (Supp Table 3 [online only]; 57.51% of the AUC), MAXTEMP (19.78% of the AUC), and the % emergence of L. huidobrensis (16.30% of the AUC) were the strongest predictors. Hence, presence of L. sativae may be affected by L. huidobrensis in the cooler highlands, and by other leafminers, such as L. trifolii in the lowlands (G. Rodriguez-Castañeda, unpublished data).

Fig. 4.

Overlap between the environmental niche of the two leafminer species. (a) Map of Guatemala and the range in which the two species co-occur in orange. (b) Unique area for Liriomyza sativae in yellow, where only 17.8% of the environmental niche of Liriomyza sativae is suitable for Liriomyza huidobrensis. The unique area for Liriomyza huidobrensis is in dark red (i.e., only 4.6% of Liriomyza huidobrensis occurs without any Liriomyza sativae).

L. huidobrensis had its highest emergence on snow peas (P. sativum; mean = 58.7 ± 7% emergence, Supp Table 4 [online only]) and L. sativae had its highest emergence on tomato (S. lycopersicum; mean = 47 ± 3% emergence, Supp Table 4 [online only]). The difference in emergence of leafminer flies in the different host plants was significant (Fig. 6; MANOVA with profile analysis, Wilk’s lambda = 3.19). The results of the profile analysis on our MANOVA show that percent emergence of the two flies is opposite in most of the 10 host plants evaluated (i.e., the significant profile analysis was emergence of L. sativae – emergence of L. huidobrensis; F9, 112 = 3.19, P < 0.002; see Supp Table 5 [online only] for full statistics). We found that temperature mediates the effect of host plant on leafminers because there was a significant difference in both host plants and the mean maximum temperatures at which the host plants were found (Figs. 2b and 5; Wilk’s lambda = 0.429 see Supp Table 5 [online only] for full statistics); hence, it is not possible to discern if leafminer flies have a host preference or do best in those crops because of the temperatures the crops are grown at.

Fig. 6.

Contrasts of the two leafminer species that either occurred individually (only one species emerged from sampled area) or co-occurred (i.e., both species of flies emerged from the same crop species in the sampled area). (a) Mean emergence ± SE of the two Liriomyza species: Liriomyza huidobrensis () and Liriomyza sativae () occurred individually or they co-occurred in area and host plant. (b) Mean maximum temperature of warmest month ± SE at which the leafminer species Liriomyza huidobrensis () and Liriomyza sativae (): occurred individually or co-occurred in space and host plant. ** When P > 0.01.

Fig. 5.

Mean ± SE adult emergence of Liriomyza huidobrensis () and Liriomyza sativae () in the different vegetable host plants. Host plants were arranged on the x-axis from the ones that grow at mean coldest temperatures (in blue) to the ones that grow at mean warmest temperatures (in red).

The analysis of covariation showed that there is a difference in emergence of L. huidobrensis whether it occurs in isolation or co-occurred with L. sativae (Fig. 6a;F5, 202 = 8.38, P = 0.0003). It also showed a difference of percent emergence with higher temperatures (Fig. 6b;F1, 202 = 14.54; P < 0.0001). Moreover, there was an interaction between areas of isolation/co-occurrence and MAXTEMP, so that the lower emergence from L. huidobrensis (observed in areas where it co-occurs with L. sativae) is not independent from the effects of MAXTEMP on L. huidobrensis in those areas (F2, 202 = 7.45. P = 0.0008). The ANCOVA analysis for the emergence of L. sativae shows that its emergence is strongly affected by whether it occurs in isolation or it co-occurs with L. huidobrensis (Fig. 6a, F2, 202 = 58.14, P < 0.0001). For L. sativae there is no effect of MAXTEMP or interactions between isolation/co-occurrence area and MAXTEMP (for full results see Supp Table 6 [online only]).

Endoparasitoids

Mean endoparasitism rates were similar between L. huidobrensis (mean = 12.88 ± 1.87, n = 78) and L. sativae (mean = 12.86 ± 1.40, n = 109). Endoparasitism rates were significant in reducing emergence for both species (r2 = 0.16, P < 0.001 for L. huidobrensis and r2 = 0.08 P < 0.01 for L. sativaeSupp Table 7 [online only]). Mean parasitism rates varied per leafminer species and per host plant (see Supp Table 4 [online only]).

We identified seven parasitoid species that emerged from L. huidobrensis, however, only two species were widespread, and comprised 84% of all instances of endoparasitism: Opius dissitus (Muesebeck) (Hymenoptera: Braconidae) and Halticoptera circulus (Walker) (Hymenoptera: Pteromalidae; Supp Fig. S1a [online only]). We identified ten species of endoparasitoids that emerged from L. sativae, but only two of these followed the species throughout its distribution: O. dissitus and Ganaspidium utilis (Beardslay) (Braconidae: Eucoilidae) which comprised 75% and 13% of the parasitism cases, respectively (Supp Fig. S1b [online only]).

We found that the same structure and directions of trophic interactions held for all host plants. However, the relative importance of maximum temperature, competition and top-down forces from parasitoids differed for all host plants (Fig. 8).

For lettuce host plants, there was strong support for MAXTEMP and parasitism rates affecting the abundance of L. huidobrensis (path coefficient PC = −0.72 and −0.68, respectively). In addition, there was a direct positive effect of MAXTEMP on the abundance of L. sativae (PC = 0.43). Temperature mediated interspecific competition between the two leafminer species where the abundance of L. huidobrensis limited the presence of L. sativae (PC = −0.46) and MAXTEMP had an indirect positive effect on L. sativae (PC = 0.33). Parasitism had a mild negative effect on L. sativae (PC = −0.15). Finally, parasitism had an indirect positive effect on the abundance of L. sativae through limiting the abundance of L. huidobrensis (PC = 0.31).

When onion host plants were studied, MAXTEMP and parasitism rates affected the abundance of L. huidobrensis (PC = −0.36 and −0.54, respectively). There was also a strong direct positive effect of MAXTEMP on the abundance of L. sativae (PC = 0.7). Temperature mediated interspecific competition between the two leafminer species where the abundance of L. huidobrensis limited the presence of L. sativae (PC = −0.55) and MAXTEMP had an indirect positive effect on L. sativae (PC = 0.3). Parasitism had a negative effect on L. sativae (PC = −0.21). Finally, parasitism had an indirect positive effect on the abundance of L. sativae through limiting the abundance of L. huidobrensis (PC = 0.2).

In tomato host plants, MAXTEMP and parasitism rates affected the abundance of L. huidobrensis (PC = −0.33 and −0.22, respectively). There was a mild positive effect of MAXTEMP on the abundance of L. sativae (PC = 0.06). Temperature mediated interspecific competition between the two leafminer species occurred so that the abundance of L. huidobrensis limited the presence of L. sativae (PC = −0.45) and MAXTEMP had an indirect positive effect on L. sativae (PC = 0.14). Parasitism had a negative effect on L. sativae (PC = −0.11). Finally, parasitism had an indirect positive effect on the abundance of L. sativae through limiting the abundance of L. huidobrensis (PC = 0.1).

For bean host plants, MAXTEMP and parasitism rates affected the abundance of L. huidobrensis (PC= −0.63 and -0.33, respectively). There was a positive effect of MAXTEMP on the abundance of L. sativae (PC = 0.22). Temperature mediated interspecific competition between the two leafminer species where the abundance of L. huidobrensis limited the presence of L. sativae (PC = −0.33) and MAXTEMP had an indirect positive effect on L. sativae (PC = 0.2). Parasitism had a strong negative effect on L. sativae (PC = −0.35). Finally, parasitism had an indirect positive effect on the abundance of L. sativae through limiting the abundance of L. huidobrensis (PC = 0.1).

Discussion

In this study, we were able to find both ecological and environmental limits to the distribution of two cosmopolitan leafminer fly species in Guatemala, which is part of their native niche.

Laboratory experiments revealed that exposure to maximum daily temperatures above 28–29°C poses a sharp physiological constraint to L. huidobrensis pupal development (Fig. 1a and b). Similarly, laboratory assays of L. huidobrensis on bean plants (P. vulgaris) found that no pupae survived at temperatures > 30°C in Korean (Maharjan and Jung 2016) and Italian populations (Lanzoni et al. 2002). Upper temperature limits were further corroborated by field collections in which we did not find L. huidobrensis in lowland sites that have maximum temperatures >26.5°C, which is the ideal developmental temperature found in the laboratory (Maharjan and Jung 2016). ENM (Fig. 3a, Supp Table S3 [online only]) and generalized linear models predicting presence of the species (Table 1) all concur on the strong effect MAXTEMP has on whether L. huidobrensis is present or absent from a vegetable plot in Guatemala.

Our finding of a narrower temperature tolerance niche for L. huidobrensis in the field than with experimental physiological limits (Figs. 1a, b and 2a) can be explained by the fact that maximum temperatures of a field site are not necessarily the same as the exact microclimatic conditions that a leafminer experiences inside the plant, which are important to consider for the risk of invasion (Capinha et al. 2014). What is interesting is that similar thresholds have been found in the native niche of Peru with no flies surviving temperatures >32°C (Mujica et al. 2017), or in the introduced ranges of Italy and Korea (i.e., no flies at temperatures >30°C; Lanzoni et al. 2002, Maharjan and Jung 2016). In contrast, growth chamber studies of physiological tolerance for L. sativae on P. vulgaris showed a broader temperature range, ranging from 10 to 40°C (Fig. 2a;Haghani et al. 2007).

Our findings support observations that L. huidobrensis is primarily a tropical montane species which has expanded its niche to its worldwide distribution. Multiple independent studies of the susceptibility to high temperatures suggest that physiological tolerance of L. huidobrensis to high temperatures is a trait with low plasticity or variability available for selection, whereas its adaptation to cooler temperatures has been documented (Chen and Kang 2004, Weintraub et al. 2017). Similar convergence of low adaptation for high temperatures and better adaptation to low temperatures has been found for the fly genus Drosophila (Hoffmann 2010). Assays of thermotolerance hardening have been performed showing it may increase heat resistance on the leafminers (Huang et al. 2007); however, it is difficult to discern from these assays how heat shocks will occur and how individual populations will respond to climatic changes in the field. This poses a question of how L. huidobrensis will respond to the future climate scenarios with temperature increases in most locations. L. huidobrensis will likely shift its range to more temperate regions, and it may be prevalent for longer periods of time in temperate zones with milder winters. In tropical and sub-tropical regions, they may be replaced by other species of leafminers that cope better with higher temperatures (such as L. sativae and L. trifolii) except for the mountain peaks which may retain suitable conditions for L. huidobrensis to thrive. Our results do not support climate envelope-based invasive scenarios modeled by Mika and Newman (2010) in which most of the areas plotted to be suitable in North America would be too hot for viable populations of L. huidobrensis to establish. This is not surprising because developmental studies of ectotherms has shown that knowledge of microclimatic conditions and individual responses to climate change is necessary in order to accurately predict future distributions (Pincebourde and Casas 2015).

Regarding international trade, a large proportion of fresh horticultural exports from Central America to the USA originate in highland regions where L. huidobrensis thrives and outcompetes L. sativae. For this reason, L. huidobrensis is the main invasive threat, especially since Scheffer and Lewis (2001) showed that L. huidobrensis is genetically distinct from its closest relatives in California (now considered Liriomyza langei, formerly huidobrensis). However, the fact that Miami, Florida is the main port of entry with maximum daily temperatures of 30°C or higher during much of the year (April–October) should preclude establishment of L. huidobrensis in south Florida in all but the coolest months. Moreover, populations that might establish when maximum daily temperatures are below 28 °C (November–March) would go extinct as soon as daily highs surpass this critical temperature.

Results from the ENM for L. sativae suggest that this species has a wider tolerance to temperatures and a larger potential distribution niche than L. huidobrensis in Guatemala (Figs. 3a, 2a). A wider range in temperature tolerance for L. sativae has also been found in studies from Iran (Haghani et al. 2007). The strongest predicting variables according to the ENM of L. sativae were PPT of the coldest quarter (Supp Table 3 [online only], which in Guatemala coincides with the dry months), PPT of the driest month and MTWQ; this suggests that desiccation may be more of a limiting factor for L. sativae than the actual environmental temperature. Interestingly the presence of L. sativae was determined by the probability of presence of L. huidobrensis.

Even though we did not experimentally test for the temperature dependent exclusion interaction that L. huidobrensis has on L. sativae, it is important to note that the elevation gradients we sampled with a wide array of temperatures (i.e., sampled MAXTEMP range 18–36°C) are part of the native niche of the two species. This acted as a “natural experiment” in which we could test if the presence of one species reduces the probability of finding the other: in our study L. huidobrensis had a small area of unique habitat in the eastern mountains of Guatemala (Fig. 4b) although L. sativae had the majority of its range free from the presence of L. huidobrensis. Hence, it is important to note that while L. sativae has the largest unshared niche, it is the species that gets displaced by L. huidobrensis, whereas for L. huidobrensis temperature is the major factor behind its distribution across a temperature gradient. This pattern was also observed in the mountains of China, where it has been reported that L. sativae was present throughout highlands of China up until 1990 when L. huidobrensis was introduced (Chen and Kang 2004, Kang et al. 2009). Now L. huidobrensis is an important pest of the temperate and highland sites of China and the distribution of L. sativae has changed to valley and coastal regions. Further, the same spatial displacement was observed in the mountains of Java, Indonesia where L. huidobrensis dominates at elevations above 1,000 m and L. sativae is abundant from 0 to 1,000 m (Tantowijoyo and Hoffman 2010). We argue that these observations of ecological interactions both in native and introduced distribution ranges point to a temperature mediated competitive exclusion in which the physiological thermotolerance constraint of one of the species (L. huidobrensis) affects the presence of another species (L. sativae). This phenomenon of competition mediated by environmental factors has, so far, only been fully documented in fish (Taniguchi et al. 1998).

L. huidobrensis and L. sativae are well-known polyphagous species, feeding on hundreds of host plant species. However, it has been recorded that their population growth varies depending on the host plant (Shepard and Braun 1998). We found support for not all host plants being equally good hosts in Guatemala for the populations of L. huidobrensis and L. sativae. While L. sativae thrived on tomato and beans, we found that L. huidobrensis does well on snowpeas. Unfortunately, our low sample size on snowpeas (n = 7) prevented us from formulating trophic cascade hypotheses. However, previous findings confirm that the highest pupae production per female and the highest pupal viability occurred when L. huidobrensis fed on snowpeas (MacVean and Pérez 1997). The emergence of L. huidobrensis does not appear to be independent of temperature as the highest emergence occurred at cultivars with the lowest maximum temperatures (21.1°C). Low parasitism rates (4.14 ± 2.17, Supp Table 1 [online only]) may also aide emergence of L. huidobrensis on snowpeas. Regardless of host plant species L. huidobrensis was limited by the MAXTEMP of the location in which it was collected. The highest pupal emergence occurred at the lowest temperatures and the extent to which temperature affected its abundance did change with host plant (path coefficients: −0.33 to −0.72; Fig. 7). Endoparasitism in L. huidobrensis was highly influenced by host plant with the highest effects on lettuce and onions (path coefficients = 0.68 and 0.54, respectively; Fig. 7). This may reflect a host plant preference by the parasitoids due to the ease of finding herbivores on these plants (Perfecto and Vet 2003).

Fig. 7.

Path analysis of the interactions between temperature, Liriomyza huidobrensis, Liriomyza sativae and parasitoids at four different host plants where they coexist and for which we had >20 samples. Solid lines depict direct effects; dashed lines depict indirect effects. Pointed arrows are for positive effects and bullet ended lines depict negative effects. Path coefficients are given for each arrow. (a) Lettuce (Lactuca sativa), (b) onion (Allium cepa), (c) tomato (Solanum lycopersicum), and (d) beans (Phaseolus vulgaris). Chi-square and P values for each model are given.

Other studies have found that despite the wide array of suitable hosts for L. sativae, tomato remains its preferred host plant (Zehnder and Trumble 1984); furthermore, olfactory assays show that preference for tomatoes is inherited and does not change when feeding on beans for generations (Zhao and Kang 2003). The degree that L. sativae benefited from higher temperatures varied greatly with a low benefit from temperature on tomato (Fig. 7). This may be due to the overlap between optimal developmental temperatures of L. sativae (Fig. 2a) and the temperatures at which tomato is grown in Guatemala (Fig. 2b). Interestingly, the effect of parasitism on L. sativae was also low when feeding on tomato host plants (mean parasitism = 23.7 ± 3). This together with the low parasitism rates on L. huidobrensis on snowpeas suggest that even a little reduction in pressure from parasitoids could have large impacts on the abundance of herbivores.

We found that across all host plants the abundance of L. huidobrensis had a strong negative effect on the abundance of L. sativae (Fig. 7, path coefficients = 0.33–0.55). Hence there was an indirect positive effect of temperature on L. sativae through its negative action on L. huidobrensis populations. Through path analysis we were able to tease apart the temperature mediated interaction between the two species and found that there are opposite direct effects of temperature on the leafminer flies (negative on L. huidobrensis and positive on L. sativae), which corresponds with their local and global distribution records (Fig. 3; Chen and Kang 2004, Tantowijoyo and Hoffman 2010). We however also found an indirect positive effect of temperature on L. sativae as competition between species is strong, such that lower emergence of L. huidobrensis’ at high temperatures allows for abundant growth in L. sativae.

The negative effect of parasitoids on L. sativae (path coefficients = −0.11 to −0.35; Fig. 7) were in general lower than on L. huidobrensis. Interestingly we also found an indirect positive effect of parasitoids on L. sativae; the strong regulating effect of parasitism on L. huidobrensis benefited abundance of L. sativae by releasing it from strong interspecific competition. These hypotheses suggest that under future climatic scenarios the distribution of these cosmopolitan species will depend on temperatures, which host plants are abundant at the locations and the success of parasitoids in regulating leafminer hosts. A possible scenario is that L. sativae will expand its range to new areas at higher elevations once they become too warm for L. huidobrensis. Populations of L. huidobrensis are likely to track their cooler climate niche and migrate to higher elevations. In the case of latitudinal shifts, L. huidobrensis could establish at latitudes that have not had invasions. In addition, populations of this leafminer that would normally die after a cold winter may now have active populations year-round.

In conclusion, we provide a strong framework to study invasion risk and study potential response scenarios of climate change by first testing physiological limits of the species instead of inferring limits by the locations where the species is found. As an important implication for trade policy regarding invasive species, we provide evidence for a low risk of establishment by L. huidobrensis in the USA through Miami, Florida based on temperature incompatibility. Second, we characterize ecological interactions between two leafminer species to address how certain ecological interactions like host plant identity, interspecific competition and top-down forces through parasitism enable or limit abundance of a cosmopolitan species in its native range. Finally, through path analysis we provide clear hypotheses scenarios of host plant, parasitoid, interspecific interaction and changes in temperature to be tested with experiments in the future.

Supplementary Data

Supplementary data are available at Journal of Insect Science online.

Click here for additional data file.