Ecology

Chai et al. Stable water isotope and surface heat flux simulation using ISOLSM: Evaluation against in-situ measurements. J. Hydrol. 523, 67–78 (2015).Article 

Google Scholar 
Brooks et al. Stable isotope estimates of evaporation: Inflow and water residence time for lakes across the United States as a tool for national lake water quality assessments. Limnol. Oceanogr. 59, 2150–2165 (2014).Article 

Google Scholar 
Good, S. P., Noone, D. & Bowen, G. Hydrologic connectivity constrains partitioning of global terrestrial water fluxes. Science 349, 175–177 (2015).CAS 
Article 

Google Scholar 
Gupta, A., Gerber, E. P. & Lauritzen, P. H. Numerical impacts on tracer transport: A proposed intercomparison test of atmospheric general circulation models. Quart. J. Roy. Meteor. Soc. 146, 3937–3964 (2020).Article 

Google Scholar 
Kanner, L. C., Buenning, N. H., Stott, L. D., Timmermann, A. & Noone, D. The role of soil processes in d18O. Global Biogeochem. Cycles 28, 239–252 (2014).CAS 
Article 

Google Scholar 
Remondi, F., Kircher, J. W., Burlando, P. & Fatichi, S. Water flux tracking with a distributed hydrologic model to quantify controls on the spatio-temporal variability of transit time distributions. Water Resour. Res. 54, 3081–3099 (2018).Article 

Google Scholar 
Abbott, B. W. et al. Using multi-tracer inference to move beyond single catchment ecohydrology. Earth-Sci. Rev. 160, 19–42 (2016).Article 

Google Scholar 
Krause, P., Boyle, D. P. & Bäse, F. Comparison of different efficiency criteria for hydrological model assessment. Adv. in Geosci. 5, 89–97 (2005).Article 

Google Scholar 
Bowen, G. J. & Good, S. P. Incorporating water isotopes in hydrological and water resource investigations. Wiley Interdiscip. Rev.: Water 2, 107–119 (2015).Article 

Google Scholar 
McGuire, K. J. & McDonnell, J. J. A review and evaluation of catchment transit time modeling. J. Hydrol. 330, 543–563 (2006).Article 

Google Scholar 
Sprenger, M. et al. The demographics of water: A review of water ages in the critical zone. Rev. Geophys. 57, 800–834 (2019).Article 

Google Scholar 
Turnadge, C. & Smerdon, B. D. A review of methods for modelling environmental tracers in groundwater: Advantages of tracer concentration simulation. J. Hydrol. 519, 3674–3689 (2014).CAS 
Article 

Google Scholar 
Fiorella, R. et al. Calibration Strategies for Detecting Macroscale Patterns in NEON Atmospheric Carbon Isotope Observations. J. Geophys. Res. Biogeosci. 126 (2021).Xiao, W., Wei, Z. & Wen, X. Evapotranspiration partitioning at the ecosystem scale using the stable isotope method—A review. Agric For Meteorol. 263, 346–361 (2018).Article 

Google Scholar 
Wu, Y. et al. Stable isotope measurements show increases in corn water use efficiency under deficit irrigation. Sci Rep 8, 14113 (2018).Article 

Google Scholar 
Al-Oqaili, F., Good, S. P., Frost, K. & Higgins, C. W. Differences in soil evaporation between row and interrow positions in furrowed agricultural fields. Vadose Zone J. 19, e20086 (2020).CAS 
Article 

Google Scholar 
Bowen, G. J., Cai, Z., Fiorella, R. P. & Putman, A. L. Isotopes in the water cycle: Regional- to global-scale patterns and applications. Annu. Rev. Earth Planet. Sci. 47, 453–479 (2019).CAS 
Article 

Google Scholar 
Lu, X. et al. Partitioning of evapotranspiration using a stable isotope technique in an arid and high temperature agricultural production system. Agric. Water Manag. 179, 103–109 (2017).Article 

Google Scholar 
Wieser, G. et al. Stable water use efficiency under climate change of three sympatric conifer species at the alpine treeline. Front. Plant Sci. 7, 799 (2016).Article 

Google Scholar 
Pataki, D. E. et al. The application and interpretation of Keeling plots in terrestrial carbon cycle research. Global Biogeochem. Cycles, 17 (2003).Miller, J. B., & Tans, P. P., Calculating isotopic fractionation from atmospheric measurements at various scales. Tellus, 55 (2003).Finkenbiner, C. E., Good, S. P., Allen, S. T., Fiorella, R. P. & Bowen, G. J. A statistical method for generating temporally downscaled geochemical tracers in precipitation. J. Hydrometeorol. 22 (2021).NEON (National Ecological Observatory Network). Precipitation (DP1.00006.001), RELEASE-2022. https://doi.org/10.48443/6wkc-1p05. Dataset accessed from https://data.neonscience.org on May 12, 2022.Lunch, C. K. & Laney, C. M. NEON (National Ecological Observatory Network). neonUtilities: Utilities for working with NEON data. R package version 1.3.4. https://github.com/NEONScience/NEON-utilities (2020).Lee, R. and S. Weintraub. NEON User Guide to Stable Isotopes in Precipitation (NEON.DPI.00038) Version B. NEON (National Ecological Observatory Network). (2021).IAEA: Global network of isotopes in precipitation. https://www.iaea.org/services/networks/gnip 2020.Allen, S. T., Kirchner, J. W. & Goldsmith, G. R. Predicting spatial patterns in precipitation isotope (δ2H and δ18O) seasonality using sinusoidal isoscapes. Geophys. Res. 45, 4859–4868 (2018).
Google Scholar 
Craig, H. Isotopic variations in meteoric waters. Science 133, 1702–1703 (1961).CAS 
Article 

Google Scholar 
Dansgaard, W. Stable isotopes in precipitation. Tellus 16, 436–468 (1964).Article 

Google Scholar 
Sklar, A. Fonctions de répartition à n dimensions et leurs marges. Publ. Inst. Stat. Univ. Paris. 8, 229–231 (1959).MATH 

Google Scholar 
NEON (National Ecological Observatory Network). Bundled data products – eddy covariance (DP4.00200.001). https://data.neonscience.org (2021).Good, S. P., Soderberg, K., Wang, L., & Caylor, K. K. Uncertainties in the assessment of the isotopic composition of surface fluxes: A direct comparison of techniques using laser‐based water vapor isotope analyzers. J. Geophys. Res. Atmos. 177 (2012).Wutzler, T. et al. Basic and extensible post-processing of eddy covariance flux data with REddyProc. Biogeosci. 15, 5015–5030 (2018).CAS 
Article 

Google Scholar 
Zobitz, J. M., Keener, J. P., Schnyder, H. & Bowling, D. R. Sensitivity analysis and quantification of uncertainty for isotopic mixing relationships in carbon cycle research. Agric For Meteorol. 136, 56–75 (2006).Article 

Google Scholar 
Wehr, R. & Saleska, S. R. An improved isotopic method for partitioning net ecosystem-atmosphere CO2 exchange. Agric For Meteorol. 214, 515–531 (2015).Article 

Google Scholar 
Bailey, A., Noone, D., Berkelhammer, M., Steen-Larsen, H. C. & Sato, P. The stability and calibration of water vapor isotope ratio measurements during long-term deployments. Atmos. Meas. Tech. 8, 4521–4538 (2015).CAS 
Article 

Google Scholar 
Rambo, J., Lai, C., Farlin, J., Schroeder, M. & Bible, K. Vapor isotope ratios using off-axis cavity-enhanced absorption spectroscopy. J Atmos. Ocean Technol. 28, 1448–1457 (2011).Article 

Google Scholar 
Finkenbiner, C. The National Ecological Observation Network Daily Isotopic Composition of Environmental Exchanges (NEON-DICEE) Dataset, HydroShare, https://doi.org/10.4211/hs.e74edc35d45441579d51286ea01b519f (2022). More

