Can Niche Dynamics and Distribution Modeling Predict the Success of Invasive Species Management Using Biocontrol? Insights From Acacia longifolia in Portugal

fevo-08-576667 October 14, 2020 Time: 16:43 # 1 ORIGINAL RESEARCH published: 19 October 2020 doi: 10.3389/fevo.2020.576667 Edited by: Rubén G. Mateo, Autonomous University of Madrid, Spain Reviewed by: Aitor Gastón, Polytechnic University of Madrid, Spain Tiago Vasconcelos, São Paulo State University, Brazil *Correspondence: Marco Dinis [email protected] Joana R. Vicente [email protected] Elizabete Marchante [email protected] Hélia Marchante [email protected] Specialty section: This article was submitted to Biogeography and Macroecology, a section of the journal Frontiers in Ecology and Evolution Received: 26 June 2020 Accepted: 25 September 2020 Published: 19 October 2020 Citation: Dinis M, Vicente JR, César de Sá N, López-Núñez FA, Marchante E and Marchante H (2020) Can Niche Dynamics and Distribution Modeling Predict the Success of Invasive Species Management Using Biocontrol? Insights From Acacia longifolia in Portugal. Front. Ecol. Evol. 8:576667. doi: 10.3389/fevo.2020.576667 Can Niche Dynamics and Distribution Modeling Predict the Success of Invasive Species Management Using Biocontrol? Insights From Acacia longifolia in Portugal Marco Dinis1,2,3*, Joana R. Vicente2*, Nuno César de Sá1,3,4, Francisco A. López-Núñez3, Elizabete Marchante3*and Hélia Marchante1,3* 1Centre for Functional Ecology - Science for People & the Planet, Instituto Politécnico de Coimbra, Escola Superior Agrária de Coimbra, Coimbra, Portugal, 2Centro de Investigação em Biodiversidade e Recursos Genéticos (CIBIO/InBIO), Universidade do Porto, Porto, Portugal, 3Centre for Functional Ecology - Science for People & the Planet, Department of Life Sciences, University of Coimbra, Coimbra, Portugal, 4Institute of Environmental Sciences, Leiden University, Leiden, Netherlands Biological invasions are among the major contemporary threats to biodiversity. Biocontrol has a long history as a safe and effective strategy for the control of invasive species in several world regions, yet the life history and ecological requirements of biocontrol agents are often poorly understood. Species distribution models and assessments of niche overlap and dynamics constitute a way to quantify and compare ecological niches and are widely used for predicting biological invasions. While these tools can also be valuable to predict the effectiveness of biocontrol programs and anticipate whether biocontrol agents can establish on areas targeted for control, they remain underutilized for this purpose. In this work, species distribution models and comparisons of niche dynamics are used to predict the success of the ongoing biocontrol program for Acacia longifolia in Portugal, and potentially along the Mediterranean Basin, using the Australian gall-forming wasp Trichilogaster acaciaelongifoliae, previously released in South Africa. Niches of both the invasive plant and the biocontrol agent were found to be highly similar and to remain conserved through the introduction process. Distribution models identify suitable climatic areas for A. longifolia in 19% of the Mediterranean Basin and predict successful establishment of T. acaciaelongifoliae in 41% of the suitable area for A. longifolia, despite excluding a few locations where very recent establishment occurred. These results allow us to quantify the risk of future A. longifolia invasion and potential success of biocontrol, as well as establish a comparative framework for similar programs being considered in other regions of the world dealing with A. longifolia invasions. Keywords: biological invasions, biocontrol, distribution modeling, niche conservatism, Acacia longifolia, Trichilogaster acaciaelongifoliae Frontiers in Ecology and Evolution | www.frontiersin.org 1October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 2 Dinis et al. Niche-Based Biocontrol Management of Acacia INTRODUCTION Biological invasions are among the major contemporary threats to biodiversity (IPBES, 2019). Specifically, invasive plants cause changes in ecosystem services and functionality, with negative impacts at different levels, e.g., at the level of plant and soil communities and nutrient dynamics (Marchante et al., 2008a,b; Le Maitre et al., 2011;Vilà et al., 2011;Simberloff et al., 2013; Zhang et al., 2019) which are often difficult to reverse or mitigate without active restoration efforts (Le Maitre et al., 2011). Control of invasive species is frequently prohibitively expensive and labor-intensive (Marais et al., 2004). Understanding biological invasion processes is essential for informing decision-making and optimizing the allocation of limited resources. Correlative approaches based on geographic occurrence data (such as species distribution models and multivariate environmental assessments of niche overlap) are valuable to quantify and compare ecological niches and are widely used to predict and improve management of biological invasions (e.g., Broennimann et al., 2007;Vicente et al., 2010, 2011, 2016;Petitpierre et al., 2012;Guisan et al., 2014). Biological control (hereafter biocontrol) has a long history (over a century) as a safe and effective strategy for the control of invasive plant species (Clewley et al., 2012) in several regions of the world, such as South Africa, Australia, or United States. Yet, the biocontrol of invasive plants was only recently (less than a decade ago) initiated in Europe (Shaw et al., 2018). Classical biocontrol consists on the deliberate introduction of a natural enemy of the target invasive species into a new region invaded by the latter with the objective of mitigating its competitive advantage. Biocontrol agents should be highly specific in order to minimize interactions with native biodiversity in the regions where they are introduced. They should also be able to establish viable populations in the invaded areas. Considering this, it follows that the chances of success of a biocontrol agent are likely to be improved when the ecological requirements for its successful establishment and long-term persistence are similar to those of its target. These requirements, due to their explicitly geographically-oriented and scenopoetic nature, fall within the Grinnellian definition of the ecological niche (sensu Peterson et al., 2011). Consequently, the similarity in ecological requirements between species can be characterized using the related concept of niche overlap (sensu Warren et al., 2008). Despite this, very few studies so far have explicitly compared biocontrol agents and their respective targets in terms of niche overlap (but see Sun et al., 2017). Acacia longifolia is a Leguminosae native to southeastern Australia which has become highly invasive in several exotic locations. Currently, invasive populations of A. longifolia are established beyond its native range in south-western Australia (Costello et al., 2000;Impson et al., 2011), as well as in New Zealand, South Africa, Chile, United States (California), Spain (Galicia), and Portugal (EFSA Panel on Plant Health, 2015). Large long-lived seed banks and the ability to germinate following disturbances (particularly fire) are common traits in the Acacia genus which account for its high invasive ability (Richardson and Kluge, 2008;Passos et al., 2017). Invasions by A. longifolia are associated to changes in soil communities and nutrient cycles (Marchante et al., 2008a,b), disturbance of plant (Costello et al., 2000;Marchante et al., 2003, 2015), and plantgall communities (López-Núñez et al., 2017), nutrient and water cycling (Werner et al., 2010), and fire regimes. Many of these effects tend to persist after removal, hindering restoration efforts (Marchante et al., 2009, 2011a). Mechanical and chemical control of invasive Acacia species is often highly demanding in terms of resources, time and labor (Marais et al., 2004). Reduction of seed banks is particularly problematic, yet essential for the longterm management of invasive plants with high seed production (Richardson and Kluge, 2008). Biological control is consequently often proposed as a valuable and indispensable asset for the management of Acacia invasions (Marais et al., 2004;Richardson and Kluge, 2008;Wilson et al., 2011). Trichilogaster acaciaelongifoliae is a gall-forming wasp endemic to Australia. Like all gall-forming insects, it is a highly specific endoparasite which requires successful interaction with its host to complete the life cycle, making it an ideal candidate for biocontrol programs due to the extremely low risk of detrimental interactions with species other than its preferred host, and its inability to expand beyond the distribution of its host (Dennill et al., 1993;Marchante et al., 2011b). Currently, only two programs using T. acaciaelongifoliae to control A. longifolia exist in the world. It was first introduced in South Africa in 1982 and 1983 (Dennill, 1985, 1988, 1990;Dennill and Donnelly, 1991; Dennill et al., 1993) where it quickly became established through most of the range of A. longifolia, successfully decreasing the reproductive potential (89–95%) and vegetative growth (53%) of its target and causing A. longifolia mortality in stressful situations (Dennill, 1985, 1988;Dennill and Gordon, 1990). Owing to the success of A. longifolia biocontrol in South Africa, a similar program was initiated in Portugal (Marchante et al., 2011b), with the first release of T. acaciaelongifoliae in the wild taking place recently, in 2015 (Marchante et al., 2017). As of 2018, four populations of T. acaciaelongifoliae established in Portugal, and it is expected that they will maintain the current trend of population growth and expansion, eventually reaching all areas of the country invaded by A. longifolia (Marchante et al., 2017). While the program was initiated in Portugal, the agent may potentially spread to other regions of the Mediterranean Basin where adequate environmental conditions are available. Understanding how the niches of a biocontrol agent and its target species relate in their native distribution may be insufficient to predict the success of biocontrol. The possibility of niche shifts during introduction into new environments (Guisan et al., 2014) may promote niche divergence between the two species, creating a mismatch in environmental space which, when transported to geographic space, can result in inability of the biocontrol agent to establish in areas targeted for control. This work aims to characterize the realized niches of A. longifolia and T. acaciaelongifoliae and test the adequacy of T. acaciaelongifoliae as a biocontrol agent for A. longifolia in Portugal as a function of niche overlap. To accomplish this, we characterize and compare the realized niches of the host plant and of the biocontrol agent in their native and introduced ranges. We also test the hypothesis of niche shifts by comparing realized niches between native and introduced ranges for both species. Frontiers in Ecology and Evolution | www.frontiersin.org 2October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 3 Dinis et al. Niche-Based Biocontrol Management of Acacia We then project the environmental affinities of A. longifolia and T. acaciaelongifoliae into geographical space using species distribution models, in order to predict the success of the biocontrol program primarily in Portugal but ultimately in the geographically and climatically proximate Mediterranean Basin. MATERIALS AND METHODS Occurrence Records A database of A. longifolia and T. acaciaelongifoliae occurrence records was compiled from the Global Biodiversity Information Facility1, Atlas of Living Australia2, South Australian Museum3 provided upon request by Doctor Peter Hudson, Australian National Insect Collection4provided upon request by Doctor Juanita Rodriguez, Invasive Plants in Portugal – Invasoras.pt5 (Marchante et al., 2017), Southern Africa Plant Invaders Atlas (Henderson, 1999), literature references (Dennill, 1985, 1987; Prinsloo and Neser, 2007;Lado, 2008;Henriksen et al., 2017, 2019); and data provided by authors and collaborators (see section “acknowledgments”). Data collection for A. longifolia and T. acaciaelongifoliae was independent and no records were shared among the two datasets. Occurrence data were considered until June 2018. A. longifolia records were included taking into account the species sensu lato (i.e., including subspecies A. longifolia longifolia,A. longifolia floribunda and A. longifolia sophorae, all considered hosts of T. acaciaelongifoliae). In June 2018, T. acaciaelongifoliae was still narrowly distributed in Portugal, and thus occurrence records in this area were not considered for subsequent analyses. Records outside the native area (southeastern Australia) and the introduced areas of interest (South Africa and Portugal) were excluded. Uncertain or duplicate records were also excluded. Records with less positional accuracy than 10 km were excluded to maximize geographic coverage while minimizing the loss of spatial precision (Gutiérrez-Rodríguez et al., 2017). To reduce potential sampling bias and spatial autocorrelation, a declustering approach (e.g., Vale et al., 2016;Dinis et al., 2019) was performed using the R package ecospat (Di Cola et al., 2017), which consisted of delimiting a 10 km buffer around each record and, when two or more buffers overlapped, removing one of the records randomly. This process was repeated until all remaining records per species had a distance of at least 10 km among them. Clustering was quantified before and after the declustering treatment to ensure an effective decrease in the level of clustering (e.g., Dinis et al., 2019). This was done by calculating the Nearest Neighbor Index in ArcMap 10.1 (Environmental Systems Research Institute, 2012). The final database consisted of 1493 and 140 occurrences for A. longifolia and T. acaciaelongifoliae, respectively (Supplementary Table 1). Occurrences were then 1https://www.gbif.org/en/ 2https://www.ala.org.au/ 3http://www.samuseum.sa.gov.au/ 4https://www.csiro.au/en/Research/Collections/ANIC 5http://www.invasoras.pt further subdivided by geographical areas of interest (southeastern Australia, South Africa, and Portugal; Figure 1). Bioclimatic Variables and Study Area The native distributions of A. longifolia and T. acaciaelongifoliae were represented by southeastern Australia. South Africa and Portugal were selected to represent areas of A. longifolia invasion because they are the areas where biocontrol programs with T. acaciaelongifoliae are ongoing. The Mediterranean Basin (slightly expanded to include also northern Spain) was selected for model projections for two reasons. First, it is the geographical area of which Portugal is part and A. longifolia occurs with invasive potential in several countries in the region (Spain, France, Italy, Lorenzo et al., 2010); at the same time, if A. longifolia is present, this creates the possibility of natural dispersion of T. acaciaelongifoliae. Second, climatic similarity with the global distribution of A. longifolia and T. acaciaelongifoliae (predominately in areas of Mediterranean climate) is expected to maximize analogy of predictors. Nineteen bioclimatic variables representing average contemporary conditions (1970–2000) were obtained from Worldclim v2.06(Hijmans et al., 2005) at a spatial resolution of 10 ×10 km, to ensure consistency with the minimum positional accuracy of the occurrence data. Other scenopoetic candidate variables, such as topography, were considered but ultimately excluded under the assumption that they would not be informative at this spatial resolution, while variables related to soil composition and vegetation cover/structure were excluded due to being dynamically linked to the presence of the target species, consequently falling outside the scope of the Grinnellian niche (Peterson et al., 2011). The temporal extent of the Worldclim dataset overlaps with the collection dates for most occurrence records used in this work, and was therefore selected as an adequate representation of climatic conditions available to the target species upon the time of observation. Background and training areas were defined by clipping bioclimatic layers to a 200 km buffer surrounding the minimum convex polygon which includes all occurrence records. Buffer size was selected taking into account coarse distributional patterns, topography, vegetation zones, and the limits of climatic regions according to the Köppen-Geiger climate classification (Köppen, 1900) and aims to represent the habitats and areas which the species may have reasonably sampled in each region (Barve et al., 2011) and prevent biases in the sampling of background data, which are known to occur when using overly broad or restrictive training areas (VanDerWal et al., 2009;Anderson and Raza, 2010). This was done separately for each species in each relevant area of distribution, i.e.: for A. longifolia in its native range in southeastern Australia as well as the invasive ranges in South Africa and Portugal; and for T. acaciaelongifoliae in its native and introduced ranges in southeastern Australia and South Africa, respectively. Analyses were also performed combining the multiple native and introduced/invaded ranges for each species (Australia, South Africa and Portugal for A. longifolia, Australia and South Africa for T. acaciaelongifoliae, 6http://www.worldclim.org Frontiers in Ecology and Evolution | www.frontiersin.org 3October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 4 Dinis et al. Niche-Based Biocontrol Management of Acacia FIGURE 1 | Study areas and occurrence records of Acacia longifolia and Trichilogaster acaciaelongifoliae used in this work. (A; top right): Portugal; (B; bottom left): South Africa; and (C; bottom right): southeastern Australia. Trichilogaster acaciaelongifoliae records from Portugal are merely indicative and were not included in any analyses due to their low number and recent establishment. Figure 2). All operations were performed in ArcMap 10.1 (Environmental Systems Research Institute, 2012). Niche Overlap Ecological niches were compared between A. longifolia and T. acaciaelongifoliae, as well as between different native and introduced areas of each species’ distribution (Figure 2). For A. longifolia, niches were compared between the native range in Australia and the two invaded ranges in South Africa and Portugal. For T. acaciaelongifoliae, comparisons took into account only the native and introduced distributions in Australia and South Africa, respectively. All three distribution areas for A. longifolia were compared, to avoid assumptions about introduction pathways. Niches were compared in multivariate environmental space using the PCA-env technique (Broennimann et al., 2012). This approach allows for pairwise comparisons of niches between groups, and consists of three steps: (1) calculation of density of occurrences along a multivariate environmental space for both target groups; (2) measurement of the overlap between density distributions of the two groups in environmental space (niche overlap), as defined by Schoener’s D metric (Schoener, 1970); and (3) statistical tests of niche equivalency (whether niche overlap is constant when randomizing occurrences among both groups) and niche similarity (whether the niche of one group is more similar to the other than would be expected by chance given the available environmental background), using a permutationbased framework (Warren et al., 2008;Broennimann et al., 2012). Representation of niches as densities along environmental space allows the quantification of differences in the densities of groups, which can be expressed as niche stability (proportion of niche B overlapping with niche A), niche expansion (proportion of niche B non-overlapping with niche A), and niche unfilling [unique proportion of niche A (i.e., non-overlapping with niche B); Guisan et al., 2014]. These metrics were developed to compare Frontiers in Ecology and Evolution | www.frontiersin.org 4October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 5 Dinis et al. Niche-Based Biocontrol Management of Acacia FIGURE 2 | Workflow for niche overlap assessments and distribution modeling. Panel (A; top): Regional species distribution datasets were compared to assess niche conservatism for Acacia longifolia and Trichilogaster acaciaelongifoliae. Niche overlap among species was quantified (overlap), tested (equivalency and similarity), and characterized (stability, expansion, and unfilling) to assess ecological overlap in environmental space. Panel (B; bottom): Distribution modeling was used to project ecological overlap into geographic space and predict the potential ability of T. acaciaelongifoliae to establish in areas susceptible to A. longifolia invasion in the Mediterranean basin. Acronyms: AU: Australia, SA: South Africa, and PT: Portugal. native and exotic niches, and represent the proportion (0–1) of the realized niche that remains constant (stability), is gained (expansion), or is lost (unfilling) during introduction/invasion. Application outside of an explicit comparison of native vs introduced ranges (such as when comparing two unrelated taxa like A. longifolia and T. acaciaelongifoliae) requires an arbitrary assignment of each group to either niche A or B, which will alter the order in which niche expansion and unfilling are calculated but will not affect the final results or interpretation. Background areas for each group were spatially delimited using the same areas used for calibration of niche models. For each pairwise comparison, background areas of the two groups under comparison were merged. Environmental space was defined by the 19 bioclimatic variables available from Worldclim v2.0 (Hijmans et al., 2005), at a resolution of 10 ×10 km. Though several of these variables are correlated in geographic space, the PCA-based approach addresses this by creating an orthogonal linear combination of the original variables. This approach allows the inclusion of all available bioclimatic variables, maximizing the environmental variance under consideration for niche assessments and minimizing assumptions about variable importance. Given what little is known about the ecology of T. acaciaelongifoliae, we opted to minimize such assumptions. Niche overlap was calculated and niche equivalency and similarity tests were performed among regions for each species, and across species. We also calculated metrics of niche stability, expansion and unfilling using the same framework. For all tests among species, A. longifolia was considered as group A and T. acaciaelongifoliae as group B. All analyses were performed in the R package ecospat (Di Cola et al., 2017). Species Distribution Models Based on the evidence for niche conservatism for both A. longifolia and T. acaciaelongifoliae (see Results section), species distribution models were calibrated using all available records for southeastern Australia, South Africa and (for A. longifolia) Portugal. This allows maximization of the number of occurrence records used to train the models and a better sampling of relevant environmental gradients. To avoid overparametrization, distribution models were developed using a subset of the bioclimatic variables used for niche overlap analysis. This approach precludes direct comparison between the two methods. Consequently, we treat the distribution models strictly as niche-based tools for predicting distributions rather than accurate representations of the niches, relying solely on the niche overlap tests for any inferences regarding the dynamics and relationships of niches. Spatial correlation between bioclimatic Frontiers in Ecology and Evolution | www.frontiersin.org 5October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 6 Dinis et al. Niche-Based Biocontrol Management of Acacia variables was assessed in ArcMap 10.1 (Environmental Systems Research Institute, 2012), and multicollinearity was assessed by calculating the Variance Inflation Factor (VIF) using the usdm package in R (Naimi et al., 2014). When two or more variables were highly correlated, preference was given to the one with the most biological sense. A set of five and six slightly correlated (Pearson’s r<0.7) and non-collinear variables (VIF <5; e.g., Vicente et al., 2013) were selected to create models for T. acaciaelongifoliae and A. longifolia, respectively, (Table 1). An ensemble modeling approach was performed in the R package biomod2 (Thuiller et al., 2009) using a combination of seven modeling techniques: Generalized Linear Models (GLM), Generalized Boosting Models (GBM), Flexible Discriminant Analysis (FDA), Classification Tree Analysis (CTA), Multivariate Adaptive Regression Splines (MARS), Random Forests (RF), and Maximum Entropy (MAXENT). Models were created using 5 datasets of 1000 pseudoabsences randomly distributed through the training areas, and 8 replicates were performed for each combination of pseudoabsence dataset and modeling algorithm, for a total of 280 models per species. Each model was calibrated setting aside 30% of available presence records for model evaluation. The resulting testing datasets were used to evaluate model discrimination by calculating the area under the curve of the receiver operating characteristic curve (AUC), and ensemble models across all modeling techniques were created for each species (e.g., Petitpierre et al., 2017), using the mean of all models where AUC >0.7. Individual models were projected for the Mediterranean Basin and binary ensemble models were created using the Maximized Sum Threshold criteria (MST; Cantor et al., 1999), implemented in biomod2 under the “binary.meth = ‘ROC”’ argument (Thuiller et al., 2009). This method selects the threshold which maximizes the sum of sensitivity and specificity. We opted to use this approach due to its good general performance in comparative studies and particularly for its superior performance with low-prevalence datasets and ability to minimize omission errors, which are generally more costly in conservation applications such as identifying areas under risk of biological invasion (Liu et al., 2005;Jiménez-Valverde and Lobo, 2007). The importance of environmental variables for each model was determined by average permutation importance (Phillips et al., 2006). RESULTS Niche Characterization Niche overlap among A. longifolia was overall low between the three regions (less than 0.5, with 0 corresponding to no overlap and 1 to complete overlap), with the greatest overlap occurring between A. longifolia in South Africa and Portugal (0.416; Table 2). Niche equivalency, i.e., whether niche overlap remains constant when randomizing occurrences among groups, was found in all pairwise tests except for A. longifolia in Australia vs. Portugal and South Africa vs. Portugal. Niche similarity, i.e., whether groups are more similar than expected considering available background, was found in all pairwise tests. Similar niche dynamics were found for all A. longifolia comparisons, with predominant niche stability (>0.89), analogous niches across regions and very little expansion into novel areas of environmental space during introduction. Niche unfilling, i.e., areas of environmental space present in original range but not after introduction, ranged between 0.2 and 0.26 for all pairs except South Africa vs. Portugal, which had very little unfilling (<0.01). For T. acaciaelongifoliae in Australia vs. South Africa, niche unfilling was of similar magnitude to A. longifolia (0.25), but niche expansion was much more pronounced (0.56), with a corresponding decrease in niche stability (0.44; Table 2). For comparisons between A. longifolia and T. acaciaelongifoliae, niche overlap ranged between 0.53 (in Australia) and 0.85 (in South Africa). The niches of the two species were found to be equivalent and similarity was not rejected in any pairwise test. Niche stability was higher than 0.96 in all comparisons, with corresponding low values of niche expansion from A. longifolia to T. acaciaelongifoliae (<0.04). Niche unfilling from A. longifolia to T. acaciaelongifoliae was lowest in South Africa (0.05) and highest in Australia (0.23; Table 2). All significant (p-value <0.01) niche similarity tests corresponded to values of overlap greater than the distribution TABLE 1 | Variable ID and minimum, maximum, and mean values across model training areas for Acacia longifolia and Trichilogaster acaciaelongifoliae of ecogeographical variables used in model creation. Acacia longifolia Trichilogaster acaciaelongifoliae Min Max Mean Min Max Mean BIO2 Mean diurnal range (C◦) 5.28 18.04 12.37 BIO3 Isothermality (%) 33.02 62.67 48.95 36.18 63.39 49.77 BIO4 Temperature seasonality 201.10 634.90 426.80 BIO6 Min temperature of coldest month (C◦)−6.15 13.42 3.18 BIO9 Mean temperature of driest quarter (C◦) 0.50 25.03 15.93 BIO13 Precipitation of wettest month (mm) 20.00 299.00 91.38 BIO14 Precipitation of driest month (mm) 1.00 119.00 26.93 2.00 119.00 25.85 BIO15 Precipitation seasonality 8.55 87.14 39.22 BIO19 Precipitation of coldest quarter (mm) 11.00 842.00 174.60 Variables are represented only for the species in which they were utilized for model creation. Frontiers in Ecology and Evolution | www.frontiersin.org 6October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 7 Dinis et al. Niche-Based Biocontrol Management of Acacia of pseudoreplicates, indicating that the groups are more similar than expected considering the available background. Species Distribution Models All 280 models per species had AUC values above 0.7 and were thus included in the ensembles. AUC values for A. longifolia models ranged between 0.79 and 0.94, with a weighted average of 0.90 ±0.02 (±SD). AUC of models for T. acaciaelongifoliae ranged between 0.71 and 0.92, with a weighted average of 0.85 ±0.04 (±SD; Table 3). Temperature seasonality (BIO4), minimum temperature of the coldest month (BIO6), and precipitation of the coldest quarter (BIO19) were the most important variables for the A. longifolia models, with marginal contributions from the remaining variables. Variable importance for T. acaciaelongifoliae was less straightforward, with all variables contributing somewhat to the models, yet isothermality (BIO3) and mean diurnal range of temperature (BIO2) were the variables with the most explanatory power (Table 3). Suitable areas for A. longifolia were predicted through the Atlantic coast of the Iberian Peninsula and Morocco, in some areas extending up to 280 km inland (Figure 3). Suitability for A. longifolia was also identified in northern Spain, southern France, some areas of the Italian and Balkan peninsulas (particularly along the Adriatic coast), most of the Mediterranean coast of North Africa and the Middle East and all Mediterranean archipelagos. T. acaciaelongifoliae is predicted to have suitable climatic conditions in most of the Iberian coastline (except for a gap of approximately 300 km in northern Portugal), the entirety of the coastline of Morocco, most of the Mediterranean areas of Libya, Egypt, Israel, and Palestine and areas in the Adriatic coast and in the islands of Cyprus, Sicily and the Aegean sea (Figure 3 and Supplementary Figure 1). In total, 19% of the total area of the Mediterranean basin was predicted to be suitable for A. longifolia, and 41% of the area identified as suitable for A. longifolia is also predicted as suitable for T. acaciaelongifoliae. Currently established population of T. acaciaelongifoliae were all correctly identified as suitable by the modeling approach used, with the exceptions of the population located in the city of Coimbra (Figure 3 and Supplementary Figure 1) and a few very recently observed populations in the northern coast (not shown). DISCUSSION Niche Dynamics of A. longifolia and T. acaciaelongifoliae While there are differences in niche dynamics between different ranges of A. longifolia and T. acaciaelongifoliae, the results are overall consistent with niche conservatism during the invasion/introduction process, respectively. Evidence for niche similarity on all cases suggests that observed differences in niches are related to differences in habitat availability and/or introduction pathways between Australia, South Africa and Portugal. Comparisons between Australia and South Africa revealed contrasting patterns for the two species, with A. longifolia invasion being associated with niche stability, and T. acaciaelongifoliae introduction being associated with strong niche expansion. It’s unlikely that this reflects actual expansion of T. acaciaelongifoliae into new environments in TABLE 2 | Niche overlap and dynamics between all pairs of native (AU: Australia) and introduced (SA: South Africa; PT: Portugal) Acacia longifolia (Al) and Trichilogaster acaciaelongifoliae (Ta). Pairs PC1 PC2 Niche overlap (D) Equivalency (sig = non-eq.) Similarity A- >B Similarity B- >A Expansion Stability Unfilling Acacia longifolia AU-PT 44.27 20.19 0.18 0.001*0.052 0.07 0.084 0.916 0.257 AU-SA 37.24 26.64 0.33 0.263 0.012 0.009*0.102 0.898 0.211 SA-PT 45.01 19.67 0.416 0.001*0.011 0.013 0.06 0.94 0.004 Trichilogaster acaciaelongifoliae AU-SA 37.24 26.64 0.157 0.987 0.141 0.134 0.561 0.439 0.254 Among species AlSA-TaSA 35.02 27.08 0.846 0.604 0.004*0.005*0.039 0.961 0.053 AlAU-TaAU 48.96 22.37 0.526 0.94 0.094 0.084 0.001 0.999 0.227 AlAll-TaAll 39.71 24.08 0.556 1 0.09 0.09 0.013 0.987 0.106 Comparisons are made within each species and between the two species considering the native (AU) and introduced (SA) ranges where both occur. Tests between the two species considering the full range where they occur (AU, SA and PT for A. longifolia; AU and SA for T. acaciaelongifoliae) are also presented (AlAll-TaAll). Percentage of explained variation of the two first principal component (PC1 and PC2), niche overlap (Schoener’s D), p-values of equivalency and similarity tests (p <0.01, significant values marked with *) and indices of niche expansion, stability and unfilling are presented. All tests of niche similarity presented null distributions lower than the observed overlap, indicating niche similarity. TABLE 3 | Number of occurrence records (N), average (and SD) Area Under the Curve (AUC) and average (and SD) permutation importance of each ecogeographical variable across all 280 species distribution model replicates for Acacia longifolia and Trichilogaster acaciaelongifoliae. ∗Indicate highest contributing variables. N AUC BIO3 BIO4 BIO6 BIO9 BIO14 BIO19 Acacia longifolia 1493 0.897 (0.022) 0.045 (0.033) 0.238* (0.089) 0.178* (0.089) 0.072 (0.052) 0.042 (0.051) 0.350* (0.094) N AUC BIO2 BIO3 BIO13 BIO14 BIO15 Trichilogaster acaciaelongifoliae 140 0.846 (0.040) 0.501* (0.109) 0.527* (0.113) 0.251 (0.153) 0.343 (0.296) 0.261 (0.236) Frontiers in Ecology and Evolution | www.frontiersin.org 7October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 8 Dinis et al. Niche-Based Biocontrol Management of Acacia FIGURE 3 | Spatial overlap between model predictions of suitability for Acacia longifolia and Trichilogaster acaciaelongifoliae in the Mediterranean region. Trichilogaster acaciaelongifoliae model was constrained to areas of predicted suitability for A. longifolia to prevent classification of suitability in areas where the species cannot persist due to lack of its host plant. Bottom: Populations of T. acaciaelongifoliae already established in Portugal as of 2018 South Africa, as its dependency on A. longifolia would require the host plant to do the same, which is rejected by the evidence. Rather, non-climatic factors such as biotic interactions are likely to constrain the distribution of T. acaciaelongifoliae in Australia, where the species is native and as such suffers higher levels of parasitism and competes with other bud-feeding insects (Neser, 1984). Despite present, the lower level of such interactions in South Africa (Manongi and Hoffmann, 1995; Seymour and Veldtman, 2010) may have allowed the species to expand into previously unoccupied areas of niche space which are nonetheless contained within the realized niche of A. longifolia. The evidence for niche expansion may also reflect the low number of occurrence records for T. acaciaelongifoliae in Australia, which thus fail to fully represent the realized niche. Gall-forming insects, particularly in native regions where they coexist with their hosts often in equilibrium, are frequently overlooked and understudied organisms, despite the importance of several species as biocontrol agents (Dennill, 1988; Harris and Shorthouse, 1996) or pests (Gil-Tapetado et al., 2018; Ferracini et al., 2019). This resulted in a relatively reduced dataset in Australia, vulnerable to spatial biases, which are known to affect niche estimates (Graham et al., 2008;Sánchez-Fernández et al., 2011). Regardless, the evidence for niche conservatism allowed us to combine the T. acaciaelongifoliae occurrence records for Australia and South Africa for the development of niche models, minimizing the effect of spatial biases. The very high niche overlap between A. longifolia and the gall-former T. acaciaelongifoliae is consistent with the high specificity of T. acaciaelongifoliae with its host plant (Dennill et al., 1993;Marchante et al., 2011a). The realized niche of T. acaciaelongifoliae is found to be almost fully contained within the niche of A. longifolia on all regions where the two species co-occur (reflected in high niche stability and near zero niche expansion). In South Africa, niches of the two species were found to be more similar than expected given the available environmental background and niche unfilling was found to be the lowest among all interspecies comparisons. This suggests a near total colonization by T. acaciaelongifoliae Frontiers in Ecology and Evolution | www.frontiersin.org 8October 2020 | Volume 8 | Article 576667 fevo-08-576667 October 14, 2020 Time: 16:43 # 9 Dinis et al. Niche-Based Biocontrol Management of Acacia of all environments where A. longifolia occurs, which has been reported by some authors (J. Hoffmann, personal communication). It also indicates that the aforementioned expansion of T. acaciaelongifoliae in South Africa in the absence of substantial negative biotic interactions took place exclusively within the range of A. longifolia, highlighting the high specificity of T. acaciaelongifoliae which is a key factor in its reliability as a biocontrol agent. Australia had the lowest niche overlap and highest niche unfilling among all interspecies comparisons, highlighting that T. acaciaelongifoliae is not known to occur (or reported) in many of the environments where A. longifolia is present in its native range. Again, this likely reflects nonclimatic influences hindering the ability of T. acaciaelongifoliae to successfully colonize all A. longifolia habitats. Implications for the Management of A. longifolia Including Biocontrol The predicted suitable area for A. longifolia in the Mediterranean basin far exceeds the current distribution. Besides environmental conditions, introduction histories are known to greatly influence species distributions (Donaldson et al., 2014). One of the main pathways for A. longifolia introduction in Portugal was for sand binding of coastal dunes (Marchante et al., 2003). This coastal introduction, associated with habitat connectivity and dispersal constraints, may account for the current species distribution which is not known to extend as far inland as predicted by our model (with a few exceptions in some particular locations). This hypothesis seems likely when considering that A. longifolia (sensu lato) is known to occur much further inland in its native distribution than in Portugal (Court et al., 2020). Likewise, in South Africa, where A. longifolia was introduced also to stabilize inland slopes (Dennill and Donnelly, 1991), the species is present as far inland as 500 km, demonstrating that it can become established far from coastal climate given the adequate means of dispersal. Regardless, we cannot dismiss the possibility that factors not included in the models, such as biotic interactions or non-climatic abiotic variables, may constrain the distribution of A. longifolia beyond the model’s prediction. The current model predicts risk of A. longifolia invasion in most of the Mediterranean and, especially, Atlantic coastline of the Mediterranean Basin region. Areas of current A. longifolia invasion in the region besides Portugal include north-western Spain (Community of Galicia), which is correctly predicted as climatically suitable. A. longifolia is currently present, even if not yet considered as invasive, in several countries which are predicted as having suitable areas for the establishment of the species, namely France, Italy, Greece and Turkey (Akanil and Middleton, 2010;EFSA Panel on Plant Health, 2015). Our results suggest that particular care should be taken to prevent A. longifolia invasion in these regions where the species is already present within an area of suitable environment. The model for A. longifolia identified a small area of suitability in the southern Pyrenees Mountains, a region characterized by montane conditions that are uncharacteristic for the species. Species distribution models have a limited ability to extrapolate beyond the environmental conditions available in the training area and may produce unexpected results when environmental homology between training and projection areas is not assured (Elith and Leathwick, 2009). Despite our efforts to maximize homology by restricting the model projections to the Mediterranean Basin, mountain ranges in the region are likely to include extreme conditions which lead to extrapolation, resulting in less reliable predictions in these areas. Trichilogaster acaciaelongifoliae is predicted to find suitable climate in a significant part of the potential A. longifolia distribution area. Current areas of A. longifolia invasion where T. acaciaelongifoliae is not predicted to find suitable conditions are located in the coastline of northern Portugal and Galicia, as well in the most inland areas of Portugal. Nevertheless, the existence of a well-established population of T. acaciaelongifoliae in the Portuguese city of Coimbra and very recent (2019; López-Núñez et al., in prep.) establishment in northern areas of Portugal (classified as unsuitable by the model) indicates that the approach used in this work may not be sufficient to fully characterize the niche of this species. Partial inconsistencies amongst modeling outputs and field observations occur in other studies (e.g., Gallien et al., 2012;Fischbein et al., 2019) for several reasons. In our study, agreement between model predictions and known distributions in the Mediterranean Basin was higher for A. longifolia than for T. acaciaelongifoliae, likely reflecting the difference in number and quality of occurrence records between the two species. Moreover, this may be explained by the fact that the biocontrol agent was only recently (2015) introduced in Portugal (Marchante et al., 2017), still not having enough time to spread and establish across the suitable area. Furthermore, the life history of T. acaciaelongifoliae is characterized by a near-total existence as an endoparasite, and as such the microenvironmental conditions and interactions with the host plant are likely to be better predictors of the species’ ability to form viable populations than macroclimatic factors. Despite the identification of ecogeographical variables that have a limiting effect on T. acaciaelongifoliae, the environmental triggers required for completion of this species’ life cycle are still poorly understood. T. acaciaelongifoliae has shown to be sensitive to such triggers, as the change in hemisphere (and consequent temporal mismatch between the species’ life cycle and seasonal queues) upon introduction to Portugal resulted in low establishment success in the first years (Marchante et al., 2017, López-Núñez et al., in prep.). Identifying these triggers is fundamental to properly assess the long-term robustness of effective biocontrol of A. longifolia mediated by T. acaciaelongifoliae. Ultimately, mechanistic approaches to modeling which integrate microhabitat and biotic interactions between host and parasite as well as an explicit understanding of how these factors contribute to the life cycle of T. acaciaelongifoliae (Palhas et al., in prep.) are expected to better approximate the fundamental niche and, when integrated with correlative approaches such as the one presented in this work, allow for more accurate predictions. Despite the limitations described above, the high niche overlap between the two species shows promise regarding the effectiveness of T. acaciaelongifoliae for biocontrol of Frontiers in Ecology and Evolution | www.frontiersin.org 9October 2020 | Volume 8 | Article 576667