InvasionWe used freshwater fish biodiversity data collated by and described in Milardi, et al.47. In summary, the dataset included 3777 sites sampled 1999–2014, recorded a total of 99 different fish species (35 of which were exotic and already established, even if some are restricted to areas with thermal springs), spanned  > 11 degrees of longitude (~ 1200 km) and 10 degrees of latitude (~ 1100 km), covering streams at altitudes -2.7–2500 m above sea level. Community turnover was not a relevant factor in our study, because fish communities are typically stable over these timescales and the data was collected in a restricted timeframe within each area29,39. Furthermore, time elapsed since last introductions was sufficient to analyze distribution patterns after major invasions had already occurred see e.g.23,48.Abundance of each species sampled during the monitoring was recorded with Moyle classes (Moyle and Nichols, 1973), which were weighted according to body-size classes in order to obtain a body-mass-corrected abundance, hereafter referred to simply as abundance. We then calculated an invasion degree, i.e. the share of introduced species in freshwater fish communities, based on the abundance of introduced and native species see e.g.9,49. A high invasion degree equals to a high share of introduced species and a low share of native species.We also selected the top 10 invasive species as further response variables, under the assumption that these would be the main components of the invasion degree, but would respond to different invasion drivers based on each species’ ecology. Invasiveness rank was defined through an index obtained by multiplying colonization (share of sites colonized) and prevalence (average relative abundance in the fish community) of each introduced species. The relative abundance of each of these species in the fish community was used as a response variable, being a measure comparable to invasion degree for single species.Invasion driversWe tested a combination of geographical, climate and anthropogenic impact factors as potential drivers of invasion. To avoid temporal mismatches, we chose time periods that overlapped as much as possible with our biological data.We used basin area, altitude and slope (derived from a seamless digital elevation model of the whole Italian territory at 10 m resolution, Tarquini, et al.50) as geographical variables.We derived climate data from available series of long-term national monitoring (http://www.scia.isprambiente.it/). We used daily air temperature (2000–2009), measured at a total of 2266 sites throughout the country, as a proxy for temperature regimes. We also used cumulated annual precipitations, number of annual dry days (precipitation  More

In this section, we provide quantitative evaluation for both spatial coverage and temporal dynamics of ReaLSAT dataset.Spatial coverageSince the dataset was created using satellite imagery analysis, it can provide more comprehensive coverage than existing datasets. However, using an automated process also has its challenges. It can invariably lead to the detection of spurious waterbodies because of issues in data (e.g., due to errors in GSW maps used as inputs in ReaLSAT).To provide more insights into the types of lakes and potential issues in the spatial coverage of ReaLSAT, we randomly sampled 5,000 lakes out of 435,717 that are only present in ReaLSAT (i.e., not available in the HydroLAKES dataset). A human annotator used Google’s satellite imagery base layer to categorize these lakes. Figure 5a shows the geographical distribution of these lakes, and Fig. 5b shows the distribution of different lake types in the sample set. Out of the 5,000 lakes, the human annotator identified 2,019 traditional lakes and reservoirs where sufficient water was visible in the satellite imagery. Another 551 lakes in the sample set showed signs of a bowl-like depression but with no (or very little) water visible in the satellite imagery and were labeled as ephemeral. There were 861 other lakes that were tagged as farm ponds because they showed geometric patterns of farming in the imagery. This diversity of waterbody types discovered by ReaLSAT that were previously unreported by HydroLakes highlights one of the strengths of our approach. In limnology, the origin/type of lake is a very important regulator of ecosystem dynamics. For instance, reservoirs will have faster water flow/lower residence time than natural lakes, and therefore nutrient and carbon processing rates will differ; floodplain lakes may dry periodically, leading to the denudation of sediments; and farm ponds will likely have much higher rates of nutrient loading and methane production than non-agriculturally influenced lakes. Hence, capturing a more comprehensive range of waterbody categories can enable various scientific studies where knowing the origin/lake type could provide a critical understanding of the process.Fig. 5(a) Geographic location of 5000 randomly selected lakes used for manual evaluation of lake type. (b) Allocation of the 5000 manually referenced lakes to specific lake types. Regular implies a traditional lake or reservoir. Unverifiable implies that the lake type could not be identified based on the available Google Earth imagery.Full size imageAlong with the lentic water types discovered in the sampled set, we also found that ReaLSAT identified 603 river segments missed by our morphological score filter. As stated earlier, this is an inherent challenge with automated approaches that use a fixed score threshold for eliminating river segments. Another 239 lakes were tagged as wetlands because of significant vegetation inside and around the lake polygon. There were also 97 lakes that were adjacent to rivers, which were labeled as riverine or floodplain lakes that were formed as a result of river channels meandering over time. Furthermore, there were 59 lakes where the polygons represented only a small portion of a larger lake and were labeled as partial. Finally, for 571 polygons, there was not enough evidence to tag them in any of the above categories. Since Google imagery represents only a single snapshot in time, these 571 waterbodies could not be definitively labeled as spurious (hence, they were labeled as unverifiable), highlighting a limitation of this evaluation pipeline. In particular, a vast majority of these waterbodies appear to be ephemeral based on their surface area timeseries (completely dry for extended periods of time). Hence, if the satellite imagery layer is from one of these timesteps, the annotator would not be able to confirm the presence of the lake.To assess whether we would obtain a similar distribution of different waterbody categories in existing datasets, we performed a similar evaluation on another 5,000 lakes sampled from ReaLSAT where each polygon has some overlap (greater than 1 pixel) with a polygon from HydroLAKES. In this sampled set, the annotator identified 4,030 lakes as traditional lakes or reservoirs, 370 as ephemeral, 138 as farm ponds, 6 as river segments, 66 as wetlands, 95 as riverine or floodplain lakes, 20 as partial, and 275 as unverifiable.Compared to previous distribution, this set of 5,000 waterbodies contains relatively fewer river segments and wetlands polygons in HydroLAKES, because these categories were manually identified and removed during HydroLAKES database creation6. Similary, this set contains relatively few farm ponds because HydroLAKES was created by manual curation of existing static databases and hence does not contain new farm ponds that got created over the years.Temporal dynamicsTo assess the quality of surface extent maps, we performed a quantitative evaluation on a random selection of extent maps. These extent maps were compared against reference maps created by a human annotator using a semi-automated pixel classification procedure. This strategy of creating reference maps is used extensively in the remote sensing literature (e.g. see36,37,38,39). Next, we describe our evaluation process in detail.Sample selectionThere are 462,574 lakes out of 681,137 total lakes where the label updates (corrections and imputations) by the ORBIT approach have trust scores within our chosen thresholds (as described in the methods section). To evaluate these candidate lakes effectively, we focus on lake extent maps where the ORBIT approach resulted in a different map than the underlying GSW extent based map. Hence, we remove maps where no updates were made by the ORBIT approach (neither corrections nor imputations) from the candidate pool of extent maps used for evaluation. We also remove maps where the percentage of missing labels was more than 90% because these maps tend to suffer from significant cloud coverage. Hence, it would be challenging to generate reference maps. Since the GSW dataset has a significant amount of missing data for most places in the world before 2000, we evaluated maps only from 2000 onwards. These three filters left us with a total of 51,077,278 water extent maps considered for selection. Figure 6a shows the distribution of percentage pixels updated made by the ORBIT approach in these water extent maps. To evaluate the robustness of our approach in comparison to GSW maps, we randomly selected 10,000 water extent maps such that extents with significant updates are given higher weight to reduce the skew in distribution towards extents with relative less updates (Fig. 6b).Fig. 6Distribution of updates made by the ORBIT approach. (a) distribution using candidate water extents (b) distribution using randomly selected 10000 water extent maps for evaluation.Full size imageSample pruningFrom these randomly selected water extent maps, we removed maps for which a reference map could not be generated due to clouds or the inability of the annotator to distinguish between land and water. A final set of 2,095 water extent maps were considered for evaluation. Figure 7a shows the distribution of percentage updates in the final set of evaluation extents and Fig. 7b shows the geographical distribution of these extent maps.Fig. 7Summary of the dataset used for evaluating water extent maps. (a) Distribution of updates made by the ORBIT approach in the water extent maps selected for evaluation. (b) Geographical location of the lakes in the evaluation set.Full size imageReference map generationFor these water extent maps, we created ground truth reference maps using a semi-automatic labeling process37,38,39. Specifically, the annotator selects land and water samples to train an SVM (Support Vector Machine) classification model for each image. The annotator keeps adding samples until a stable map is generated. As a final step, the annotator masks out pixels affected by clouds, cloud shadows, and any other region where the annotator is not confident about the accuracy of the reference labels. This process enables a quick and robust generation of reference maps. Supplementary Fig. S7 shows one of the reference maps in the evaluation set. While this strategy of comparing maps is different from the traditional approach of comparing pixels (often selected using stratified sampling), it provides a much more exhaustive evaluation of surface extent maps. The reference maps used for evaluation in this study are also available as part of the dataset.ComparisonTo compare the extent maps generated by ReaLSAT with the reference maps, we used accuracy as the evaluation metric, a widely used metric to measure the quality of classification maps. Accuracy is simply defined as the ratio of pixels with correct labels over a total number of pixels. Specifically, we assign 1 to water pixels and 0 to land pixels. Since GSW based extent maps contain missing labels, they are assigned a value of 0.5 to reflect the uncertainty between land and water. Accuracy is then calculated as follows:$$Accuracy=1-frac{1}{Rast C}mathop{sum }limits_{i=1}^{R}mathop{sum }limits_{j=1}^{C}left|ReferenceMap[i,j]-PredictedMap[i,j]right|$$
(2)
where, R is the number of rows and C is the number of columns of the map.When the accuracy of RealSAT and GSW labels are compared, a vast majority of points lie above the diagonal 1:1 line, which implies that ReaLSAT labels were more accurate overall (Fig. 8a). In Fig. 8 the points are colored based on % of pixels where GSW labels were missing. To better show the improvement in RealSAT labeling, we plot the distribution of the difference in accuracy values between the two datasets as shown in Fig. 8b. A positive value indicates that the surface extent map from the ReaLSAT dataset had better accuracy than the map from the GSW dataset and vice versa. For ease of visualization, we plot this distribution after excluding cases where the accuracy from both datasets was equal. The positively skewed distribution demonstrates the efficacy of the ORBIT approach.Fig. 8Comparison of accuracy values using GSW labels vs ReaLSAT labels. (a) Scatter plot of accuracy values using GSW labels vs ReaLSAT labels. (b) Histogram of difference in accuracy between ReaLSAT labels vs GSW labels. Positive value represents cases where ReaLSAT labels were more accurate than GSW labels. (c) Histogram of difference in accuracy values for the scenario where pixels labelled as land by both products as well as ground truth were removed to reduce the skew of surrounding land pixel on the accuracy values.Full size imageNote that the shape of a lake will influence the number of land pixels surrounding it, which might bias the accuracy values. For example, the reference map shown in Supplementary Fig. S7 contains more than 70% of land pixels. To address this bias, we also calculated accuracy values after removing pixels that were labeled as land by both datasets as well as the ground truth. This variation allows a more strict evaluation of water extent maps. Figure 8c shows the distribution of the difference in accuracy values under this scenario (after excluding cases with equal accuracy). As shown, a vast majority of the distribution is still towards positive values. Furthermore, the distribution has a larger spread towards high positive values, suggesting significant improvement made by the ORBIT approach.From Fig. 8, we can see that for some cases ReaLSAT based extent maps are less accurate relative to GSW. As described earlier, violation of assumptions made by the ORBIT approach could lead to the observed poor performance. Out of 2,095 extent maps, GSW labels show better accuracy than ReaLSAT for 323 of them. On visual analysis of errors in these maps, we found that 165 maps are slightly different only at the lake’s boundary. We categorized the remaining extent maps based on the reason behind the observed poor performance. In particular, 45 maps have poor performance due to occlusion of water surface by algae, 18 maps contain farm ponds, 8 contain mining lakes, 27 maps have unreliable bathymetry, 30 maps have issues due to the weighting factor used by ORBIT approach, and 30 maps have class conditional missing data. All the reference maps and corresponding maps from GSW and ReaLSAT are provided with the dataset.Next, we describe some of these cases in detail.Impact of algae: It can be difficult to visually differentiate surface algae or floating aquatic plants from terrestrial vegetation40, as they have similar reflectance spectra. Therefore, surface algal blooms often get incorrectly labeled as land in the reference maps. However, in most cases, the appearance and disappearance of algae on a lake are independent of the bathymetry. Thus, algae pixels get detected as physically inconsistent by the ORBIT approach, and consequently, these pixels are updated based on the labels of other pixels without algae. In many cases, while the accuracy with respect to the reference map is poor (because algae get labeled as land), ReaLSAT based extent maps are closer to the true extent of the lake. For example, Supplementary Fig. S8 illustrates the impact of algae on the extent mapping of Center Lake, Texas. In this example, the bimodal distribution of fraction values (either low or high) reveals high confidence in lake persistence (Supplementary Fig. S8b). On Oct 22, 2008, false-color composite processing of LANDSAT-5 imagery reveals a strong vegetative signal on the west side of the lake (Supplementary Fig. S8c). Since we know that this is a lake, we can assume that the west side of the lake is experiencing a large surface algal bloom with a similar reflectance to the surrounding terrestrial landscape. Because of the strong vegetative reflectance signal, the semi-automated reference mapping labels the west side of the lake as land (Supplementary Fig. S8d), as does most GSW labels (Supplementary Fig. S8e). Conversely, the ReaLSAT extent map labels the west side of the lake as water (Supplementary Fig. S8f). However, we calculate accuracy based on the semi-automated reference map (Supplementary Fig. S8d). Due to this, the GSW extent map is considered more accurate than the ReaLSAT map, even though this is not true because the reference map is incorrectly labeled. Therefore, some negative accuracy values may be a misrepresentation of reality due to surface algal blooms.Impact of variable bathymetry: Even though we tried to remove lakes with unreliable bathymetry by using score-based filters defined in an earlier section, not all cases were removed. For example, agricultural ponds often have small sections that are connected and change shape based on agricultural needs. Supplementary Fig. S9 highlights an example of labeling issues on agricultural ponds in Mexico. In this area, satellite imagery and the GSW fraction map confirm the presence of agricultural ponds (Supplementary Fig. S9a,b). These individual ponds are filled and drained based on operational decisions and do not follow a consistent pattern of growing or shrinking. Thus, the ORBIT approach can introduce spurious updates in water extent maps for these farms. In the Landsat-5 imagery from 2009–10–08, some of the ponds are dry, while others are filled (Supplementary Fig. S9c). This distribution of water is evident from a visual inspection and is confirmed in the semi-automated reference map (Supplementary Fig. S9d). Due to the similar elevations between the individual pond sections, the ORBIT approach spuriously fills the remaining sections with water based on the incorrectly learned bathymetry (Supplementary Fig. S9f). While quantification of such uncertainties is outside the scope of this paper, we hope that the wider research community can use RealSAT to address such questions. In particular, changes in bathymetry of a lake can be identified using spatial-temporal patterns in the label corrections. Specifically, if the elevation of some pixels in a lake increases after a certain time (e.g., sediment deposits leading to increase the elevation of a pixel), they will appear as physically inconsistent to the ORBIT framework, and hence the labels for these locations will be changed from land to water much more frequently after this increase in elevation.Impact of bias in errors and missing data: As mentioned earlier in the methods section, based on our observation, the confidence of water labels is higher than land labels in the GSW dataset. To account for this bias, we used a weighting factor of 3 for the water class. While this weighting factor improves the ORBIT approach’s performance in most cases, this assumption leads to an overestimation of water for some lakes. For example, Supplementary Fig. S10 compares the water extent maps with and without the weighting factor for a small reservoir in eastern Brazil. As we can see, the GSW labels contain false positives, and due to the weighting factor of 3, ORBIT prefers to update the land labels to water which further increases the number of false positives, as shown in Supplementary Fig. S10e. However, if we use a weighting factor of 1 for this example, the ORBIT approach can effectively remove many of the false positives in the GSW map, as shown in Supplementary Fig. S10f.Similarly, apart from missing data due to clouds in the GSW dataset, there can also be missing values on pixels where the GSW classification model is not confident. Hence, for some water extent maps, class-dependent missing data (compared to missing data which is class independent) adversely impact the ORBIT approach. For example, Supplementary Fig. S11 shows a water extent map for Zhongleng Reservoir in China, where missing data along the eastern edges is not independent but has resulted from ambiguous pixels around the lake where the GSW’s approach was not confident. In such a scenario, the ORBIT approach heavily relies on information from nearby timesteps to infer labels for missing pixels, leading to errors in ReaLSAT maps if there is a significant variation in lake extent in nearby timesteps, as shown in Supplementary Fig. S11e. More

Aissi A, Beghami Y, Heuertz M (2019) Le chêne faginé (Quercus faginea, Fagaceae) en Algérie: potentiel germinatif et variabilité morphologique des glands et des semis. Plant Ecol Evol 152:437–449Article 

Google Scholar 
Aissi A, Beghami Y, Lepais O, Véla E (2021) Morphological and taxonomic analysis of Quercus faginea (Fagaceae) complex in Algeria. Botany 99:99–113Article 

Google Scholar 
Alberto F, Niort J, Derory J, Lepais O, Vitalis R, Galop D et al. (2010) Population differentiation of sessile oak at the altitudinal front of migration in the French Pyrenees. Mol Ecol 19:2626–2639CAS 
PubMed 
Article 

Google Scholar 
Allendorf FW, Luikart G (2007) Conservation and the genetics of populations. Blackwell PubDe Barba M, Miquel C, Lobréaux S, Quenette PY, Swenson JE, Taberlet P (2016) High-throughput microsatellite genotyping in ecology: Improved accuracy, efficiency, standardization and success with low-quantity and degraded DNA. Mol Ecol Resour 17:492–507PubMed 
Article 
CAS 

Google Scholar 
Barthe S, Gugerli F, Barkley NA, Maggia L, Cardi C, Scotti I (2012) Always look on both sides: phylogenetic information conveyed by simple sequence repeat allele sequences. PLoS One 7:e40699CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Beaumont MA, Zhang W, Balding DJ (2002) Approximate Bayesian computation in population genetics. Genetics 162:2025–35PubMed 
PubMed Central 
Article 

Google Scholar 
Belkhir K, Borsa P, Chikhi L, Raufaste N, Bonhomme F (2004) GENETIX 4.05, logiciel sous Windows TM pour la génétique des populations. Laboratoire Génome, Populations, Interactions, CNRS UMR 5171, Université de Montpellier II, Montpellier (France)
Google Scholar 
Bradbury IR, Wringe BF, Watson B, Paterson I, Horne J, Beiko R et al. (2018) Genotyping-by-sequencing of genome-wide microsatellite loci reveals fine-scale harvest composition in a coastal Atlantic salmon fishery. Evol Appl 11:918–930CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Buschbom J, Yanbaev Y, Degen B (2011) Efficient long-distance gene flow into an isolated relict oak stand. J Hered 102:464–472PubMed 
Article 

Google Scholar 
Chapuis M, Raynal L, Plantamp C, Meynard CN, Blondin L, Marin J et al. (2020) A young age of subspecific divergence in the desert locust inferred by ABC Random Forest. Mol Ecol 29:4542–4558PubMed 
Article 

Google Scholar 
Cornuet J-M, Ravigné V, Estoup A (2010) Inference on population history and model checking using DNA sequence and microsatellite data with the software DIYABC (v1.0). BMC Bioinforma 11:401Article 
CAS 

Google Scholar 
Crow JF, Aoki K (1984) Group selection for a polygenic behavioral trait: Estimating the degree of population subdivision. Proc Natl Acad Sci USA 81:6073–6077CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Curto M, Winter S, Seiter A, Schmid L, Scheicher K, Barthel LMF et al. (2019) Application of a SSR-GBS marker system on investigation of European Hedgehog species and their hybrid zone dynamics. Ecol Evol 9:2814–2832PubMed 
PubMed Central 
Article 

Google Scholar 
Darby BJ, Erickson SF, Hervey SD, Ellis-Felege SN (2016) Digital fragment analysis of short tandem repeats by high-throughput amplicon sequencing. Ecol Evol 6:4502–4512PubMed 
PubMed Central 
Article 

Google Scholar 
Dickey AM, Hall PM, Shatters RG, Mckenzie CL (2013) Evolution and homoplasy at the Bem6 microsatellite locus in three sweetpotato whitefly (Bemisia tabaci) cryptic species. BMC Res Notes 6:249PubMed 
PubMed Central 
Article 

Google Scholar 
Do C, Waples RS, Peel D, Macbeth GM, Tillett BJ, Ovenden JR (2013) NeEstimator v2: re-implementation of software for the estimation of contemporary effective population size (Ne) from genetic data. Mol Ecol Resour 14:209–214PubMed 
Article 

Google Scholar 
Durand J, Bodenes C, Chancerel E, Frigerio JM, Vendramin G, Sebastiani F et al. (2010) A fast and cost-effective approach to develop and map EST-SSR markers: oak as a case study. BMC Genomics 11:570PubMed 
PubMed Central 
Article 

Google Scholar 
Earl DA, vonHoldt BM (2012) STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method. Conserv Genet Resour 4:359–361Article 

Google Scholar 
Estoup A, Jarne P, Cornuet J-M (2002) Homoplasy and mutation model at microsatellite loci and their consequences for population genetics analysis. Mol Ecol 11:1591–1604CAS 
PubMed 
Article 

Google Scholar 
Estoup A, Raynal L, Verdu P, Marin J-M (2018) Model choice using Approximate Bayesian Computation and Random Forests: analyses based on model grouping to make inferences about the genetic history of Pygmy human populations. J la Société Fr Stat 159:167–190
Google Scholar 
Excoffier L, Dupanloup I, Huerta-Sánchez E, Sousa VC, Foll M (2013) Robust demographic inference from genomic and SNP data. PLoS Genet 9:e1003905PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
Excoffier L, Lischer HEL (2010) Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol Ecol Resour 10:564–567PubMed 
Article 

Google Scholar 
Falush D, Stephens M, Pritchard JK (2003) Inference of population structure using multilocus genotype data: linked loci and correlated allele frequencies. Genetics 164:1567–87CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Feliner GN (2014) Patterns and processes in plant phylogeography in the Mediterranean Basin. A review. Perspect Plant Ecol, Evol Syst 16:265–278Article 

Google Scholar 
Flagel L, Brandvain Y, Schrider DR (2019) The unreasonable effectiveness of convolutional neural networks in population genetic inference. Mol Biol Evol 36:220–238CAS 
PubMed 
Article 

Google Scholar 
Foll M, Gaggiotti O (2008) A genome-scan method to identify selected loci appropriate for both dominant and codominant markers: A Bayesian perspective. Genetics 180:977–993PubMed 
PubMed Central 
Article 

Google Scholar 
Gaggiotti OE, Chao A, Peres-Neto P, Chiu CH, Edwards C, Fortin MJ et al. (2018) Diversity from genes to ecosystems: A unifying framework to study variation across biological metrics and scales. Evol Appl 11:1176–1193PubMed 
PubMed Central 
Article 

Google Scholar 
García Murillo P., Harvey-Brown Y (2017) Quercus canariensis. In: The IUCN Red List of Threatened Species,, p e.T78809256A80570536GBIF Secratariat (2021a) Quercus faginea Lam. In: GBIF Backbone Taxonomy. Checklist dataset https://doi.org/10.15468/39omeiaccessed via GBIF.org on 2022-05-10GBIF Secratariat (2021b). Quercus canariensis Willd. In: GBIF Backbone Taxonomy. Checklist dataset https://doi.org/10.15468/39omeiaccessed via GBIF.org on 2022-05-10Gómez A, Lunt DH (2007) Refugia within refugia: Patterns of phylogeographic concordance in the Iberian peninsula. In: Weiss S, Ferrand N (eds) Phylogeography of Southern European Refugia. Springer Netherlands, Dordrecht, p 155–188Chapter 

Google Scholar 
Gorener V, Harvey-Brown Y, Barstow M (2017) Quercus canariensis. IUCN red List Threat species e.T7880925Goudet J (1995) FSTAT (Version 1.2): A computer program to calculate F-statistics. J Hered 86:485–486Article 

Google Scholar 
Hampe A, Petit RJ (2005) Conserving biodiversity under climate change: the rear edge matters. Ecol Lett 8:461–7PubMed 
Article 

Google Scholar 
Hardy OJ, Charbonnel N, Fréville H, Heuertz M (2003) Microsatellite allele sizes: a simple test to assess their significance on genetic differentiation. Genetics 163:1467–82CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Hardy OJ, Vekemans X (2002) SPAGeDi: a versatile computer program to analyse spatial genetic structure at the individual or population levels. Mol Ecol Notes 2:618–620Article 
CAS 

Google Scholar 
Harvey-Brown Y, García Murillo PG, Buira A (2017) Quercus faginea. IUCN Red List Threat Species: e.T78916251A80570540.Henriques R, von der Heyden S, Matthee CA (2016) When homoplasy mimics hybridization: a case study of Cape hakes (Merluccius capensis and M. paradoxus). PeerJ 4:e1827PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
Hewitt GM (1999) Post-glacial re-colonization of European biota. Biol J Linn Soc 68:87–112Article 

Google Scholar 
Hey J, Won YJ, Sivasundar A, Nielsen R, Markert JA (2004) Using nuclear haplotypes with microsatellites to study gene flow between recently separated Cichlid species. Mol Ecol 13:909–919CAS 
PubMed 
Article 

Google Scholar 
Hoban S, Bruford M, D’Urban Jackson J, Lopes-Fernandes M, Heuertz M, Hohenlohe PA et al. (2020) Genetic diversity targets and indicators in the CBD post-2020 Global Biodiversity Framework must be improved. Biol Conserv 248:108654Article 

Google Scholar 
Hoogenboom J, de Knijff P, Laros JFJ, de Leeuw RH, van der Gaag KJ, Sijen T (2016) FDSTools: A software package for analysis of massively parallel sequencing data with the ability to recognise and correct STR stutter and other PCR or sequencing noise. Forensic Sci Int Genet 27:27–40PubMed 
Article 
CAS 

Google Scholar 
Hothorn T, Hornik K, van de Wiel MA, Zeileis A (2008) Implementing a class of permutation tests: the coin package. J Stat Softw 28:1–23Article 

Google Scholar 
Jamieson IG, Allendorf FW (2012) How does the 50/500 rule apply to MVPs? Trends Ecol Evol 27:578–584PubMed 
Article 

Google Scholar 
Jerome D, Vasquez F (2018) Quercus faginea. IUCN Red List Threat Species e.T7891625Kalinowski ST (2005) HP-RARE 1.0: A computer program for performing rarefaction on measures of allelic richness. Mol Ecol Notes 5:187–189CAS 
Article 

Google Scholar 
Kampfer S, Lexer C, Glössl J, Steinkellner H (1998) Characterization of (GA)n microsatellite loci from Quercus robur. Hereditas 129:183–186CAS 
Article 

Google Scholar 
Kivelä M, Arnaud-Haond S, Saramäki J (2015) EDENetworks: A user-friendly software to build and analyse networks in biogeography, ecology and population genetics. Mol Ecol Resour 15:117–122PubMed 
Article 

Google Scholar 
Kopelman NM, Mayzel J, Jakobsson M, Rosenberg NA, Mayrose I (2015) Clumpak: A program for identifying clustering modes and packaging population structure inferences across K. Mol Ecol Resour 15:1179–1191CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Layton KKS, Dempson B, Snelgrove PVR, Duffy SJ, Messmer AM, Paterson IG et al. (2020) Resolving fine‐scale population structure and fishery exploitation using sequenced microsatellites in a northern fish. Evol Appl: eva.12922.Lepais O, Chancerel E, Boury C, Salin F, Manicki A, Taillebois L et al. (2020) Fast sequence-based microsatellite genotyping development workflow. PeerJ 8:e9085PubMed 
PubMed Central 
Article 

Google Scholar 
Lepais O, Leger V, Gerber S (2006) Short note: high throughput microsatellite genotyping in oak species. Silvae Genet 55:238Article 

Google Scholar 
Lepais O, Muller SD, Ben Saad-Limam S, Benslama M, Rhazi L, Belouahem-Abed D et al. (2013) High genetic diversity and distinctiveness of rear-edge climate relicts maintained by ancient tetraploidisation for Alnus glutinosa. PLoS One 8:e75029CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Leroy T, Roux C, Villate L, Bodénès C, Romiguier J, Paiva JAP et al. (2017) Extensive recent secondary contacts between four European white oak species. N. Phytol 214:865–878CAS 
Article 

Google Scholar 
Lye GC, Lepais O, Goulson D (2011) Reconstructing demographic events from population genetic data: the introduction of bumblebees to New Zealand. Mol Ecol 20:2888–900CAS 
PubMed 
Article 

Google Scholar 
Magri D, Fineschi S, Bellarosa R, Buonamici A, Sebastiani F, Schirone B et al. (2007) The distribution of Quercus suber chloroplast haplotypes matches the palaeogeographical history of the western Mediterranean. Mol Ecol 16:5259–66CAS 
PubMed 
Article 

Google Scholar 
Marin J, Pudlo P, Estoup A, Robert C (2018) Likelihood-free model choice. In: Sisson S A, Fan Y, Beaumont M (eds) Handbook of Approximate Bayesian Computation, CRC Press, pp. 153.Médail F, Diadema K (2009) Glacial refugia influence plant diversity patterns in the Mediterranean Basin. J Biogeogr 36:1333–1345Article 

Google Scholar 
Moracho E, Moreno G, Jordano P, Hampe A (2016) Unusually limited pollen dispersal and connectivity of Pedunculate oak (Quercus robur) refugial populations at the species’ southern range margin. Mol Ecol 25:3319–3331CAS 
PubMed 
Article 

Google Scholar 
Mountain JL, Knight A, Jobin M, Gignoux C, Miller A, Lin AA et al. (2002) SNPSTRs: Empirically derived, rapidly typed, autosomal haplotypes for inference of population history and mutational processes. Genome Res 12:1766–1772CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Muir G, Lowe AJ, Fleming CC, Vogl C (2004) High nuclear genetic diversity, high levels of outcrossing and low differentiation among remnant populations of Quercus petraea at the margin of its range in Ireland. Ann Bot 93:691–697CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Neophytou C, Gärtner SM, Vargas-Gaete R, Michiels H-G (2015) Genetic variation of Central European oaks: shaped by evolutionary factors and human intervention? Tree Genet Genomes 11:1–15Article 

Google Scholar 
Payseur BA, Cutter AD (2006) Integrating patterns of polymorphism at SNPs and STRs. Trends Genet 22:424–429CAS 
PubMed 
Article 

Google Scholar 
Petit RJ, Brewer S, Bordács S, Burg K, Cheddadi R, Coart E et al. (2002) Identification of refugia and post-glacial colonisation routes of European white oaks based on chloroplast DNA and fossil pollen evidence. Ecol Manag 156:49–74Article 

Google Scholar 
Petit RJ, Hampe A, Cheddadi R (2005) Climate changes and tree phylogeography in the Mediterranean. Taxon 54:877–885Article 

Google Scholar 
Press MO, Hall AN, Morton EA, Queitsch C (2019) Substitutions are boring: Some arguments about parallel mutations and high mutation rates. Trends Genet 35:253–264CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Press MO, Mccoy RC, Hall AN, Akey JM, Queitsch C (2018) Massive variation of short tandem repeats with functional consequences across strains of Arabidopsis thaliana. Genome Res 28:1169–1178CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Pritchard JK, Stephens M, Donnelly P (2000) Inference of population structure using multilocus genotype data. Genetics 155:945–59CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Pudlo P, Marin J-M, Estoup A, Cornuet J-M, Gautier M, Robert CP (2016) Reliable ABC model choice via random forests. Bioinformatics 32:859–866CAS 
PubMed 
Article 

Google Scholar 
Ramakrishnan U, Mountain JL (2004) Precision and accuracy of divergence time estimates from STR and SNPSTR variation. Mol Biol Evol 21:1960–1971CAS 
PubMed 
Article 

Google Scholar 
Raynal L, Marin J-M, Pudlo P, Ribatet M, Robert CP, Estoup A (2019) ABC random forests for Bayesian parameter inference. Bioinformatics 35:1720–1728CAS 
PubMed 
Article 

Google Scholar 
Rodríguez-Sánchez F, Hampe A, Jordano P, Arroyo J (2010) Past tree range dynamics in the Iberian Peninsula inferred through phylogeography and palaeodistribution modelling: A review. Rev Palaeobot Palynol 162:507–521Article 

Google Scholar 
Šarhanová P, Pfanzelt S, Brandt R, Himmelbach A, Blattner FR (2018) SSR-seq: Genotyping of microsatellites using next-generation sequencing reveals higher level of polymorphism as compared to traditional fragment size scoring. Ecol Evol 8:10817–10833PubMed 
PubMed Central 
Article 

Google Scholar 
Scotti-Saintagne C, Mariette S, Porth I, Goicoechea PG, Barreneche T, Bodénès C et al. (2004) Genome scanning for interspecific differentiation between two closely related oak species [Quercus robur L. and Q. petraea (Matt.) Liebl.]. Genetics 168:1615–26CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
Vartia S, Villanueva-Cañas JL, Finarelli J, Farrell ED, Collins PC, Hughes GM et al. (2016) A novel method of microsatellite genotyping-by-sequencing using individual combinatorial barcoding. R Soc Open Sci 3:150565PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
Viruel J, Haguenauer A, Juin M, Mirleau F, Bouteiller D, Boudagher-Kharrat M et al. (2018) Advances in genotyping microsatellite markers through sequencing and consequences of scoring methods for Ceratonia siliqua (Leguminosae). Appl Plant Sci 6:e01201PubMed 
PubMed Central 
Article 

Google Scholar 
Wang J (2016) Individual identification from genetic marker data: developments and accuracy comparisons of methods. Mol Ecol Resour 16:163–175CAS 
PubMed 
Article 

Google Scholar 
Waples RS, Do C (2010) Linkage disequilibrium estimates of contemporary Ne using highly variable genetic markers: A largely untapped resource for applied conservation and evolution. Evol Appl 3:244–262PubMed 
Article 

Google Scholar 
Xie KT, Wang G, Thompson AC, Wucherpfennig JI, Reimchen TE, MacColl ADC et al. (2019) DNA fragility in the parallel evolution of pelvic reduction in stickleback fish. Science (80-) 363:81–84CAS 
Article 

Google Scholar  More

In our simulations for the autotrophs, we varied two of the three trade-off properties (level of defence, plasticity costs and defence costs; see Fig. 1b) at a time and kept the third one constant. This results in three constellations reflecting three different trade-offs between these properties (Table 1):

parallel: trade-off between defence and plasticity costs;

crossing: trade-off between defence costs and plasticity costs;

angle: trade-off between defence and defence costs.

Table 1 Description of the three constellations parallel, crossing, and angle defining the position of the four phenotypes in the trait space of defence and growth rate.Full size tableIn all three constellations, the autotrophic species B spanned the entire defence range, i.e. it had a completely undefended phenotype Bu and a maximally defended phenotype Bd. A either had a more limited defence range (in constellations parallel and angle) or spanned the entire range as well (in constellation crossing), representing three distinct ways that the trade-off between defence, growth rate, and plasticity range may play out. For each constellation, we varied the maximum switching rate χmax over 5 orders of magnitude to investigate the effect of plasticity (Table 1, middle row). This parameter determines how rapidly a species can switch between phenotypes (see “Methods”, “Exchange rates”); higher values indicate faster adaptation. These results were also compared with a non-plastic baseline scenario where both phenotypes of each species are presented but χmax = 0 (Table 1, upper row), as well as a rigid scenario where the species have only a single phenotype (Table 1, bottom row). All parameters and their values can be found in Supplementary Table S1.In the following, we give a detailed description of the results for constellation parallel, where the autotroph species A and B have the same defence costs resulting in parallel trade-off lines between defence and growth rate, while varying the level of defence for A and varying the plasticity costs for B (Table 1, left column). We start with examining patterns for the phenotype biomasses, coexistence and community stability in the non-plastic baseline scenario “parallel 0”, and then compare the corresponding scenarios with a low exchange rate (“parallel 0.01”) and a high exchange rate (“parallel 1”). We next discuss the other two constellations (crossing and angle, Table 1) more briefly. Finally, we generalize across all scenarios and focus on the coexistence, the degree of maladaptive switching, and the consumer and total autotroph biomasses.Non-plastic baseline dynamics: scenario parallel 0
In this scenario, four single phenotypes unconnected by exchange compete with each other. Thus, species coexistence here depends entirely on phenotype coexistence: the trade-offs have to be such that for each species, at least one phenotype is a good enough competitor to survive. Which phenotypes survive depend on the two trade-off parameters, defence of the defended phenotype of species A (dAu) and plasticity costs for species B (pcB), which thus determine whether coexistence is possible.The defence costs were kept constant at an intermediate value of 0.3 for both species, resulting in parallel trade-off lines (Table 1, scenario “parallel 0”). The undefended phenotype of A, Au, is a growth-specialist with the highest growth rate of all phenotypes. The defended phenotype of the same species, Ad, has a defence between 0 and 0.9 and a relatively high growth rate, and can be viewed as a generalist. Species B has variable plasticity costs that lower the growth rate of both phenotypes. The defended phenotype of species B, Bd, has the lowest growth rate of all phenotypes but is very well-defended, and thus a defence-specialist. Its undefended phenotype, Bu, is as undefended as Au but has a lower growth rate; it is thus always an inferior competitor and inevitably goes extinct (Fig. 2c).Figure 2Biomasses, coexistence and trait space for scenario parallel 0. Biomasses of the four autotrophic phenotypes (a–d), their coexistence patterns (e), the consumer biomass (f) and the autotrophs’ trait values (g–j) (higher biomasses are shown by darker colours). Lines in (a–f) separate the regions I–III of different coexistence patterns. Note that in (a–f), the y-axis is reversed to show increasing fitness along all axes. An exemplary trait combination for every region is shown in (g–j); larger symbols indicate the surviving phenotypes. Shaded areas in (e) depict oscillating systems (quarter-lag predator–prey cycles in dense shading, antiphase cycles in loose shading).Full size imageAs Bu never survives, coexistence of the autotroph species requires the survival of defence-specialist Bd. Bd can only survive if Ad is not too defended, because Ad has a higher growth rate than Bd and will outcompete Bd in the “defended” niche otherwise (region Ib; Fig. 2d,h). A second criterion is that the plasticity costs for B must not be too high, because then the benefits of the defence of Bd no longer outweigh the costs, and it will go extinct even if there are no other highly defended phenotypes around (region Ia; Fig. 2d,g). In the regions where Bd goes extinct, species coexistence is not possible (Fig. 2e). The generalist Ad either survives by itself (region Ia in Fig. 2b,g) if its defence is low to intermediate, or together with the growth-specialist Au if its defence is high (region Ib in Fig. 2a,b,h). In the regions II and III where Bd survives, it never survives on its own, but always together with one of the phenotypes of A. It coexists with the growth-specialist Au if the plasticity costs are very low (region II in Fig. 2,i), and together with Ad if they are low to intermediate (region III in Fig. 2,j). These two regions do support species coexistence (Fig. 2e).In three of the four regions (Ib, II and III in Fig. 2f), consumer biomass is low, because the final community always contains a well-defended phenotype (Ad in region Ib, and Bd in regions II and III); the overall level of defence of the community is relatively high in these regions (Supplementary Figure S1). Conversely, consumer biomass is relatively high in region Ia, because the only surviving autotroph phenotype is relatively fast-growing and fairly undefended (Fig. 2f,g). The regions where a well-defended phenotype survives often show antiphase cycles (Ib, II and III in Fig. 2e). These cycles do not occur in the region where only Ad survives (Ia in Fig. 2e); but regular quarter-lag predator–prey cycles can be found here if Ad is almost entirely undefended.While the community defence (i.e. mean defence of the autotroph community) depends strongly on the coexisting phenotypes, the community growth rate is roughly constant because over the entire trait space, at least one phenotype with a high growth rate always survives (Supplementary Figure S1). The standing variance of the community defence was high when two phenotypes coexist as they occupy different niches along the defence axis (Fig. 2h–j). In contrast, the variance of the community growth rate was very low and almost constant across all regions.Effect of phenotypic plasticityEven a little bit of plasticity in the scenario parallel 0.01 (χmax = 0.01) can change the above patterns for coexistence, stability, and average consumer biomass (Fig. 3a–d). While the autotrophs are intuitively expected to benefit from being plastic, the effect of plasticity on consumer biomass always turned out to be positive (Fig. 3a). This may be explained by the fact that switching was always, on average, maladaptive (Fig. 3c,d), measured by the adaptation index Φ (see Eqs. (11–13) in “Methods”). This index combines information on the net “flow” of individuals due to switching (i.e. whether more undefended individuals switch to defended or vice versa) with the fitness difference between the two phenotypes, and thus measures whether overall, more individuals switch from a low-fitness to a high-fitness phenotype (adaptive) or the reverse (maladaptive). This index can approach zero, but is always negative at equilibrium (see Appendix B), indicating maladaptive switching.Figure 3Consumer biomass, autotroph coexistence and maladaptive switching for the scenarios parallel 0.01 (a–d) and parallel 1 (e–h). Consumer biomass (a,e), the autotroph coexistence patterns (b,f), and the autotrophs’ maladaptive switching Φ (c,d,g,h) (higher biomasses or more intensive maladaptive switching are shown by darker colours). Lines separate the regions I–III of different autotroph coexistence. The y-axis is reversed to follow the pattern of increasing fitness. Grey areas in (c,d,g,h) depict areas where the species was extinct. Shaded areas in b and f depict oscillating systems (quarter-lag predator–prey cycles in dense shading, antiphase cycles in loose shading).Full size imageThe most striking effect of plasticity was on coexistence, which was affected both positively and negatively by plasticity in different regions of the parameter space (Fig. 3b, Supplementary Figure S4a–d). A negative effect on coexistence is seen in region II, where the autotroph species previously coexisted (Fig. 2e), while with plasticity, B outcompeted A (Fig. 3b). Without plasticity, coexistence was possible in this region because Au and Bd survived; importantly, Au outcompeted Bu due to its higher growth rate, even though the difference between their growth rates is very small in this region (Fig. 2i). Plasticity reverses the competitive exclusion pattern between the two undefended phenotypes: Bu receives a constant flow of biomass from the well-defended Bd, which compensates for its slightly lower growth rate and allows it to outcompete Au. Thus, coexistence is reduced as a direct consequence of maladaptive switching.Plasticity can also promote coexistence, as the coexistence region now extends into former region Ib where the generalist Ad is highly defended (Fig. 2b, Supplementary Figure S1a). This is also an effect of maladaptive switching, though in this case the effect is indirect, mediated through the effect of plasticity on consumer biomass. Without plasticity, coexistence was impossible in region Ib because Bd was always outcompeted by Ad: even though the latter had a slightly lower level of defence, this was outweighed by its higher growth rate, making Ad the superior competitor over Bd. However, plasticity changes this because maladaptive switching increases the consumer biomass, which in turn alters the cost/ benefit balance of defence: Bd derives a stronger benefit from its high level of defence, which now outweighs the cost and allows it to survive. Coexistence through this mechanism is not possible when the plasticity costs for B are too high or when Ad is too well-defended, explaining the narrowing of the coexistence “tail” for high defence of Ad (Fig. 3b).While the patterns of coexistence changed when allowing for plasticity, the patterns in the trait values were nearly indistinguishable from the previous scenario (Supplementary Figure S1, S2). Finally, plasticity had a strong impact on the community dynamics, as most of the antiphase cycles were stabilized (Ib, II, III in Fig. 3b). Their area decreased sharply as these cycles were characterized by asynchronous dynamics between the two prey phenotypes, which were reduced by plasticity. In contrast, the area of the quarter-lag predator–prey cycles remained unaffected by plasticity.All the above patterns were found to a far stronger degree with a higher amount of plasticity (χmax = 1; Fig. 3e–h, Supplementary Figure S4e–h). Consumer biomass increased strongly everywhere (cf. Fig. 3a,e), reflecting the strong increase in the degree of maladaptive switching (cf. Fig. 3c,d,g,h). The higher exchange rates led to more synchronization between the phenotypes, extinguishing the antiphase cycles completely (Fig. 3f). It also decreased the biomass of both defended phenotypes (cf. Supplementary Figure S4b,d,f,h). This in turn led to a lower community defence and a higher community growth rate (Supplementary Figure S3) both contributing to a higher consumer biomass. Finally, there was a sharp decrease in the coexistence region for high plasticity (Fig. 3e). Region II, where B outcompetes A through maladaptive switching, doubled in size due to the much higher degree of maladaptive switching (Fig. 3g,h). Region I, where A outcompetes B, now also increased, when the level of defence of Ad is relatively low (Fig. 3e). This is again an indirect effect of maladaptive switching causing a strong increase in consumer biomass, affecting the cost/ benefit balance of defence: while Bd derives a strong benefit from its high level of defence, Bu is completely undefended, and is at an extra disadvantage because of its low growth rate. Thus, while Bd would have been able to survive by itself, the high exchange rate causes a strong source-sink dynamic that drives B extinct.Effect of plasticity in constellations crossing and angle
In constellation crossing the trade-off lines of both species cross in the trait space, as the level of defence is the same for both defended phenotypes; species B has a lower growth rate for its undefended phenotype than species A due to plasticity costs, while its defence costs are low and thus the growth rate of its defended phenotype is higher than for species A (Table 1, Supplementary Figure S5). Without plasticity the crossing trade-off lines lead to coexistence of both species in all simulations as Au and Bd were always the only survivors, mostly showing antiphase oscillations (Supplementary Figure S5).Allowing for phenotypic plasticity has the same results as were observed for constellation parallel: consumer biomass sharply increases (Fig. 4a,e); antiphase cycles are dampened or absent; and the area of coexistence decreases (Fig. 4b,f). All these changes are more pronounced for higher exchange rates (cf. Fig. 4a,b,e,f). Again, the biomass of the defended phenotypes decreased for high exchange rates (Supplementary Figure S6). Switching was always maladaptive for high exchange rates (Fig. 4g,h), and mostly maladaptive for low exchange rates (Fig. 4c,d). As was seen for constellation parallel, maladaptive switching was the reason for the decrease in coexistence. B can outcompete A when B has low plasticity costs. Bd has a much higher growth rate than Ad, while the undefended phenotypes have similar growth rates. The direction of competitive exclusion between Au and Bu is thus easily reversed by Bd donating biomass to the sink Bu, allowing B to occupy both niches and outcompete A (region II in Fig. 4b,f). The same mechanism happens in reverse for high plasticity and defence costs of B: the differences in growth rate for the undefended phenotypes are high, while the defended phenotypes have very similar growth rates. Au can support Ad, and A outcompetes B (region III in Fig. 4b,f).Figure 4Coexistence and maladaptive switching for scenario crossing 0.01 (a-d) and crossing 1 (e–h). Consumer biomass (a,e), the autotroph coexistence patterns (b,f), and the autotrophs’ maladaptive switching Φ (c,d,g,h) (higher biomasses or more intensive maladaptive switching are shown by darker colours). Lines separate the regions I–III of different autotroph coexistence. The x- and y-axis are reversed to follow the pattern of increasing fitness. Shaded areas in (b) depict antiphase cycles. Grey areas in (c,d,g,h) depict areas where the species was extinct. Shaded grey areas depict areas without simulations (cf. “Methods”). Note that (c,d,g,h) have each a different colour scale.Full size imageIn constellation angle there are no plasticity costs, and thus the undefended phenotypes Au and Bu have identical growth rates. The defended phenotypes take the same places in trait space as in the parallel constellation: Ad is a generalist, with a lower level of defence and a relatively high growth rate due to low defence costs, whereas Bd is a defence-specialist with a high level of defence but a low growth rate. This leads to the trade-off lines forming an angle (see Table 1). Without phenotypic plasticity, the coexistence patterns are the same as in constellation parallel, except that no competitive exclusion occurs between the undefended phenotypes; instead, they (neutrally) coexist in regions Ib, II and III (Supplementary Figure S7; cf. Fig. 2).With plasticity, neutral coexistence vanished: the defended phenotype that survived (Ad in region Ib, Bd in region III) could support the undefended phenotype of its own species, driving the other species extinct (Fig. 5b,f). As in the other constellations, the area of coexistence and the biomasses of the defended phenotypes decreased and antiphase cycles vanished with increasing χmax (Fig. 5b,f, Supplementary Figure S8), while maladaptive switching and the consumer biomass increased (Fig. 5).Figure 5Coexistence and maladaptive switching for scenario angle 0.01 (a–d) and angle 1 (e–h). Consumer biomass (a,e), the autotroph coexistence patterns (b,f), and the autotrophs’ maladaptive switching Φ (c,d,g,h) (higher biomasses or more intensive maladaptive switching are shown by darker colours). Lines separate the regions I–III of different autotroph coexistence. The y-axis is reversed to follow the pattern of increasing fitness. Shaded areas in (b) depict antiphase cycles. Grey areas in (c,d,g,h) depict areas where the species was extinct. Note that (c,d,g,h) have each a different colour scale.Full size imageGeneral resultsAs plasticity had very similar effects across all three constellations, we here generalize our results: we compare the three constellations for exchange rates over 5 orders of magnitude, as well as the non-plastic scenario and the rigid scenario (Table 1). That is, all simulations from one scenario (e.g. parallel 0) were summarized into one bar respective point in Fig. 6.Figure 6General patterns for coexistence, maladapative switching and biomasses. Share of surviving species in percent (A, B, coexistence or neutral coexistence) (a–c), median absolute value of maladaptive switching Φ (d–f) and median of total autotroph biomass (A + B), median consumer biomass C and share of defended phenotypes ((Ad + Bd)/(A + B)) (g–i) for the three constellations and increasing maximum exchange rates χmax. χmax = 0 denotes the non-plastic scenarios; *denotes the rigid scenarios. Maladaptive switching and the share of defended phenotypes do not apply for the rigid scenarios.Full size imageFor all constellations, the fraction of simulation runs leading to coexistence was highest in the non-plastic scenario and decreased with increasing χmax (Fig. 6a–c). In constellation parallel the share of coexistence for increasing χmax continuously decreased from 51 to 3% (Fig. 6a). In crossing, the share decreased from full to no coexistence (Fig. 6b). In angle, the share of coexistence was 88% in the non-plastic scenario when taking also neutral coexistence into account (Fig. 6c). Its share decreased to 9% for a χmax of 10 and increased again to 25% for the rigid scenario. Maladaptive switching increased for both species and all constellations for increasing χmax (Fig. 6d–f). The increased plasticity led to a lower total autotroph biomass and a lower share of defended phenotypes (Fig. 6g–i), which resulted in higher consumer biomass (Fig. 6g–i).Interestingly, and counterintuitively, the above patterns show that increasing the speed of plasticity (by increasing χmax) makes the system behave more like the rigid system. The coexistence patterns in scenarios with high χmax approach those of the rigid scenarios in two of the constellations (Fig. 6a,b). Similarly, the total autotroph and consumer biomasses approach the ones in the rigid scenarios (Fig. 6g–i). Thus, we found the higher χmax make the autotrophs not more adaptive, but behave more like non-adaptive species. More

Insects and plantsLarvae of fall armyworm (FAW, Spodoptera frugiperda) were cultured at the Department of Entomology at the Max Planck Institute for Chemical Ecology, and reared on a semi-artificial diet based on pinto bean59, and maintained under controlled light and temperature conditions (12:12 h light/dark, 21 °C).Feeding experiments3rd–4th instar FAW larvae were utilized for all experiments. Insects were starved overnight prior to feeding experiments. The following day insects were fed with a semi-artificial, pinto bean-based diet or put on maize leaves in small plastic cups and allowed to feed on the respective diets for a day. Insects were dissected in cold phosphate buffered saline (PBS, pH = 7.4) to harvest larval tissues (guts, Malphigian tubules, fat bodies, cuticle), which were stored at − 80 °C until further use. For droplet feeding, 12.5 mM DIMBOA was prepared by dissolving the compound in DMSO. This DIMBOA solution was further diluted in 10% aqueous sucrose solution. The larvae were stimulated with forceps to encourage regurgitation, and 2 μL DIMBOA-sucrose solution was administered directly to the larval mouthparts. Insects were then fed on semi-artificial diet for up to 6 h; following which gut tissue was dissected using cold phosphate buffer and the tissue samples were stored at − 80 °C until further use.Insect cell culturesSpodoptera frugiperda Sf9 cells and Trichoplusia ni Hi5 cells were cultured in Sf-900 II serum-free medium (Gibco) and ExpressFive serum-free medium (Gibco), respectively. Adherent cultures were maintained at 27 °C, and sub-cultured every 3–4 days.Cell treatmentsInsect cells were seeded in 6 well culture plates (Corning) and left at 27 °C overnight. For transcript stability tests, a fresh cycloheximide (CHX) stock (50 mg/mL) was prepared in ethanol and added to the cultured cells at a concentration of 50 µg/mL. Incubations with CHX were performed up to 6 h. For testing substrate specificity, cells were then treated with the following compounds for 1 h—DIMBOA (25–100 μM), indole (50–100 μM), quercetin (50–100 μM), and esculetin (50–100 μM). All the stocks were prepared in DMSO and cells treated with the corresponding volume of pure DMSO served as a control. The range of concentrations used for the substrates was based on previous work38.RNA extraction, reverse transcription and real time-PCR analysisTissue samples from the larvae were homogenized and total RNA extracted using the innuPREP RNA Mini Kit (Analytik Jena). Cell cultures used for RNA extraction were obtained during sub-culturing at full confluency, and centrifuged at 500×g for 5 min. The culture medium was discarded, and the fresh pellets were directly used for RNA extraction. RNA concentrations were measured with the NanoDrop 2000 UV–Vis Spectrophotometer (Thermo Scientific). First strand cDNA was synthesized from 1 μg total RNA using SuperScript III Reverse Transcriptase and OligodT primers from Invitrogen. Sequences were successfully amplified using Phusion High Fidelity DNA Polymerase (New England Biolabs) (PCR protocol: 30 s at 98 °C; 35 cycles of 10 s at 98 °C, 20 s at 55 °C, 45 s at 72 °C; and 5 min at 72 °C). The PCR products were purified with a PCR cleanup kit (Qiagen) and cloned into pCR-Blunt II-TOPO vector (Life Technologies) and transformed into NEB cells (Life Technologies), which were plated on selective LB agar medium containing 100 μg/mL ampicillin and incubated overnight at 37 °C. Positive colonies were identified by PCR using vector-specific M13 primers. Positive clones were confirmed by sequencing. Real time PCR analyses were carried out using Brilliant III SYBR Master Mix, employing SYBR Green chemistry. Relative quantification of the transcript levels was done using the 2−∆∆Ct method60. SfRPL10 was used as reference gene for all analyses. The primer pairs used for distinguishing between the variants are listed in Supplementary Table 1. As the expression of full-length and variants of SfUGT33F28 differed according to the strains, tissues, and treatments being analyzed, variant expression is reported as ratios relative to the canonical transcript to facilitate comparisons.Preparation of minigenes for alternative splicing studiesGenomic DNA was isolated from S. frugiperda larvae using the cetyl trimethyl ammonium bromide (CTAB) protocol61. DNA concentration was measured with the NanoDrop 2000 UV–Vis Spectrophotometer (Thermo Scientific). The minigene was amplified using Phusion High Fidelity DNA Polymerase (New England Biolabs) (PCR protocol: 30 s at 98 °C; 35 cycles of 10 s at 98 °C, 30 s at 55–60 °C, 1 min 30 s at 72 °C; and 10 min at 72 °C), cloned into a pCR-Blunt II-TOPO vector (Life Technologies) and sequenced using M13 primers. The confirmed sequence was eventually cloned into a pIB/V5-His-TOPOvector (Life Technologies) and transformed into NEB cells (Life Technologies). Positive colonies were identified by colony PCR using vector-specific OpIE2 primers, sub-cultured overnight at 37 °C in liquid LB medium containing 100 μg/mL ampicillin and used for plasmid DNA purification with the NucleoSpin Plasmid kit (Macherey-Nagel). Concentration and purity of the obtained construct was assessed by the NanoDrop 2000 UV–Vis Spectrophotometer (Thermo Scientific) and the correct orientation of the PCR products was confirmed by DNA sequencing.Nuclear protein isolationNuclear proteins were isolated from insect cells62 using the protocol originally described with few modifications. Cells grown to concentrations of up to 1 × 106 cells/well were harvested and washed with PBS (pH 7.4). The extracts were centrifuged at 12,000×g for 10 min and pellets were re-suspended in 400 μL cell lysis buffer (10 mM HEPES, pH 7.5, 10 mM KCl, 0.1 mM EDTA pH 8.0, 1 mM DTT, 0.5% Nonidet-40 and 10 μL protease inhibitor cocktail). Cells were allowed to swell on ice for 20 min with intermittent mixing. Suspensions were vortexed to disrupt the cell membranes and then centrifuged at 12,000×g for 10 min at 4 °C. Pelleted nuclei were washed thrice with cell lysis buffer, re-suspended in 50 μL nuclear extraction buffer (20 mM HEPES pH 7.5, 400 mM KCl, 1 mM EDTA pH 8.0, 1 mM DTT, 10% glycerol and protease inhibitor) and incubated on ice for 30 min. Nuclear fractions were collected by centrifugation at 12,000g for 15 min at 4 °C. Protein concentrations were measured by Bradford and extracts were stored at − 80 °C until further use.Electrophoretic mobility shift assay (EMSA)EMSA was performed using the LightShift Chemiluminescent EMSA kit (Thermo Scientific) following the manufacturer’s instructions. Genomic DNA fragments of 20–25 bp corresponding to the 5′ flanking region of UGT33F28 exon 1 (with and without AhR-ARNT motif deletion) were synthesized with covalently linked biotin (Sigma). The DNA probes used in the experiment are listed in Supplementary Table 6. EMSA was performed in 20 µL reactions containing 20 fmol biotinylated DNA probe with 3.5–4 µg nuclear protein extracted from insect cells, according to manufacturer’s instructions. A reaction comprising the above along with the excess of unlabeled canonical DNA probe (200 molar excess) was further employed as a control. The reaction was assembled at room temperature and incubated for 30 min. The reactions were separated on a 5% TBE gel in 0.5X TBE at 100 V for 60 min. The samples were then transferred to a positively charged nylon membrane (Hybond N+, Amersham Bioscience) using semi-dry transfer at 15 V for 30 min. The membrane was cross-linked for 1 min using the auto cross-link function on the UV cross-linker (Stratagene). The biotinylated DNA–protein complex was detected by the streptavidin–horseradish peroxidase conjugated antibody provided in the kit. The membrane was washed and incubated with the chemiluminescence substrate for 5 min and the signals were developed by exposing the membrane to an X-ray film for 1 min.Streptavidin affinity purificationStreptavidin agarose (Sigma-Aldrich) was employed for protein purification. Briefly, 50–100 μL of agarose was packed into a 1.5 mL Eppendorf tube for each sample. The agarose was allowed to settle with a short centrifugation (500×g, 5 min) and the supernatant was discarded. The agarose was washed 4–5 times with binding buffer (PBS containing 1 mM EDTA, 1 mM DTT, 4 µg poly dI. dC as non-specific competitor DNA and protease inhibitor). Simultaneously, the binding reaction with the nuclear protein fraction and the DNA probe was assembled as described above. A 100 μg amount of total nuclear protein was incubated with 4 μg of biotinylated DNA probe at room temperature for 20 min. The reaction was loaded onto the streptavidin column equilibrated with the binding buffer and incubated for another 1 h at room temperature with gentle shaking. Subsequently, the agarose was washed 4–5 times with the binding buffer. After the final wash, the supernatant was aspirated and 10 μL was left above the beads. For protein separation, 20–30 μL pf the SDS loading buffer was added onto the agarose, boiled at 95 °C for 5 min and the sample thus obtained was utilized for electrophoresis.Deletion mutagenesisFor deletion mutagenesis, a pair of primers flanking the sequence to be deleted (non-overlapping) was designed. The pCR-Blunt II-TOPO vector (Life Technologies) clone for the SfUGT33F28 exon 1–2 minigene was utilized as a template. Sequence was successfully amplified using Phusion High Fidelity DNA Polymerase (New England Biolabs) (PCR protocol: 30 s at 98 °C; 20 cycles of 10 s at 98 °C, 30 s at 55–60 °C, 4 min at 72 °C; and 10 min at 72 °C). A DpnI digest was performed to remove the background DNA, followed by ligation and transformation into fresh cells. The sequence of the mutant TOPO clone was then confirmed and utilized as a template for cloning into pIB/V5-His-TOPO vector (Life Technologies) for transfection into insect cells.Cloning and heterologous expression of SfUGTsSequences were amplified from S. frugiperda gut cDNA samples using Phusion High Fidelity DNA Polymerase (New England Biolabs) (PCR protocol: 30 s at 98 °C; 35 cycles of 10 s at 98 °C, 20 s at 55–60 °C, 45 s at 72 °C; and 5 min at 72 °C). The resulting amplified products were purified with a PCR cleanup kit (Qiagen) and incubated with GoTaq DNA polymerase (Promega) for 15 min at 72 °C in order to add A overhangs. The products were cloned into the pIB/V5-His-TOPO vector (Life Technologies) and transformed into NEB cells (Life Technologies), which were plated on selective LB agar medium containing 100 μg/mL ampicillin and incubated overnight at 37 °C. Positive colonies were identified by PCR using vector-specific OpIE2 primers, sub-cultured overnight at 37 °C in liquid LB medium containing 100 μg/mL ampicillin and used for plasmid DNA purification with the NucleoSpin Plasmid kit (Macherey-Nagel). Concentration and purity of the obtained constructs were assessed by NanoDrop 2000 UV–Vis Spectrophotometer (Thermo Scientific) and the correct orientation of the PCR products was confirmed by DNA sequencing.Insect cell transfectionFor transfection, Sf9 cells and Hi5 cells were sub-cultured at full confluency in a 6-well plate in a 1:3 dilution and left overnight to adhere to the flask surface. The medium was replaced, and transfections were carried out using FuGENE HD Transfection Reagent (Promega) in a 1:3 plasmid/lipid ratio (1.7 μg plasmid and 5.0 μL lipid for 3 mL medium). Cells were incubated for 48–72 h at 27 °C and re-suspended in fresh medium containing 50 μg/mL blasticidin for 2 weeks. Stable cell cultures were subsequently maintained at 10 μg/mL blasticidin.Cell lysate preparationCells were obtained from cultures 2 weeks post transfection growing stably on 50 μg/mL blasticidin. A 1 mL quantity of cells was harvested for each construct and re-suspended into 100 µL buffer. Protein concentrations were measured using the Bradford reagent, and 1–2 μg of the cell lysate was used for enzyme assays.Microsome preparationFor microsome extraction, confluent, stably transfected cells from five T-75 flasks (10 mL culture) per recombinant plasmid were harvested by scraping the cells off the bottom using a sterile cell scraper (Sarstedt AG, Nuembrecht, Germany). The obtained cell suspensions were combined into a 50 mL falcon tube and centrifuged at 1000×g for 15 min at 4 °C (AvantiTM J-20 XP Centrifuge, Beckman Coulter, Krefeld, Germany). The supernatant was discarded, the cells were washed twice with ice-cold PBS buffer (pH 7.4) and centrifuged at 1000×g for 15 min. The resulting cell pellet was re-suspended in 10 mL hypotonic buffer (20 mM Tris, 5 mM EDTA, 1 mM DTT, 20% glycerol, pH 7.5), containing 0.1% BenzonaseR nuclease and 100 μL Protease Inhibitor Cocktail (Serva) followed by incubation on ice for 30 min. After cell lysis, the cells were homogenized by 20–30 strokes in a Potter–Elvehjem tissue grinder (Kontes Glass Co., Vineland, USA) and were subsequently mixed with an equal volume of sucrose buffer (20 mM Tris, 5 mM EDTA, 1 mM DTT, 500 mM sucrose, 20% glycerol, pH 7.5). The homogenate was centrifuged at 1200×g and 4 °C for 10 min (AvantiTM J-20 XP Centrifuge, Beckman Coulter), and the supernatant was transferred into Beckman polycarbonate ultracentrifugation bottles (25 × 89 mm) (Beckman Coulter) and centrifuged at 100,000×g and 4 °C for 1.5 h in a fixed angle Type 70 Ti rotor (OptimaTM L-90K Ultracentrifuge, Beckman Coulter). After ultracentrifugation, the clear supernatant, containing the cytosolic fraction, was aliquoted into 1.5 mL Eppendorf tubes. The pellet, containing the microsomal fractions, was re-suspended in 1 mL of phosphate buffer (100 mM K2HPO4, pH 7.0), containing 10 μL Protease Inhibitor Cocktail (Serva) and stored at − 80 °C until further use. Typically, 5–10 μg of the microsome fraction so obtained was utilized for the enzyme assays.Cross-linking assaysCross-linking assays were performed using dimethyl suberimidate (DMS) as the cross-linking agent. A fresh stock of DMS (5 mg/mL) was prepared in 0.2 M triethanolamine (pH 8.0) at the start of each assay. DMS was added to a final concentration of 2.5 mg/mL to insect cell microsomes with gentle shaking up to 3 h, and samples were subsequently stored at − 20 °C until further use. All protein samples were electrophoresed using a 12% Mini-PROTEAN tris glycine gel, blotted onto PVDF membrane using wet transfer at 70 V for 30–45 min, followed by detection using the V5-HRP conjugate.V5-based affinity purificationAnti-V5 agarose affinity gel (Sigma-Aldrich) was employed for protein purification. Briefly, 50–75 μL of the agarose was packed into a 1.5 mL Eppendorf tube for each sample. The agarose was allowed to settle with a short centrifugation and the supernatant was discarded. The agarose was washed 4–5 times with PBS (pH 7.4). Samples to be purified were incubated with 5% digitonin on ice for 20 min and subject to centrifugation at 16,000×g for 30 min. Clarified cell lysate or microsomal extract was added onto the resin (up to 200 μL, volume adjusted by addition of PBS) and incubated for 1.5 h on a shaker. Subsequently, the agarose was washed 4–5 times with PBS. After the final wash, the supernatant was aspirated and 10 μL was left above the beads. This fraction was used for both protein electrophoresis and enzyme assays (separate purifications). For SDS-PAGE, 20–30 μL pf the SDS loading buffer was added onto the agarose, boiled at 95 °C for 5 min and sample thus obtained was utilized for electrophoresis.LC–MS/MS peptide analysisProtein bands of Coomassie Brilliant blue R250 stained gels were cut from the gel matrix and tryptic digestion was carried out63. For LC–MS/MS analysis of the resulting peptides, samples were reconstituted in 20 μL aqueous 1% formic acid, and 1 μL was injected onto an UPLC M-class system (Waters, Manchester, UK) coupled to a Synapt G2-si mass spectrometer (Waters, Manchester, UK). Samples were first pre-concentrated and desalted using a Symmetry C18 trap column (100 Å, 180 µm × 20 mm, 5 µm particle size) at a flow rate of 15 µL/min (0.1% aqueous formic acid). Peptides were eluted onto a ACQUITY UPLC HSS T3 analytical column (100 Å, 75 µm × 200 mm, 1.8 µm particle size) at a flow rate of 350 nL/min with the following gradient: 3–15% over 3 min, 15–20% B over 7 min, 20–40% B over 30 min, 40–50% B over 5 min, 50–70% B over 5 min, 70–95% B over 3 min, isocratic at 95% B for 1 min, and a return to 1% B over 1 min. Phases A and B were composed of 0.1% formic acid and 100% acetonitrile in 0.1% formic acid, respectively). The analytical column was re-equilibrated for 10 min prior to the next injection. The eluted peptides were transferred into the mass spectrometer operated in V-mode with a resolving power of at least 20,000 full width at half height FWHM. All analyses were performed in a positive ESI mode. A 100 fmol/μL sample of human Glu-Fibrinopeptide B in 0.1% formic acid/acetonitrile (1:1 v/v) was infused at a flow rate of 1 μL/min through the reference sprayer every 45 s to compensate for mass shifts in MS and MS/MS fragmentation mode. Data were acquired using data-dependent acquisition (DDA). The acquisition cycle for DDA analysis consisted of a survey scan covering the range of m/z 400–1800 Da followed by MS/MS fragmentation of the ten most intense precursor ions collected at 0.5 s intervals in the range of 50–2000 m/z. Dynamic exclusion was applied to minimize multiple fragmentations for the same precursor ions. MS data were collected using MassLynx v4.1 software (Waters, Manchester, UK).Data processing and protein identificationDDA raw data were processed and searched against a sub-database containing common contaminants (human keratins and trypsin) using ProteinLynx Global Server (PLGS) version 2.5.2 (Waters, Manchester, UK). Spectra remaining unmatched by database searching were interpreted de novo to yield peptide sequences and subjected to homology-based searching using the MS BLAST program64 installed on a local server. MS BLAST searches were performed against a Spodoptera frugiperda database obtained by in silico translation of the S. frugiperda transcriptome37 and against arthropoda database (NCBI). PKL-files of MS/MS spectra were generated and searched against Spodoptera frugiperda database combined with NCBI nr (downloaded on May 24, 2020) using MASCOT software version 2.6.2. The following searching parameters were applied: fixed precursor ion mass tolerance of 15 ppm for the survey peptide, fragment ion mass tolerance of 0.1 Da, 1 missed cleavage, fixed carbamidomethylation of cysteines and possible oxidation of methionine.Enzymatic assaysFor UGT assays, samples from insect cell cultures (transient or stable) were prepared in phosphate buffer (pH 7.0, 100 mM). Typical enzyme reactions included 5–10 µg cell microsomal extracts, 2 μL of 12.5 mM DIMBOA in DMSO (25 nmol), 4 μL of 12.5 mM UDP-glucose in water (50 nmol), and phosphate buffer (pH 7.0, 100 mM) to give an assay volume of 50 μL. Controls containing either boiled enzymatic preparation, or only the protein suspension and buffer were included. After incubation at 30 °C for 60 min, the enzyme reactions were interrupted by adding 50 μL of 1:1 (v:v) methanol/formic acid solution. For enzyme assays involving resin purified microsomal extracts, equal amounts of extracts were employed for resin purification and the enzyme assay (buffer + substrate) was pipetted directly onto the resin. Post incubation, samples were centrifuged, supernatant was collected, and reaction was stopped by addition of methanol/formic acid solution. Assays were centrifuged at 5000g for 5 min and the obtained supernatant was collected and analyzed by LC–MS/MS.Chromatographic methodsFor all analytical chromatography procedures, formic acid (0.05%) in water and acetonitrile were used as mobile phases A and B, respectively, and the column temperature was maintained at 25 °C. Analyses of enzymatic assays and plant samples used a Zorbax Eclipse XDB-C18 column (50 × 4.6 mm, 1.8 μm, Agilent Technologies) with a flow rate of 1.1 mL/min and with the following elution profile: 0–0.5 min, 95% A; 0.5–6 min, 95–67.5% A; 6.02–7 min, 100% B; 7.1–9.5 min, 95% A. LC–MS/MS analyses were performed on an Agilent 1200 HPLC system (Agilent Technologies) coupled to an API 6500 tandem spectrometer (AB Sciex) equipped with a turbospray ion source operating in negative ionization mode. Multiple reaction monitoring (MRM) was used to monitor analyte parent ion to product ion conversion with parameters from the literature for DIMBOA65 and DIMBOA-Glc16. Analyst (version 1.6.3, Applied Biosystems) software was used for data acquisition and processing.Statistical analysisAll statistical analyses were carried out using SigmaPlot 12.0 and R studio (version 3.6.3). Data were tested for homogeneity of variance and normality and were appropriately transformed to meet these criteria where required. The specific statistical method used for each data set is described in the figure legends. More

Wink, M. Evolution of secondary metabolites from an ecological and molecular phylogenetic perspective. Phytochemistry 64, 3–19 (2003).CAS 
PubMed 
Article 

Google Scholar 
Khare, S. et al. Plant secondary metabolites synthesis and their regulations under biotic and abiotic constraints. J. Plant Biol. 63, 203–216 (2020).CAS 
Article 

Google Scholar 
Gajger, I. T. & Dar, S. A. Plant allelochemicals as sources of insecticides. Insects 12, 189 (2021).Article 

Google Scholar 
Vig, A. P., Rampal, G., Thind, T. S. & Arora, S. Bio-protective effects of glucosinolates: A review. LWT Food Sci. Technol. 42, 1561–1572 (2009).CAS 
Article 

Google Scholar 
Sikorska-Zimny, K. & Beneduce, L. The glucosinolates and their bioactive derivatives in Brassica: A review on classification, biosynthesis and content in plant tissues, fate during and after processing, effect on the human organism and interaction with the gut microbiota. Crit. Rev. Food Sci. Nutr. 61, 2544–2571. https://doi.org/10.1080/10408398.2020.1780193 (2020).CAS 
Article 
PubMed 

Google Scholar 
Radojčić Redovniković, I., Glivetić, T., Delonga, K. & Vorkapić-Furač, J. Glucosinolates and their potential role in plant. Period. Biol. 110, 297–309 (2008).
Google Scholar 
Wittstock, U., Kliebenstein, D. J., Lambrix, V., Reichelt, M. & Gershenzon, J. Chapter five glucosinolate hydrolysis and its impact on generalist and specialist insect herbivores. Recent Adv. Phytochem. 37, 101–125 (2003).CAS 
Article 

Google Scholar 
Noret, N. et al. Palatability of Thlaspi caerulescens for snails: Influence of zinc and glucosinolates. New Phytol. 165, 763–772 (2005).CAS 
PubMed 
Article 

Google Scholar 
Hopkins, R. J., Van Dam, N. M. & Van Loon, J. J. A. Role of glucosinolates in Insect-plant relationships and multitrophic interactions. Annu. Rev. Entomol. 54, 57–83 (2008).Article 
CAS 

Google Scholar 
Guleria, S. & Tiku, A. K. Botanicals in pest management: Current status and future perspectives. Integr. Pest Manag. 1, 317–329 (2009).
Google Scholar 
Clay, N. K., Adio, A. M., Denoux, C., Jander, G. & Ausubel, F. M. Glucosinolate metabolites required for an Arabidopsis innate immune response. Science 323, 95–101 (2009).ADS 
CAS 
PubMed 
Article 

Google Scholar 
Yang, B. et al. Inhibitory effect of allyl and benzyl isothiocyanates on ochratoxin a producing fungi in grape and maize. Food Microbiol. 100, 103865 (2021).CAS 
PubMed 
Article 

Google Scholar 
Agrawal, A. A. & Kurashige, N. S. A role for isothiocyanates in plant resistance against the specialist herbivore Pieris rapae. J. Chem. Ecol. 296(29), 1403–1415 (2003).Article 

Google Scholar 
Müller, C. et al. The role of the glucosinolate-myrosinase system in mediating greater resistance of Barbarea verna than B. vulgaris to Mamestra brassicae Larvae. J. Chem. Ecol. 44, 1190–1205 (2018).PubMed 
Article 
CAS 

Google Scholar 
Kegley, S. E., Hill, B. R., Orme, S. & Choi, A. H. PAN Pesticide Database (Pesticide Action Network, 2000).
Google Scholar 
Worfel, R. C., Schneider, K. S. & Yang, T. C. S. Suppressive effect of allyl isothiocyanate on populations of stored grain insect pests. J. Food Process. Preserv. 21, 9–19 (1997).CAS 
Article 

Google Scholar 
Wu, H., Zhang, G. A., Zeng, S. & Lin, K. C. Extraction of allyl isothiocyanate from horseradish (Armoracia rusticana) and its fumigant insecticidal activity on four stored-product pests of paddy. Pest Manag. Sci. 65, 1003–1008 (2009).CAS 
PubMed 
Article 

Google Scholar 
Bhushan, S., Gupta, S., Kaur Sohal, S., Arora, S. & Saroj Arora, C. Assessment of insecticidal action of 3-Isothiocyanato-1-propene on the growth and development of Spodoptera litura (Fab.) (Lepidoptera: Noctuidae). J. Entomol. Zool. Stud. 4, 1068–1073 (2016).
Google Scholar 
Dhillon, M. K., Naresh, J. S., Singh, R. & Sharma, N. K. Reaction of different bitter gourd (Momordica charantia L.) genotypes to melon fruit fly, Bactrocera cucurbitae (Coquillett). Int. J. Plant Prot. 33, 55–59 (2005).
Google Scholar 
Ekesi, S., Nderitu, P. W. & Chang, C. L. Adaptation to and small-scale rearing of invasive fruit fly Bactrocera invadens (Diptera: Tephritidae) on artificial diet. Ann. Entomol. Soc. Am. 100, 562–567 (2007).Article 

Google Scholar 
Jakhar, S. et al. Estimation losses due to fruit fly, Bactrocera cucurbitae (Coquillett) on long melon in semi-arid region of Rajasthan. J. Entomol. Zool. Stud. 8, 632–635 (2020).MathSciNet 

Google Scholar 
Ladania, M. S. Physiological Disorders and Their Management. Citrus Fruit: Biology, Technology and Evaluation 451–463 (Academic press, 2008).
Google Scholar 
Du, Y., Grodowitz, M. J. & Chen, J. Insecticidal and enzyme inhibitory activities of isothiocyanates against red imported fire ants, Solenopsis invicta. Biomolecules 10, 716 (2020).CAS 
PubMed Central 
Article 

Google Scholar 
Tsao, R., Reuber, M., Johnson, L. & Coats, J. R. Insecticidal toxicities of glucosinolate· containing extracts from crambe seeds. J. Agric. Urban Entomol. 13, 109–120 (1996).CAS 

Google Scholar 
Li, Q., Eigenbrode, S. D., Stringam, G. R. & Thiagarajah, M. R. Feeding and growth of Plutella xylostella and Spodoptera eridania on Brassica juncea with varying glucosinolate concentrations and myrosinase activities. J. Chem. Ecol. 26, 2401–2419 (2000).CAS 
Article 

Google Scholar 
Noble, R. R., Harvey, S. G. & Sams, C. E. Toxicity of Indian mustard and allyl isothiocyanate to masked chafer beetle larvae. Plant Health Prog. 3, 9 (2002).Article 

Google Scholar 
Sousa, A. H., Faroni, L. R. A., Pimentel, M. A. G. & Freitas, R. S. Relative toxicity of mustard essential oil to insect-pests of stored products. Rev. Caatinga 27, 222–226 (2014).
Google Scholar 
de Souza, L. P., Faroni, L. R. D. A., Lopes, L. M., de Sousa, A. H. & Prates, L. H. F. Toxicity and sublethal effects of allyl isothiocyanate to Sitophilus zeamais on population development and walking behavior. J. Pest Sci. 91, 761–770 (2018).Article 

Google Scholar 
Freitas, R. C. P., Faroni, L. R. D. A., Haddi, K., Jumbo, L. O. V. & Oliveira, E. E. Allyl isothiocyanate actions on populations of Sitophilus zeamais resistant to phosphine: Toxicity, emergence inhibition and repellency. J. Stored Prod. Res. 69, 257–264 (2016).Article 

Google Scholar 
Jabeen, A., Zaitoon, A., Lim, L. T. & Scott-Dupree, C. Toxicity of five plant volatiles to adult and egg stages of Drosophila suzukii matsumura (Diptera: Drosophilidae), the spotted-wing Drosophila. J. Agric. Food Chem. 69, 9511–9519 (2021).CAS 
PubMed 
Article 

Google Scholar 
Wu, H., Liu, X. R., Yu, D. D., Zhang, X. & Feng, J. T. Effect of allyl isothiocyanate on ultra-structure and the activities of four enzymes in adult Sitophilus zeamais. Pestic. Biochem. Physiol. 109, 12–17 (2014).CAS 
PubMed 
Article 

Google Scholar 
Zhang, C., Wu, H., Zhao, Y., Ma, Z. & Zhang, X. Comparative studies on mitochondrial electron transport chain complexes of Sitophilus zeamais treated with allyl isothiocyanate and calcium phosphide. Pestic. Biochem. Physiol. 126, 70–75 (2016).CAS 
PubMed 
Article 

Google Scholar 
Jeschke, V. et al. How glucosinolates affect generalist lepidopteran larvae: Growth, development and glucosinolate metabolism. Front Plant Sci. 8, 1995. https://doi.org/10.3389/fpls.2017.01995 (2017).Article 
PubMed 
PubMed Central 

Google Scholar 
Agnihotri, A. R., Hulagabali, C. V., Adhav, A. S. & Joshi, R. S. Mechanistic insight in potential dual role of sinigrin against Helicoverpa armigera. Phytochemistry 145, 121–127. https://doi.org/10.1016/j.phytochem.2017.10.014 (2018).CAS 
Article 
PubMed 

Google Scholar 
Jeschke, V. et al. So much for glucosinolates: A generalist does survive and develop on Brassicas, but at what cost?. Plants 10, 962. https://doi.org/10.3390/plants10050962 (2021).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
Benrey, B. & Denno, R. F. The slow-growth-high-mortality hypothesis: A test using the cabbage butterfly. Ecology 78, 987–999 (1997).
Google Scholar 
Shroff, R., Vergara, F., Muck, A., Svatoš, A. & Gershenzon, J. Nonuniform distribution of glucosinolates in Arabidopsis thaliana leaves has important consequences for plant defense. Proc. Natl. Acad. Sci 05, 6196–6201 (2008).ADS 
Article 
CAS 

Google Scholar 
Bai, P. P. et al. Inhibition of phenoloxidase activity delays development in Bactrocera dorsalis (Diptera: Tephritidae). Fla. Entomol. 97, 477–485. https://doi.org/10.1653/024.097.0218 (2014).Article 

Google Scholar 
Datta, R., Kaur, A., Saraf, I., Singh, I. P. & Kaur, S. Effect of ethyl acetate extract and purified compounds of Alpinia galanga (L.) on Immune Response of a Polyphagous Lepidopteran pest, Spodoptera litura (Fabricius). Asian J. Adv. Basic Sci. 6, 16–21 (2018).CAS 

Google Scholar 
Hartzer, K. L., Zhu, K. Y. & Baker, J. E. Phenoloxidase in larvae of Plodia interpunctella (Lepidoptera: Pyralidae): Molecular cloning of the proenzyme cDNA and enzyme activity in larvae paralyzed and parasitized by Habrobracon hebetor (Hymenoptera: Braconidae). Arch. Insect Biochem. Physiol. 59, 67–79 (2005).CAS 
PubMed 
Article 

Google Scholar 
Silva, C. J. M. et al. Immune response triggered by the ingestion of polyethylene microplastics in the dipteran larvae Chironomus riparius. J. Hazard. Mater. 414, 125401. https://doi.org/10.1016/j.jhazmat.2021.125401 (2021).CAS 
Article 
PubMed 

Google Scholar 
Aucoin, R. R., Philogène, B. J. R. & Arnason, J. T. Antioxidant enzymes as biochemical defenses against phototoxin induced oxidative stress in three species of herbivorous Lepidoptera. Arch. Insect Biochem. Physiol. 16, 139–152 (1991).CAS 
Article 

Google Scholar 
Wang, Y., Branicky, R., Noë, A. & Hekimi, S. Superoxide dismutases: Dual roles in controlling ROS damage and regulating ROS signaling. Int. J. Cell Biol. 217, 1915–1928. https://doi.org/10.1083/jcb.201708007 (2018).CAS 
Article 

Google Scholar 
Zhang, C., Ma, Z., Zhang, X. & Wu, H. Transcriptomic alterations in Sitophilus zeamais in response to allyl isothiocyanate fumigation. Pest. Biochem. Physiol. 137, 62–70. https://doi.org/10.1016/j.pestbp.2016.10.001 (2017).CAS 
Article 

Google Scholar 
Felton, G. W. & Summers, C. B. Antioxidant systems in insects. Arch. Insect Biochem. Physiol. 29, 187–197. https://doi.org/10.1002/arch.940290208 (1995).CAS 
Article 
PubMed 

Google Scholar 
Cadenas, E. Mechanisms of oxygen activation and reactive oxygen species detoxification. In Oxidative Stress and Antioxidant Defenses in Biology (ed. Ahmad, S.) 1–46 (Chapman & Hall, 1995). https://doi.org/10.1007/978-1-4615-9689-9_1.Chapter 

Google Scholar 
Schramm, K., Vassão, D. G., Reichelt, M., Gershenzon, J. & Wittstock, U. Metabolism of glucosinolate- derived isothiocyanates to glutathione conjugates in generalist lepidopteran herbivores. Insect Biochem. Mol. Biol. 42, 174–182. https://doi.org/10.1016/j.ibmb.2011.12.002 (2012).CAS 
Article 
PubMed 

Google Scholar 
Falk, K. L. et al. The role of glucosinolates and the jasmonic acid pathway in resistance of Arabidopsis thaliana against molluscan herbivores. Mol. Ecol. 23, 1188–1203. https://doi.org/10.1111/mec.12610 (2014).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
Gloss, A. D. et al. Evolution in an ancient detoxification pathway is coupled with a transition to herbivory in the Drosophilidae. Mol. Biol. Evol. 31, 2441–3245. https://doi.org/10.1093/molbev/msu201 (2014).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
Bhatt, P., Zhou, X., Huang, Y., Zhang, W. & Chen, S. Characterization of the role of esterases in the biodegradation of organophosphate, carbamate, and pyrethroid pesticides. J. Hazard. Mater. 1, 125026. https://doi.org/10.1016/j.jhazmat.2020.125026 (2021).CAS 
Article 

Google Scholar 
Murfadunnisa, S. et al. Larvicidal and enzyme inhibition of essential oil from Spheranthus amaranthroids (Burm.) against lepidopteran pest Spodoptera litura (Fab.) and their impact on non-target earthworms. Biocatal. Agric. Biotechnol. 21, 101324. https://doi.org/10.1016/j.bcab.2019.101324 (2019).Article 

Google Scholar 
Sengottayan, S. N. Physiological and biochemical effect of neem and other Meliaceae plants secondary metabolites against Lepidopteran insects. Front. Physiol. 4, 359. https://doi.org/10.3389/fphys.2013.00359 (2013).Article 

Google Scholar 
Augustyniak, M., Gladysz, M. & Dziewięcka, M. The Comet assay in insects: Status, prospects and benefits for science. Mutat. Res. Rev. Mutat. Res. 767, 67–76. https://doi.org/10.1016/j.mrrev.2015.09.001Get (2016).CAS 
Article 
PubMed 

Google Scholar 
Foster, E. R. & Downs, J. A. Histone H2A phosphorylation in DNA double strand break repair. FEBS J. 272, 3231–3240. https://doi.org/10.1111/j.1742-4658.2005.04741.x (2005).CAS 
Article 
PubMed 

Google Scholar 
Porichha, S. K., Sarangi, P. K. & Prasad, R. Genotoxic effect of chlorpyrifosin Channa punctatus. Cytol. Genet. 9, 631–638 (1998).
Google Scholar 
Kalita, M. K., Haloi, K. & Devi, D. Larval exposure to chlorpyrifos affects nutritional physiology and induces genotoxicity in silkworm Philosamia ricini (Lepidoptera: Saturniidae). Front. physiol. 7, 1–14. https://doi.org/10.3389/fphys.2016.00535 (2016).Article 

Google Scholar 
Datta, R. et al. Assessment of genotoxic and biochemical effects of purified compounds of Alpinia galanga on a polyphagous lepidopteran pest Spodoptera litura (Fabricius). Phytoparasitica 48, 501–511. https://doi.org/10.1007/s12600-020-00813-8 (2020).CAS 
Article 

Google Scholar 
Afify, A. & Negm, A. A. K. H. Genotoxic effect of insect growth regulators on different stages of peach fruit fly, Bactrocera zonata (Saunders)(Diptera: Tephritidae). Afr. Entomol. 26, 154–161 (2018).Article 

Google Scholar 
Gupta, J. N., Verma, A. N. & Kashyap, R. K. An improved method for mass rearing for melon fruit fly Dacus cucurbitae Coquillett. Indian J. Entomol. 40, 470–471 (1978).
Google Scholar 
Srivastava, B. G. A chemically defined diet for Dacus cucurbitae (Coq.) larvae under aseptic conditions. Entomol. News Lett. 5, 24 (1975).
Google Scholar 
Kumar, A., Sood, S., Mehta, V., Nadda, G. & Shanker, A. Biology of Thysanoplusia orichalcea (Fab.) in relation to host preference and suitability for insect culture and bioefficacy. Indian J. Appl. Entomol. 18, 16–21 (2004).
Google Scholar 
Martinez, S. S. & Emden, H. F. V. Growth disruption, abnormalities and mortality of Spodoptera littoralis (Boisduval) (Lepidoptera: Noctuidae) caused by azadirachtin. Neotrop. Entomol. 30, 113–125 (2001).CAS 
Article 

Google Scholar 
Khan, Z. R. & Saxena, R. C. Behavioural and physiological responses of Sogatella furcifera (Homoptera: Delphacidae) to selected resistant and susceptible rice cultivars. J. Econ. Entomol. 78, 1280–1286 (1985).Article 

Google Scholar 
Zimmer, M. Phenol oxidation. In Methods to Study Litter Decomposition (eds Graça, M. A. et al.) (Springer, 2005).
Google Scholar 
Kono, Y. Generation of superoxide radical during auto-oxidation of hydroxylamine and an assay for superoxide dismutase. Arch. Biochem. Biophys. 186, 189–195. https://doi.org/10.1016/0003-9861(78)90479-4 (1978).CAS 
Article 
PubMed 

Google Scholar 
Bergmeyer, H. U. Reagents for enzymatic analysis. In Methods of Enzymatic Analysis (eds Bergmeyer, H. U. & Gawehn, K.) 438 (Verlag Chemie, 1974).
Google Scholar 
Chien, C. & Dauterman, W. C. Studies on glutathione S-transferases in Helicoverpa (=Heliothis) zea. Insect Biochem. 21, 857–864. https://doi.org/10.1016/0020-1790(91)90092-S (1991).CAS 
Article 

Google Scholar 
Katzenellenbogen, B. & Kafatos, F. C. General esterases of silk worm moth moulting fluid: Preliminary characterization. J. Insect Physiol. 17, 1139–1151. https://doi.org/10.1016/0022-1910(71)90016-3 (1971).CAS 
Article 

Google Scholar 
Mac Intyre, R. J. A method for measuring activities of acid phosphatases separated by acrylamide gel electrophoresis. Biochem. Genet. 5, 45–56 (1971).CAS 
Article 

Google Scholar 
Singh, N. P., McCoy, M. T., Tice, R. R. & Schneider, E. L. A simple technique for quantitation of low levels of DNA damage in individual cells. Exp. Cell Res. 175, 184–191 (1988).CAS 
PubMed 
Article 

Google Scholar  More

Infection properties of clade A and clade B T7-like cyanophagesWe set out to test the hypothesis that the phylogenetic separation of T7-like cyanophages into two major clades reflects differences in their infection physiology. To do this we investigated a suite of infection properties of three pairs of clade A and B phages, each pair infecting the same Synechococcus host (Table 1) to allow us to control for variability in host genetics and physiology. These six cyanophages are representatives of 3 clade A and 2 clade B cyanophage subclades (SI Appendix, Table S1).Table 1 Summary of infection physiology of three pairs of clade A and clade B cyanophages infecting the same Synechococcus hosts.Full size tableWe began by investigating adsorption kinetics and the length of time taken to produce new phages in the infection cycle, the latent period, from phage growth curve experiments. In all three pairs of phages, adsorption was 7–15-fold more rapid in the clade A phage versus the clade B phage (Fig. 1, Table 1). Furthermore, the clade A phage had a faster infection cycle with a latent period that was 3-5-fold shorter than the clade B phage on the same host (Fig. 1a–c) (Table 1). To determine how representative these findings are for a greater diversity of T7-like cyanophages we report the latent period of nine additional non-paired phages that infect a variety of hosts and span the diversity of this cyanophage genus, measured here and taken from the literature (SI Appendix, Table S1). These phages showed the same pattern as observed between phage pairs, although one clade A phage had a relatively long latent period (see SI Appendix, Table S1). Overall, the 5 clade A phages representative of 5 subclades had a significantly shorter latent period (3.3 ± 3.6 h, n = 5 phages (mean ± SD) than the 10 clade B phages from 7 subclades (7.7 ± 2.0 h, n = 10 phages) (Kruskal-Wallis: χ2 = 4.72, df = 1; p = 0.029, n = 15). No significant differences in the length of the latent period were found for clade B phages that infected Synechococcus and Prochlorococcus (Kruskal-Wallis: χ2 = 1.13, df = 1; p = 0.29, n = 10).Fig. 1: Comparison of the infection physiology between pairs of clade A and clade B T7-like cyanophage infecting the same Synechococcus host.a–c Cyanophage growth curves, d–f burst sizes, g–i virulence as the percentage of lysed host cells, j–l decay as loss of infectivity, m–o plaque sizes. a, d, g, j, m Clade A Syn5 phage and clade B S-TIP37 phage infecting WH8109. b, e, h, k, n Clade A S-CBP42 phage and clade B S-RIP2 phage infecting WH7803. c, f, i, l, o Clade A S-TIP28 phage and clade B S-TIP67 phage infecting CC9605. The host strain is shown at the right of the panels. Red and blue lines or bars show results for clade A and clade B phages, respectively. a–c, g–I Error bars indicate standard deviations. d–f Burst size results are for single cells. j–l The solid line shows the fitted multi-level linear model. m–o The time after infection at which plaques were photographed appears above the images. *p value  More