Ecology

Ethics statementAll procedures were conducted in accordance with the current laws in Thailand on experimental animals and were approved by the safety management committee for experiments of the Laboratory Animal Center, Chiang Mai University (Project Number 2561/FA-0001). The study also followed the recommendations in the ARRIVE guidelines.Species-specific primer designAll the DNA tissue analysed originated from the mucus of the individual giant snakehead. Total DNA was extracted from the mucus sample using the Qiagen DNeasy Blood and Tissue Kit (Qiagen, Valencia, CA). Extracted DNA was used as a template for qPCR assay together with synthetic fragments. DNA samples were quantified using a Qubit fluorometer (Life Technologies) calibrated with the Quant-iT dsDNA HS Assay following the manufacturer’s instructions. For each replicate, 3 µL volumes were measured.Species-specific primers and a minor-groove binding (MGB) probe incorporating a 5′ FAM reporter dye and a 3′ non-fluorescent quencher were designed to amplify an 127 bp targeting within the 16S region for the giant snakehead (C. micropeltes), using Primer Express (V3.0, Life Technologies; Table 3). Probe and primer sequences were matched against the National Centre for Biotechnology Information (NCBI, http://www.ncbi.nlm.nih.gov/) nucleotide database with BLASTn (Basic Local Alignment Search Tool) to confirm the species’ specificity for the giant snakehead in silico assays.Table 3 Details of species-specific primers and the probe designed to amplify a 127 bp fragment of the 16S region of Channa micropeltes (Cuvier, 1831).Full size tableTo ensure that the assay only amplified the giant snakehead, it was deployed on a closely related species commonly found in Thai freshwater environments using conventional PCR amplification and visualization on a 1.5% agarose gel stained with SYBR Safe DNA Gel Stain (Life Technologies).qPCR assayThe qPCR assay was deployed using Environmental Master Mix (Applied Biosystems) on mucus samples from the giant snakehead and related species to ensure the species specificity to the qPCR assay. In addition, eDNA qPCR assay for the giant snakehead, a water sample collected from tank at Phayao Freshwater Aquarium (Phayao Inland Fisheries Research and Development Center) was known to have only the giant snakehead was included as a positive control for the presence of amplifiable eDNA in water samples. The tank contains around 4.5 m3 of water with one individual of giant snakehead resides in the tank (the fish is about 60–70 cm in length).All eDNA qPCR amplifications were performed in three replicates in a final volume of 20 µL, using 10.0 µL of 2 × TaqMan Environmental Master Mix 2.0 (Thermo Fisher Scientific), 2.0 µL of DNA template, 900 nM each of the F/R primers, and 125 nM of the probe. Samples were run under the following conditions: an initial 10 min incubation at 95 °C followed by 50 cycles of denaturation at 95 °C for 15 s and annealing/extension at 60 °C for 1 min. Negative controls with all PCR reagents but no template (three replicates) were run in parallel to assess potential contamination. The quantification cycle (Cq) was converted to quantities per unit volume using the linear regression obtained from the synthesized target gene standard curve (Integrated DNA Technologies Pte. Ltd., Singapore). The giant snakehead eDNA concentrations were then reported as copies/mL. The limit of detection (LOD) and the limit of quantification (LOQ) were also measured using the standard dilution series of synthesized target gene fragment with known copy numbers. A dilution series containing 1.5 × 101 to 1.5 × 104 copies per PCR tube were prepared and used as quantification standards. The calculation of LOD and LOQ was done using published R script by Klymus et al.26.DNA extraction from the filtersDNA trapped on the filters obtained from the aquarium experiments and field collections were extracted using Qiagen DNeasy Blood and Tissue Kit (Qiagen, Hilden, Germany) using a protocol modified from the manufacturer’s protocol with the following changes: the DNA from all samples were eluted twice with 25 µL AE buffer, in a total volume of 50 µL to obtain a more concentrated eDNA solution. The volume of ATL buffer (360 µL), Proteinase K (40 µL), AL buffer (400 µL) and Ethanol (400 µL) were doubled.Aquarium experimentAn aquarium experiment was used to test the extent to which qPCR of water samples can detect eDNA of giant snakehead at low simulated densities. The juvenile giant snakehead was obtained from the fish store and transported to a laboratory at Chiang Mai University. The giant snakeheads were then held in separate 120 L plastic holding containers in which the water was continuously filtered. The fish were fed frozen shrimp/commercially available flake fish food three times a week, and were held at 23 ± 1 °C.The sensitivity of eDNA detection in the aquaria was evaluated by conducting three aquarium experiments using plastic tanks (30 × 45 × 25 cm) filled with 120 L of aged-tap water. The water in the tanks was continuously aerated through a filter. In each experiment, the giant snakeheads were randomly assigned to the tanks (10 individuals per tank). The average size of the snakeheads was 9.7 cm (body length ranging from 9.1 to 10.6 cm). The average weight was 8.15 g (ranging from 6.7 to 10.6 g). The water in the tanks was maintained at 23 ± 1 °C. A 300 mL water sample from each tank was collected at each time point (0, 3, 6, 12, 24, 48, 72, 96, 120, 144, and 168 after removal of the fishes from the tanks) in triplicate. Collected water was filtered on a GF/F filter (0.7 μm Whatman International Ltd., Maidstone, UK). The eDNA from each sample solution was extracted using a Qiagen DNeasy Blood and Tissue Kit (Qiagen, Hilden, Germany) in a final volume of 50 µL, detailed in DNA extraction from the filters. To confirm the absence of the giant snakehead eDNA in the water prior to the experiments, three tanks without giant snakehead were prepared and water sample was collected and treated as described above.Real-time PCR was performed with the species-specific primers and probe set using a Rotor-Gene Q system (Qiagen, Hilden, Germany). The reaction conditions were the same as described in qPCR assay. Three replicates were conducted for each sample including the negative PCR control and positive control.eDNA field collectionWater samples were collected at 6 points within Kwan Payao according to the survey locations of the Inland Fisheries Research and Development Center. Additional water samples were collected from 11 and 6 locations in Ing River where water flowed into and out of Kwan Payao, respectively (Fig. 1). To avoid contamination, all field equipment was sterilized using 10% bleach, UV-Crosslinker or autoclaved and sealed prior to transport to the study site, and a separate pair of nitrile disposable gloves were used for each sample. At each site, water samples were collected three replicate in bucket that had been previously decontaminated with a 10% bleach rinse followed by two distilled water rinses.In total, water samples were collected from 6 sites (in Kwan Phayao) and from 17 sites (in the Ing River) from 15th February to 5th March 2019, the middle of the dry season. Each site was sampled in triplicate, 300 mL samples of water were collected and filtered on GF/F filter (0.7 μm Whatman International Ltd., Maidstone, UK). In total, 306 water samples were collected from the surface water of lakes and rivers. For every sampling day, deionised water (300 mL) was filtrated as a negative control. The water samples and real-time PCR was processed as described above in qPCR assay. More

1.Heger, T. & Jeschke, J. M. The enemy release hypothesis as a hierarchy of hypotheses. Oikos 123, 741–750 (2014).Article 

Google Scholar 
2.Elton, C. S. The Ecology of Invasions by Animals and Plants (Springer, 1958).3.Keane, R. M. & Crawley, M. J. Exotic plant invasions and the enemy release hypothesis. Trends Ecol. Evol. 17, 164–170 (2002).Article 

Google Scholar 
4.Mitchell, C. E. & Power, A. G. Release of invasive plants from fungal and viral pathogens. Nature 421, 625–627 (2003).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
5.Colautti, R. I., Ricciardi, A., Grigorovich, I. A. & MacIsaac, H. J. Is invasion success explained by the enemy release hypothesis? Ecol. Lett. 7, 721–733 (2004).Article 

Google Scholar 
6.Liu, H. & Stiling, P. Testing the enemy release hypothesis: a review and meta-analysis. Biol. Invasions 8, 1535–1545 (2006).Article 

Google Scholar 
7.Meijer, K., Schilthuizen, M., Beukeboom, L. & Smit, C. A review and meta-analysis of the enemy release hypothesis in plant-herbivorous insect systems. PeerJ 4, e2560v1 (2016).Article 

Google Scholar 
8.Jeschke, J. M. & Heger, T. (eds) Invasion Biology: Hypotheses and Evidence (CABI, 2018).9.Levine, J. M., Adler, P. B. & Yelenik, S. G. A meta-analysis of biotic resistance to exotic plant invasions. Ecol. Lett. 7, 975–989 (2004).Article 

Google Scholar 
10.Maron, J. L. & Vilà, M. When do herbivores affect plant invasion? Evidence for the natural enemies and biotic resistance hypotheses. Oikos 95, 361–373 (2001).Article 

Google Scholar 
11.Callaway, R. M. & Ridenour, W. M. Novel weapons: invasive success and the evolution of increased competitive ability. Front. Ecol. Environ. 2, 436–443 (2004).Article 

Google Scholar 
12.Cappuccino, N. & Arnason, J. T. Novel chemistry of invasive exotic plants. Biol. Lett. 2, 189–193 (2006).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
13.Bezemer, T. M., Harvey, J. A. & Cronin, J. T. Response of native insect communities to invasive plants. Ann. Rev. Entomol. 59, 119–141 (2014).CAS 
Article 

Google Scholar 
14.Keeler, M. S. & Chew, F. S. Escaping an evolutionary trap: preference and performance of a native insect on an exotic invasive host. Oecologia 156, 559–568 (2008).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
15.Eckberg, J. O., Tenhumberg, B. & Louda, S. M. Insect herbivory and propagule pressure influence Cirsium vulgare invasiveness across the landscape. Ecology 93, 1787–1794 (2012).PubMed 
Article 
PubMed Central 

Google Scholar 
16.Bürki, C. & Nentwig, W. Comparison of herbivore insect communities of Heracleum sphondylium and H. mantegazzianum in Switzerland (Spermatophyta: Apiaceae). Entomol. Gen. 22, 147–155 (1997).Article 

Google Scholar 
17.Cincotta, C. L., Adams, J. M. & Holzapfel, C. Testing the enemy release hypothesis: a comparison of foliar insect herbivory of the exotic Norway maple (Acer platanoides L.) and the native sugar maple (A. saccharum L.). Biol. Invasions 11, 379–388 (2008).Article 

Google Scholar 
18.Cronin, J. T., Bhattarai, G. P., Allen, W. J. & Meyerson, L. A. Biogeography of a plant invasion: plant-herbivore interactions. Ecology 96, 1115–1127 (2015).PubMed 
Article 
PubMed Central 

Google Scholar 
19.Hu, X.-T. & Dong, B.-C. Herbivory and nitrogen availability affect performance of an invader Alternanthera philoxeroides and its native congener A. sessilis. Flora 257, 151412 (2019).Article 

Google Scholar 
20.Agrawal, A. A. & Kotanen, P. M. Herbivores and the success of exotic plants: a phylogenetically controlled experiment. Ecol. Lett. 6, 712–715 (2003).Article 

Google Scholar 
21.Agrawal, A. A. et al. Enemy release? An experiment with congeneric plant pairs and diverse above- and belowground enemies. Ecology 86, 2979–2989 (2005).Article 

Google Scholar 
22.Parker, J. D. & Hay, M. E. Biotic resistance to plant invasions? Native herbivores prefer non-native plants. Ecol. Lett. 8, 959–967 (2005).Article 

Google Scholar 
23.Parker, J. D., Burkepile, D. E. & Hay, M. E. Opposing effects of native and exotic herbivores on plant invasions. Science 311, 1459–1461 (2006).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
24.Parker, I. M. & Gilbert, G. S. When there is no escape: the effects of natural enemies on native, invasive, and non-native plants. Ecology 88, 1210–1224 (2007).PubMed 
Article 
PubMed Central 

Google Scholar 
25.Dostál, P. et al. Enemy damage of exotic plant species is similar to that of natives and increases with productivity. J. Ecol. 101, 388–399 (2013).Article 

Google Scholar 
26.Meijer, K. et al. Phytophagous insects on native and non-native host plants: combining the community approach and the biogeographical approach. PLoS ONE 10, e0125607 (2015).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
27.Schultheis, E. H., Berardi, A. E. & Lau, J. A. No release for the wicked: enemy release is dynamic and not associated with invasiveness. Ecology 96, 2446–2457 (2015).PubMed 
Article 
PubMed Central 

Google Scholar 
28.Beckstead, J. & Parker, I. M. Invasiveness of Ammophila arenaria: release from soil-borne pathogens? Ecology 84, 2824–2831 (2003).Article 

Google Scholar 
29.van Kleunen, M., Weber, E. & Fischer, M. A meta-analysis of trait differences between invasive and non-invasive plant species. Ecol. Lett. 13, 235–245 (2010).PubMed 
Article 
PubMed Central 

Google Scholar 
30.Ashton, I. W. & Lerdau, M. T. Tolerance to herbivory, and not resistance, may explain differential success of invasive, naturalized, and native North American temperate vines. Divers. Distrib. 14, 169–178 (2008).Article 

Google Scholar 
31.Seebens, H. et al. No saturation in the accumulation of alien species worldwide. Nat. Commun. 8, 14435 (2017).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
32.Eppinga, M. B., Rietkerk, M., Dekker, S. C., De Ruiter, P. C. & van der Putten, W. H. Accumulation of local pathogens: a new hypothesis to explain exotic plant invasions. Oikos 114, 168–176 (2006).Article 

Google Scholar 
33.Bufford, J. L. et al. Taxonomic similarity, more than contact opportunity, explains novel plant–pathogen associations between native and alien taxa. N. Phytol. 212, 657–667 (2016).CAS 
Article 

Google Scholar 
34.White, E. M., Wilson, J. C. & Clarke, A. R. Biotic indirect effects: a neglected concept in invasion biology. Divers. Distrib. 12, 443–455 (2006).Article 

Google Scholar 
35.Allen, W. J. in Plant Invasions: The Role of Species Interactions (CABI Publishing, 2020).36.Holt, R. D. Predation, apparent competition, and the structure of prey communities. Theor. Popul. Biol. 12, 197–229 (1977).MathSciNet 
CAS 
PubMed 
Article 

Google Scholar 
37.Holt, R. D. & Bonsall, M. B. Apparent competition. Annu. Rev. Ecol. Evol. Syst. 48, 447–471 (2017).Article 

Google Scholar 
38.Sessions, L. & Kelly, D. Predator-mediated apparent competition between an introduced grass, Agrostis capillaris, and a native fern, Botrychium australe (Ophioglossaceae), in New Zealand. Oikos 96, 102–109 (2002).Article 

Google Scholar 
39.Dangremond, E. M., Pardini, E. A. & Knight, T. M. Apparent competition with an invasive plant hastens the extinction of an endangered lupine. Ecology 91, 2261–2271 (2010).PubMed 
PubMed Central 
Article 

Google Scholar 
40.Bhattarai, G. P., Meyerson, L. A. & Cronin, J. T. Geographic variation in apparent competition between native and invasive Phragmites australis. Ecology 98, 349–358 (2017).PubMed 
Article 
PubMed Central 

Google Scholar 
41.Simberloff, D. & Von Holle, B. Positive interactions of nonindigenous species: invasional meltdown? Biol. Invasions 1, 21–32 (1999).Article 

Google Scholar 
42.Allen, W. J. et al. Community-level direct and indirect impacts of an invasive plant favour exotic over native species. J. Ecol. 108, 2499–2510 (2020).Article 

Google Scholar 
43.Morris, R. J., Lewis, O. T. & Godfray, C. J. Experimental evidence for apparent competition in a tropical forest food web. Nature 428, 310–313 (2004).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
44.Tack, A. J. M., Gripenberg, S. & Roslin, T. Can we predict indirect interactions from quantitative food webs? – an experimental approach. J. Anim. Ecol. 80, 108–118 (2011).PubMed 
Article 

Google Scholar 
45.Frost, C. M. et al. Apparent competition drives community-wide parasitism rates and changes in host abundance across ecosystem boundaries. Nat. Commun. 7, 12644 (2016).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
46.Bardgett, R. D. & Wardle, D. A. Aboveground–Belowground Linkages: Biotic Interactions, Ecosystem Processes and Global Change (Oxford University Press, 2010).47.Heinen, R., Biere, A., Harvey, J. A. & Bezemer, T. M. Effects of soil organisms on aboveground plant-insect interactions in the field: patterns, mechanisms and the role of methodology. Front. Ecol. Evol. 6, 106 (2018).Article 

Google Scholar 
48.Bever, J. D., Westover, K. M. & Antonovics, J. Incorporating the soil community into plant population dynamics: the utility of the feedback approach. J. Ecol. 85, 561–573 (1997).Article 

Google Scholar 
49.Kulmatiski, A., Beard, K. H., Stevens, J. R. & Cobbold, S. M. Plant–soil feedbacks: a meta‐analytical review. Ecol. Lett. 11, 980–992 (2008).PubMed 
Article 

Google Scholar 
50.Levine, J. M., Pachepsky, E., Kendall, B. E., Yelenik, S. G. & Lambers, J. H. Plant-soil feedbacks and invasive spread. Ecol. Lett. 9, 1005–1014 (2006).PubMed 
Article 
PubMed Central 

Google Scholar 
51.Suding, K. N., Harpole, W. S., Fukami, T., Kulmatiski, A., MacDougall, A. S., Stein, C. & van der Putten, W. H. Consequences of plant–soil feedbacks in invasion. J. Ecol. 101, 298–308 (2013).Article 

Google Scholar 
52.Crawford, K. M. et al. When and where plant-soil feedback may promote plant coexistence: a meta-analysis. Ecol. Lett. 22, 1274–1284 (2019).PubMed 
PubMed Central 

Google Scholar 
53.Cornelissen, T., Fernandes, G. W. & Vasconcellos-Neto, J. Size does matter: variation in herbivory between and within plants and the plant vigor hypothesis. Oikos 117, 1121–1130 (2008).Article 

Google Scholar 
54.Price, P. W. The plant vigor hypothesis and herbivore attack. Oikos 62, 244–251 (1991).Article 

Google Scholar 
55.Waller, L. P. et al. Biotic interactions drive ecosystem responses to plant invaders. Science 368, 967–972 (2020).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
56.Kozlov, M. V., Lanta, V., Zverev, V. & Zvereva, E. L. Global patterns in background losses of woody plant foliage to insects. Glob. Ecol. Biogeogr. 24, 1126–1135 (2015).Article 

Google Scholar 
57.Zas, R., Moreira, X. & Sampedro, L. Tolerance and induced resistance in a native and an exotic pine species: relevant traits for invasion ecology. J. Ecol. 99, 1316–1326 (2011).Article 

Google Scholar 
58.Croy, J. R., Meyerson, L. A., Allen, W. J., Bhattarai, G. P. & Cronin, J. T. Lineage and latitudinal variation in Phragmites australis tolerance to herbivory: implications for invasion success. Oikos 129, 1341–1357 (2020).Article 

Google Scholar 
59.Liu, G., Huang, Q.-Q., Lin, Z.-G., Huang, F.-F., Liao, H.-X. & Peng, S.-L. High tolerance to salinity and herbivory stresses may explain the expansion of Ipomoea cairica to salt marshes. PLoS ONE 7, e48829 (2012).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
60.Paynter, Q. et al. Why did specificity testing fail to predict the field host-range of the gorse pod moth in New Zealand. Biol. Control 46, 453–462 (2008).Article 

Google Scholar 
61.Groenteman, R., Fowler, S. V. & Sullivan, J. J. St. John’s wort beetles would not have been introduced to New Zealand now: a retrospective host range test of New Zealand’s most successful weed biocontrol agents. Biol. Control 57, 50–58 (2011).Article 

Google Scholar 
62.Blossey, B. & Nötzold, R. Evolution of increased competitive ability in invasive nonindigenous plants: a hypothesis. J. Ecol. 83, 887–889 (1995).Article 

Google Scholar 
63.Felker-Quinn, E., Schweitzer, J. A. & Bailey, J. K. Meta-analysis reveals evolution in invasive plant species but little support for evolution of increased competitive ability (EICA). Ecol. Evol. 3, 739–751 (2013).PubMed 
PubMed Central 
Article 

Google Scholar 
64.Rotter, M. C. & Holeski, L. M. A meta-analysis of the evolution of increased competitive ability hypothesis: genetic-based trait variation and herbivory resistance trade-offs. Biol. Invasions 20, 2647–2660 (2018).Article 

Google Scholar 
65.Shelby, N. et al. No difference in the competitive ability of introduced and native Trifolium provenances when grown with soil biota from their introduced and native ranges. AoB Plants 8, plw016 (2016).PubMed 
PubMed Central 
Article 

Google Scholar 
66.Harvey, J. A., Bukovinszky, T. & van der Putten, W. H. Interactions between invasive plants and insect herbivores: a plea for a multitrophic perspective. Biol. Conserv. 143, 2251–2259 (2010).Article 

Google Scholar 
67.Allen, W. J. et al. Multitrophic enemy release of invasive Phragmites australis and its introduced herbivores in North America. Biol. Invasions 17, 3419–3432 (2015).Article 

Google Scholar 
68.Kim, T. N. & Underwood, N. Plant neighborhood effects on herbivory: damage is both density and frequency dependent. Ecology 96, 1431–1437 (2015).PubMed 
Article 

Google Scholar 
69.Bartomeus, I., Vilà, M. & Santamaría, L. Contrasting effects of invasive plants in plant-pollinator networks. Oecologia 155, 761–770 (2008).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
70.Lekberg, Y., Gibbons, S. M., Rosendahl, S. & Ramsey, P. W. Severe plant invasions can increase mycorrhizal fungal abundance and diversity. ISME J. 7, 1424–1433 (2013).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
71.Fernandez-Conradi, P., Jactel, H., Robin, C., Tack, A. J. M. & Castagneyrol, B. Fungi reduce preference and performance of insect herbivores on challenged plants. Ecology 99, 300–311 (2018).PubMed 
Article 
PubMed Central 

Google Scholar 
72.Reinhart, K. O. & Callaway, R. M. Soil biota and invasive plants. N. Phytol. 170, 445–457 (2006).Article 

Google Scholar 
73.Gioria, M. & Osborne, B. A. Resource competition in plant invasions: emerging patterns and research needs. Front. Plant Sci. 5, 501 (2014).PubMed 
PubMed Central 
Article 

Google Scholar 
74.Frost, C. M. et al. Using ecological network theory to predict biological invasions. Trends Ecol. Evol. 34, 831–843 (2019).PubMed 
Article 
PubMed Central 

Google Scholar 
75.Sauve, A. M. C., Thébault, E., Pocock, M. J. O. & Fontaine, C. How plants connect pollination and herbivory networks and their contribution to community stability. Ecology 97, 908–917 (2016).PubMed 
PubMed Central 

Google Scholar 
76.Pilosof, S., Porter, M. A., Pascual, M. & Kéfi, S. The multilayer nature of ecological networks. Nat. Ecol. Evol. 1, 0101 (2017).Article 

Google Scholar 
77.Weir, B. S., Turner, S. J., Silvester, W. B., Park, D. C. & Young, J. M. Unexpectedly diverse Mesorhizobium strains and Rhizobium leguminosarum nodulate native legume genera of New Zealand, while introduced legume weeds are nodulated by Bradyrhizobium species. Appl. Environ. Microbiol. 70, 5980–5987 (2004).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
78.Cappuccino, N. & Carpenter, D. Invasive exotic plants suffer less herbivory than non-invasive exotic plants. Biol. Lett. 1, 435–438 (2005).PubMed 
PubMed Central 
Article 

Google Scholar 
79.Liu, H., Stiling, P. & Pemberton, R. W. Does enemy release matter for invasive plants? evidence from a comparison of insect herbivore damage among invasive, non-invasive and native congeners. Biol. Invasions 9, 773–781 (2007).Article 

Google Scholar 
80.Howell, C. Consolidated List of Environmental Weeds in New Zealand. DOC Research & Development Series 292 (Department of Conservation, 2008).81.Ghanizadeh, H. & Harrington, K. C. Weed management in New Zealand pastures. Agronomy 9, 448 (2019).CAS 
Article 

Google Scholar 
82.Kos, M., Tuijl, M. A. B., de Roo, J., Mulder, P. P. J. & Bezemer, T. M. Species-specific plant–soil feedback effects on aboveground plant-insect interactions. J. Ecol. 103, 904–914 (2015).CAS 
Article 

Google Scholar 
83.Heinen, R., Biere, A. & Bezemer, T. M. Plant traits shape soil legacy effects on individual plant–insect interactions. Oikos 129, 261–273 (2020).CAS 
Article 

Google Scholar 
84.Bezemer, T. M et al. Above‐and below‐ground herbivory effects on below‐ground plant–fungus interactions and plant–soil feedback responses. J. Ecol. 101, 325–333 (2013).Article 
CAS 

Google Scholar 
85.Heinze, J., Wacker, A. & Kulmatiski, A. Plant–soil feedback effects altered by aboveground herbivory explain plant species abundance in the landscape. Ecology 101, e03023 (2020).PubMed 
Article 
PubMed Central 

Google Scholar 
86.Müller, C. B., Adriaanse, I. C. T., Belshaw, R. & Godfray, H. C. J. The structure of an aphid-parasitoid community. J. Anim. Ecol. 68, 346–370 (1999).Article 

Google Scholar 
87.R Core Team. R: a language and environment for statistical computing. Version 3.6.1. R Foundation for Statistical Computing http://www.R-project.org (2019).88.Bates, D. et al. lme4: linear mixed-effects models using ‘Eigen’ and S4. R package version 1.1-21 http://CRAN.R-project.org/package=lme4 (2019).89.Lenth, R., Singmann, H., Love, J., Buerkner, P. & Herve, M. emmeans: estimated marginal means, aka least-squares means. R package version 1.3.5.1 http://CRAN.R-project.org/package=emmeans (2019).90.Dormann, C. F., Fruend, J. & Gruber, B. bipartite: visualising bipartite networks and calculating some (ecological) indices. R package version 2.13 http://CRAN.R-project.org/package=bipartite (2019). More

1.Gómez, C., White, J. C. & Wulder, M. A. Optical remotely sensed time series data for land cover classification: A review. ISPRS J. Photogramm. Remote. Sens. 116, 55–72 (2016).ADS 
Article 

Google Scholar 
2.Dale, V. H. The relationship between land-use change and climate change. Ecol. Appl. 7, 753–769 (1997).Article 

Google Scholar 
3.Dai, A. Increasing drought under global warming in observations and models. Nat. Clim. Chang. 3, 52–58 (2013).ADS 
Article 

Google Scholar 
4.Houghton, R. A. et al. Carbon emissions from land use and land-cover change. Biogeosciences 9, 5125–5142 (2012).ADS 
CAS 
Article 

Google Scholar 
5.Lambin, E. F. et al. The causes of land-use and land-cover change: Moving beyond the myths. Glob. Environ. Chang. 11, 261–269 (2001).Article 

Google Scholar 
6.Song, X. P. et al. Global land change from 1982 to 2016. Nature 560, 639–643 (2018).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
7.Grimm, N. B. et al. Global change and the ecology of cities. Science 319, 756–760 (2008).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
8.Utkin, L. V. et al. A weighted random survival forest. Knowl. Based Syst. 177, 136–144 (2019).Article 

Google Scholar 
9.Gašparović, M., Zrinjski, M. & Gudelj, M. Automatic cost-effective method for land cover classification (ALCC). Comput. Environ. Urban Syst. 76, 1–10 (2019).Article 

Google Scholar 
10.Hu, Y., Dong, Y. & Batunacun. ,. An automatic approach for land-change detection and land updates based on integrated NDVI timing analysis and the CVAPS method with GEE support. ISPRS J. Photogramm. Remote Sens. 146, 347–359 (2018).ADS 
Article 

Google Scholar 
11.ESA. Sentinel-2 Mission Requirements Document. Earth (2007).12.Main-Knorn, M. et al. Sen2Cor for Sentinel-2. In 3 (2017). https://doi.org/10.1117/12.2278218.13.Gong, P. et al. Stable classification with limited sample: Transferring a 30-m resolution sample set collected in 2015 to mapping 10-m resolution global land cover in 2017. Sci. Bull. 64, 370–373 (2019).Article 

Google Scholar 
14.Homer, C., Huang, C., Yang, L., Wylie, B. & Coan, M. Development of a 2001 National Land-Cover Database for the United States. Photogramm. Eng. Remote. Sens. 70, 829–840 (2004).Article 

Google Scholar 
15.Fry, J. A. et al. Completion of the 2006 national land cover database for the conterminous united states. Photogramm. Eng. Remote Sens. 77, 858–864 (2011).
Google Scholar 
16.Homer, C. et al. Completion of the 2011 national land cover database for the conterminous United States—Representing a decade of land cover change information. Photogramm. Eng. Remote Sens. 81, 345–354 (2015).
Google Scholar 
17.Yang, L. et al. A new generation of the United States National Land Cover Database: Requirements, research priorities, design, and implementation strategies. ISPRS J. Photogramm. Remote Sens. 146, 108–123 (2018).ADS 
Article 

Google Scholar 
18.Chen, J. et al. Global land cover mapping at 30 m resolution: A POK-based operational approach. ISPRS J. Photogramm. Remote Sens. 103, 7–27 (2015).ADS 
Article 

Google Scholar 
19.Hoang, T. T., Truong, V. T., Hayashi, M., Tadono, T. & Nasahara, K. N. New JAXA high-resolution land use/land cover map for Vietnam aiming for natural forest and plantation forest monitoring. Remote Sens. 12, 2707 (2020).20.Phan, D. C., Trung, T. H., Nasahara, K. N. & Tadono, T. JAXA high-resolution land use/land cover map for Central Vietnam in 2007 and 2017. Remote Sens. 10, 1406 (2018).ADS 
Article 

Google Scholar 
21.Nemani, R. Nasa earth exchange: Next generation earth science collaborative. ISPRS Int. Arch. Photogramm. Remote Sens. Spat. Inf. Sci. XXXVIII-8/, 17–17 (2012).22.Mugiraneza, T., Nascetti, A. & Ban, Y. Continuous monitoring of urban land cover change trajectories with landsat time series and landtrendr-google earth engine cloud computing. Remote Sens. 12, 2883 (2020).23.Jin, Y., Liu, X., Yao, J., Zhang, X. & Zhang, H. Mapping the annual dynamics of cultivated land in typical area of the Middle-lower Yangtze plain using long time-series of Landsat images based on Google Earth Engine. Int. J. Remote Sens. 41, 1625–1644 (2020).ADS 
Article 

Google Scholar 
24.Hu, Y. & Hu, Y. Detecting forest disturbance and recovery in Primorsky Krai, Russia, using annual landsat time series and multi-source land cover products. Remote Sens. 12, 129 (2020).Article 

Google Scholar 
25.Miettinen, J., Shi, C. & Liew, S. C. 2015 Land cover map of Southeast Asia at 250 m spatial resolution. Remote Sens. Lett. 7, 701–710 (2016).Article 

Google Scholar 
26.Ghorbanian, A. et al. Improved land cover map of Iran using Sentinel imagery within Google Earth Engine and a novel automatic workflow for land cover classification using migrated training samples. ISPRS J. Photogramm. Remote Sens. 167, 276–288 (2020).ADS 
Article 

Google Scholar 
27.Huang, H. et al. The migration of training samples towards dynamic global land cover mapping. ISPRS J. Photogramm. Remote Sens. 161, 27–36 (2020).ADS 
Article 

Google Scholar 
28.Bagan, H. & Yamagata, Y. Land-cover change analysis in 50 global cities by using a combination of Landsat data and analysis of grid cells. Environ. Res. Lett. 9, 064015 (2014).29.Radoux, J. et al. Automated training sample extraction for global land cover mapping. Remote Sens. 6, 3965–3987 (2014).ADS 
Article 

Google Scholar 
30.Tran, H., Tran, T. & Kervyn, M. Dynamics of land cover/land use changes in the Mekong Delta, 1973–2011: A Remote sensing analysis of the Tran Van Thoi District, Ca Mau Province, Vietnam. Remote Sens. 7, 2899–2925 (2015).ADS 
Article 

Google Scholar 
31.Chi, V. K. et al. Land transitions in Northwest Vietnam: An integrated analysis of biophysical and socio-cultural factors. Hum. Ecol. 41, 37–50 (2013).Article 

Google Scholar 
32.Villamor, G. B., Catacutan, D. C., Truong, V. A. T. & Thi, L. D. Tree-cover transition in Northern Vietnam from a gender-specific land-use preferences perspective. Land Use Policy 61, 53–62 (2017).Article 

Google Scholar 
33.Truong, V. T. et al. JAXA annual forest cover maps for Vietnam during 2015–2018 Using ALOS-2/PALSAR-2 and auxiliary data. Remote Sens. 11, 2412 (2019).ADS 
Article 

Google Scholar 
34.And, R. D. of A. Vietnam’S Modified Submission on Refreence Levels for Redd+ Results Based Payments Under Unfccc. https://redd.unfccc.int/files/vietnam_frl_modified__submission_final_for_posting.pdf (2016).35.Xu, X., Jain, A. K. & Calvin, K. V. Quantifying the biophysical and socioeconomic drivers of changes in forest and agricultural land in South and Southeast Asia. Glob. Chang. Biol. 25, 2137–2151 (2019).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
36.Hansen, M. C. et al. High-resolution global maps of 21st-century forest cover change. Science (80- ) 342, 850–853 (2013).ADS 
CAS 
Article 

Google Scholar 
37.Di Gregorio, A., and Jansen, L. J. M. Land Cover Classification System (LCCS): Classification Concepts and User Manual. Fao http://www.fao.org/3/x0596e/x0596e00.htm (2000).38.Van Thinh, T., Phan, D. C., Nasahara, K. N. & Tadono, T. How does land use/land cover map’s accuracy depend on number of classification classes? Sci. Online Lett. Atmos. 15, 28–31 (2019).
Google Scholar 
39.Talukdar, S. et al. Land-use land-cover classification by machine learning classifiers for satellite observations—a review. Remote Sens. 12, 1135 (2020).ADS 
Article 

Google Scholar 
40.Song, S., Gong, W., Zhu, B. & Huang, X. Wavelength selection and spectral discrimination for paddy rice, with laboratory measurements of hyperspectral leaf reflectance. ISPRS J. Photogramm. Remote Sens. 66, 672–682 (2011).ADS 
Article 

Google Scholar 
41.Szantoi, Z., Smith, S. E., Strona, G., Koh, L. P. & Wich, S. A. Mapping orangutan habitat and agricultural areas using Landsat OLI imagery augmented with unmanned aircraft system aerial photography. Int. J. Remote Sens. 38, 2231–2245 (2017).ADS 
Article 

Google Scholar 
42.Hill, R. A. & Thomson, A. G. Mapping woodland species composition and structure using airborne spectral and LiDAR data. Int. J. Remote Sens. 26, 3763–3779 (2005).ADS 
Article 

Google Scholar 
43.Kontgis, C., Schneider, A. & Ozdogan, M. Mapping rice paddy extent and intensification in the Vietnamese Mekong River Delta with dense time stacks of Landsat data. Remote Sens. Environ. 169, 255–269 (2015).ADS 
Article 

Google Scholar 
44.Kontgis, C. et al. Monitoring peri-urbanization in the greater Ho Chi Minh City metropolitan area. Appl. Geogr. 53, 377–388 (2014).Article 

Google Scholar 
45.D’Amour, C. B. et al. Future urban land expansion and implications for global croplands. Proc. Natl. Acad. Sci. USA 114, 8939–8944 (2017).Article 
CAS 

Google Scholar 
46.Ha, T. V., Tuohy, M., Irwin, M. & Tuan, P. V. Monitoring and mapping rural urbanization and land use changes using Landsat data in the northeast subtropical region of Vietnam. Egypt. J. Remote Sens. Sp. Sci. 23, 11–19 (2020).
Google Scholar 
47.Nhan, T. Q., Van Ly, L. T. & Tan, L. V. How much do rice farmers earn from their crops? Evidence from a rice-exporting country. J. Agric. Stud. 8, 302 (2020).
Google Scholar 
48.Truong, T. D. & Do, L. H. Mangrove forests and aquaculture in the Mekong river delta. Land Use Policy 73, 20–28 (2018).Article 

Google Scholar 
49.Lam-Dao, N., Pham-Bach, V., Nguyen-Thanh, M., Pham-Thi, M.-T. & Hoang-Phi, P. Change detection of land use and riverbank in Mekong Delta, Vietnam using time series remotely sensed data. J. Resour. Ecol. 2, 370–374 (2011).
Google Scholar 
50.Ha, T. T. P., van Dijk, H. & Visser, L. Impacts of changes in mangrove forest management practices on forest accessibility and livelihood: A case study in mangrove-shrimp farming system in Ca Mau Province, Mekong Delta, Vietnam. Land Use Policy 36, 89–101 (2014).Article 

Google Scholar 
51.Le, T. N., Bregt, A. K., van Halsema, G. E., Hellegers, P. J. G. J. & Nguyen, L. D. Interplay between land-use dynamics and changes in hydrological regime in the Vietnamese Mekong Delta. Land Use Policy 73, 269–280 (2018).Article 

Google Scholar 
52.Khuc, Q. V., Tran, B. Q., Meyfroidt, P. & Paschke, M. W. Drivers of deforestation and forest degradation in Vietnam: An exploratory analysis at the national level. For. Policy Econ. 90, 128–141 (2018).Article 

Google Scholar 
53.Meyfroidt, P., Vu, T. P. & Hoang, V. A. Trajectories of deforestation, coffee expansion and displacement of shifting cultivation in the Central Highlands of Vietnam. Glob. Environ. Chang. 23, 1187–1198 (2013).Article 

Google Scholar 
54.Ngo-Duc, T., Kieu, C., Thatcher, M., Nguyen-Le, D. & Phan-Van, T. Climate projections for Vietnam based on regional climate models. Clim. Res. 60, 199–213 (2014).Article 

Google Scholar 
55.Lindesay, J. A. et al. International geosphere-biosphere programme/international global atmospheric chemistry SAFARI-92 field experiment: Background and overview. J. Geophys. Res. Atmos. 101, 23521–23530 (1996).ADS 
CAS 
Article 

Google Scholar 
56.Klemas, V. V., Dobson, J. E., Ferguson, R. L. & Haddad, K. D. A coastal land cover classification system for the NOAA Coastwatch Change Analysis Project. J. Coast. Res. 9, 862–872 (1993).
Google Scholar 
57.Saah, D. et al. Primitives as building blocks for constructing land cover maps. Int. J. Appl. Earth Obs. Geoinf. 85, 101979 (2020).Article 

Google Scholar 
58.Keys, R. G. Cubic convolution interpolation for digital image processing. IEEE Trans. Acoust. 29, 1153–1160 (1981).MathSciNet 
MATH 
Article 

Google Scholar 
59.Farr, T. G. et al. The shuttle radar topography mission. Rev. Geophys. 45, RG2004 (2007).60.Filipponi, F. Sentinel-1 GRD preprocessing workflow. Proceedings 18, 11 (2019).Article 

Google Scholar 
61.Soenen, S. A., Peddle, D. R. & Coburn, C. A. SCS+C: A modified sun-canopy-sensor topographic correction in forested terrain. IEEE Trans. Geosci. Remote Sens. 43, 2148–2159 (2005).ADS 
Article 

Google Scholar 
62.Saleous, N. & Kutler, J. LEDAPS Calibration, Reflectance, Atmospheric Correction Preprocessing Code, Version 2 Summary: Data Citation: Model Product Description: Oak Ridge National Laboratory Distributed Active Archive Center http://daac.ornl.gov/cgi-bin/dsviewer.pl?ds_id=1146 (2012) 10.3334/ORNLDAAC/1146.63.Vermote, E., Justice, C., Claverie, M. & Franch, B. Preliminary analysis of the performance of the Landsat 8/OLI land surface reflectance product. Remote Sens. Environ. 185, 46–56 (2016).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
64.Foga, S. et al. Cloud detection algorithm comparison and validation for operational Landsat data products. Remote Sens. Environ. 194, 379–390 (2017).ADS 
Article 

Google Scholar 
65.Louis, J. et al. Sentinel-2 Global Surface Reflectance Level-2a Product Generated with Sen2Cor. In International Geoscience and Remote Sensing Symposium (IGARSS) 8522–8525 (2019). https://doi.org/10.1109/IGARSS.2019.8898540.66.Roy, D. P. et al. A general method to normalize Landsat reflectance data to nadir BRDF adjusted reflectance. Remote Sens. Environ. 176, 255–271 (2016).ADS 
Article 

Google Scholar 
67.Roy, D. P. et al. Examination of Sentinel-2A multi-spectral instrument (MSI) reflectance anisotropy and the suitability of a general method to normalize MSI reflectance to nadir BRDF adjusted reflectance. Remote Sens. Environ. 199, 25–38 (2017).ADS 
Article 

Google Scholar 
68.Lee, J. S., Ainsworth, T. L., Wang, Y. & Chen, K. S. Polarimetric SAR speckle filtering and the extended sigma filter. IEEE Trans. Geosci. Remote Sens. 53, 1150–1160 (2015).ADS 
Article 

Google Scholar 
69.Roy, D. P. et al. Characterization of Landsat-7 to Landsat-8 reflective wavelength and normalized difference vegetation index continuity. Remote Sens. Environ. 185, 57–70 (2016).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
70.Claverie, M. et al. The Harmonized Landsat and Sentinel-2 surface reflectance data set. Remote Sens. Environ. 219, 145–161 (2018).ADS 
Article 

Google Scholar 
71.Flood, N. Seasonal composite landsat TM/ETM+ Images using the medoid (a multi-dimensional median). Remote Sens. 5, 6481–6500 (2013).ADS 
Article 

Google Scholar 
72.Kaufman, Y. J. & Tanré, D. Atmospherically Resistant Vegetation Index (ARVI) for EOS-MODIS. IEEE Trans. Geosci. Remote Sens. 30, 261–270 (1992).ADS 
Article 

Google Scholar 
73.Saah, D. et al. Land cover mapping in data scarce environments: Challenges and opportunities. Front. Environ. Sci. 7, 150 (2019).Article 

Google Scholar 
74.Li, C., Wang, J., Wang, L., Hu, L. & Gong, P. Comparison of classification algorithms and training sample sizes in urban land classification with landsat thematic mapper imagery. Remote Sens. 6, 964–983 (2014).ADS 
Article 

Google Scholar 
75.Richards, J. A. Remote sensing digital image analysis: An introduction. Remote Sensing Digital Image Analysis: An Introduction vol. 9783642300 (2013).76.Kruse, F. A. et al. The spectral image processing system (SIPS)-interactive visualization and analysis of imaging spectrometer data. Remote Sens. Environ. 44, 145–163 (1993).ADS 
Article 

Google Scholar 
77.Breiman, L. Random forests. Mach. Learn. 45, 5–32 (2001).MATH 
Article 

Google Scholar 
78.He, Y., Lee, E. & Warner, T. A. A time series of annual land use and land cover maps of China from 1982 to 2013 generated using AVHRR GIMMS NDVI3g data. Remote Sens. Environ. 199, 201–217 (2017).ADS 
Article 

Google Scholar 
79.Zhao, H. & Chen, X. Use of normalized difference bareness index in quickly mapping bare areas from TM/ETM+. In International Geoscience and Remote Sensing Symposium (IGARSS) vol. 3 1666–1668 (IEEE, 2005).80.Baloloy, A. B., Blanco, A. C., Raymund Rhommel, R. R. C. & Nadaoka, K. Development and application of a new mangrove vegetation index (MVI) for rapid and accurate mangrove mapping. ISPRS J. Photogramm. Remote Sens. 166, 95–117 (2020).ADS 
Article 

Google Scholar 
81.García, M. J. L. & Caselles, V. Mapping burns and natural reforestation using thematic mapper data. Geocarto Int. 6, 31–37 (1991).Article 

Google Scholar 
82.Wright, C. & Gallant, A. Improved wetland remote sensing in Yellowstone National Park using classification trees to combine TM imagery and ancillary environmental data. Remote Sens. Environ. 107, 582–605 (2007).ADS 
Article 

Google Scholar 
83.Hutchinson, C. F. Techniques for combining Landsat and ancillary data for digital classification improvement. Photogramm. Eng. Remote Sens. 48, 123–130 (1982).
Google Scholar 
84.Tadono, T. et al. Precise global DEM generation by ALOS PRISM. ISPRS Ann. Photogramm. Remote Sens. Spat. Inf. Sci. II–4, 71–76 (2014).Article 

Google Scholar 
85.Jokar Arsanjani, J., Zipf, A., Mooney, P. & Helbich, M. OpenStreetMap in GIScience. OpenStreetMap in GIScience: Experiences, Research, Applications (2015). https://doi.org/10.1007/978-3-319-14280-7.86.Open Development Mekong. OpenDevelopmentMekong. East-West Management Institute https://data.opendevelopmentmekong.net/organization/6f37a27d-2790-4b9a-8570-a36cb1d8108f?res_format=KML (2015).87.Truong, V. T. et al. JAXA annual forest cover maps for Vietnam during 2015–2018 Using ALOS-2/PALSAR-2 and auxiliary data. Remote Sens. 11 (2019).88.Olofsson, P. et al. Good practices for estimating area and assessing accuracy of land change. Remote Sens. Environ. 148, 42–57 (2014).ADS 
Article 

Google Scholar 
89.Congalton, R. G. A review of assessing the accuracy of classifications of remotely sensed data. Remote Sens. Environ. 37, 35–46 (1991).ADS 
Article 

Google Scholar 
90.Schmidt, M. The Sankey diagram in energy and material flow management—part II: Methodology and current applications. J. Ind. Ecol. 12, 173–185 (2008).Article 

Google Scholar 
91.Zha, Y., Gao, J. & Ni, S. Use of normalized difference built-up index in automatically mapping urban areas from TM imagery. Int. J. Remote Sens. 24, 583–594 (2003).ADS 
Article 

Google Scholar 
92.Lacaux, J. P., Tourre, Y. M., Vignolles, C., Ndione, J. A. & Lafaye, M. Classification of ponds from high-spatial resolution remote sensing: Application to Rift Valley Fever epidemics in Senegal. Remote Sens. Environ. 106, 66–74 (2007).ADS 
Article 

Google Scholar 
93.Tucker, C. J. Red and photographic infrared linear combinations for monitoring vegetation. Remote Sens. Environ. 8, 127–150 (1979).ADS 
Article 

Google Scholar 
94.McFeeters, S. K. The use of the Normalized Difference Water Index (NDWI) in the delineation of open water features. Int. J. Remote Sens. 17, 1425–1432 (1996).ADS 
Article 

Google Scholar 
95.Huete, A. R. A soil-adjusted vegetation index (SAVI). Remote Sens. Environ. 25, 295–309 (1988).ADS 
Article 

Google Scholar 
96.Penuelas, J., Baret, F. & Filella, I. Semi-empirical indices to assess carotenoids/chlorophyll a ratio from leaf spectral reflectance. Photosynthetica 31, 221–230 (1995).CAS 

Google Scholar 
97.Jayamanna, S., Kawamura, M. & Tsujiko, Y. Relation between social and environmental conditions in colombo, sri lanka and the urban index estimated by satellite remote sensing data. Int. Arch. Photogram. Remote Sens. 31, 321–326 (1996).
Google Scholar 
98.Shen, L. & Li, C. Water body extraction from Landsat ETM+ imagery using adaboost algorithm. In 2010 18th International Conference on Geoinformatics, Geoinformatics 2010 (2010). https://doi.org/10.1109/GEOINFORMATICS.2010.5567762.99.Perry, C. R. & Lautenschlager, L. F. Functional equivalence of spectral vegetation indices. Remote Sens. Environ. 14, 169–182 (1984).ADS 
Article 

Google Scholar 
100.As-syakur, A. R., Adnyana, I. W. S., Arthana, I. W. & Nuarsa, I. W. Enhanced built-UP and bareness index (EBBI) for mapping built-UP and bare land in an urban area. Remote Sens. 4, 2957–2970 (2012).ADS 
Article 

Google Scholar 
101.Liu, H. Q. & Huete, A. Feedback based modification of the NDVI to minimize canopy background and atmospheric noise. IEEE Trans. Geosci. Remote Sens. 33, 457–465 (1995).ADS 
Article 

Google Scholar 
102.Gitelson, A. A., Viña, A., Ciganda, V., Rundquist, D. C. & Arkebauer, T. J. Remote estimation of canopy chlorophyll content in crops. Geophys. Res. Lett. 32, 1–4 (2005).Article 
CAS 

Google Scholar  More

Our monetary compensation programme for poor rural communities in Mexico’s Sierra de Manantlán Biosphere Reserve encouraged them to forgo exploitation of their natural resources to provide ecosystem services for the city of Colima (see Nature 591, 178; 2021). But focusing solely on the monetary value of ecosystems isn’t enough.The National Forestry Commission of Mexico and the Fund for Natural Protected Areas have made compensatory payments of almost US$1 million to communities in the Cerro Grande region of the reserve since 2003, backed since 2013 by a local fiduciary fund of voluntary contributions from Colima’s citizens. However, the voluntary payments are minimal because most people don’t acknowledge the importance of the mountain forest that provides 90% of their water.As well as monetary schemes, the cutural traditions and the land-tenure rights of communal organizations must be recognized. They should be empowered to draw up contracts between owners of natural resources and urban beneficiaries that will promote their common social, economic and livelihood interests. To increase productivity and family income for impoverished small landowners, payment for ecosystem services could be implemented by using diverse marketing approaches that include sustainable agroforestry and livestock production. More

Insect culturePhyllotreta armoraciae beetles were reared on potted Brassica juncea cultivar “Bau Sin” plants or on potted Brassica rapa cultivar “Yo Tsai Sum” plants (Known-You Seed Co., Ltd) in mesh cages (Bugdorm, MegaView Science Co., Ltd) in a controlled environment chamber at 24 °C, 60% relative humidity, and a 16-h photoperiod. Food plants were cultivated in a growth chamber at 24 °C, 55% relative humidity, and a 14-h photoperiod. Beetles were provided with 3–4-week-old plants once per week, and plants with eggs were kept separately for larval development. Larvae were allowed to pupate in the soil, and after 3 weeks, the soil with pupae was transferred into plastic boxes (9 L, Lock&Lock) until the new generation of beetles emerged.RNA isolation, RNA sequencing, and de novo transcriptome assemblyTotal RNA was extracted from dissected foregut, midgut, hindgut, and Malpighian tubule tissue of newly emerged P. armoraciae beetles that were reared on B. juncea using the innuPREP RNA Mini Kit (Analytik Jena). Tissues from at least ten beetles were pooled per sample. RNA integrity was verified using an Agilent Technologies 2100 Bioanalyzer with the RNA 6000 Nano Kit (Agilent Technologies). RNA quantity was determined using a Nanodrop ND-1000 spectrophotometer (PEQlab Biotechnologie GmbH). One set of RNA samples was sequenced by GATC Biotech on the HiSeq 2500 System from Illumina in Rapid Run mode, using the paired-end (2 × 125 bp) read technology at a depth of 15–25 million reads for each sample. For a second set of samples consisting of four biological replicates per tissue, we additionally performed an on-column DNA digestion with the innuPREP DNase I Digest Kit (Analytik Jena) according to the manufacturer’s instructions. RNA samples were poly(A)-enriched, fragmented, and sequenced at the Max Planck Genome Centre Cologne on the HiSeq 3000 Sequencing System from Illumina, using the paired-end (2 × 150 bp) read technology at a depth of 22 million reads for each sample. Sequencing reads were filtered to remove bad-quality reads based on fastq file scores and trimmed based on read length using CLC Genomics Workbench software version 10.1. With a randomly sampled set of 420 million reads from the two sets of sequencing data, a transcriptome was assembled de novo with the following parameters: nucleotide mismatch cost = 1; insertion = deletion costs = 2; length fraction = 0.6; similarity = 0.9. Conflicts among the individual bases were resolved by voting for the base with the highest frequency. After removing contigs shorter than 250 bp, the final assembly contained 36,445 contigs with an N50 contig size of 2115 bp.Identification of coleopteran MFS transportersWe predicted a protein dataset for P. armoraciae by translating each contig of the gut and Malpighian tubule-specific transcriptome into all six reading frames. After removing sequences shorter than 50 amino acids, we submitted the protein dataset (267,568 sequences) to the TransAAP hosted at the TransportDB 2.0 web portal35. This initial annotation predicted a total of 1401 putative transporter sequences and revealed the MFS and the ABC transporters to be the largest transporter families (Supplementary Data 1). We focused on MFS transporters, which are classified into >80 families53. We used one protein sequence from each MFS family as a query to search candidate MFS transporters in the protein dataset from P. armoraciae using Blastp (E value threshold of 10−5), and assigned each candidate to an MFS family based on sequence similarity to transporter sequences deposited in TCDB. Additional candidates were identified by repeating the search procedure with an extended dataset including the candidate MFS transporters from P. armoraciae. The number of TMDs for each candidate was predicted using the TMHMM Server v.2.054. Partial sequences encoding less than six predicted TMDs were removed from the dataset. The same strategy was used to identify putative MFS transporters in protein datasets that were predicted from the genomes of Leptinotarsa decemlineata, genome annotations v.0.5.355, Anoplophora glabripennis, assembly Agla_1.056, and Tribolium castaneum, assembly Tcas3.057, respectively. The predicted protein sequences are provided in Supplementary Data 2.Digital gene expression analysisDigital gene expression analysis of putative MFS transporters identified in the P. armoraciae transcriptome was carried out using CLC Genomics Workbench v.10.1 by mapping the Illumina reads from the second set of samples onto the reference transcriptome, and counting the reads to estimate gene expression levels. For the cloned MFS genes, complete open-reading frames (ORFs) were used as reference sequences for mapping. For read alignment, we used the following parameters: nucleotide mismatch cost = 2; insertion = deletion costs = 3; length fraction = 0.6; similarity fraction = 0.9; maximum number of hits for a read = 15. Each pair of reads was counted as two. Biases in the sequence datasets and different transcript sizes were corrected using the TPM (transcripts per kilobase million) normalization method to obtain correct estimates for relative expression levels between samples.Phylogenetic analyses of coleopteran MFS transportersWe inferred the lineage-specific diversification patterns of putative MFS transporters from P. armoraciae, L. decemlineata, A. glabripennis, and T. castaneum in phylogenetic analyses with two different datasets, one containing all identified MFS transporters (867 sequences), the other containing a subset of putative sugar porters (120 sequences) from the above four species and 35 sugar porters from C. populi52. The corresponding protein sequences were aligned using the MUSCLE algorithm58 implemented in MEGA 7 with default parameters. The alignments were trimmed manually and the best substitution models were determined using ProtTest 3.4.259. Maximum-likelihood phylogenetic trees were constructed in IQ-TREE version 1.6.060 using the VT + G + F substitution model with 1000 ultrafast bootstrap replicates for the full dataset, and the LG + G + F substitution model with 1000 bootstrap replicates for the subset of putative sugar porters.Identification and sequencing of candidate transportersBased on our phylogenetic analysis, we selected the largest clade of putative MFS transporters that was specifically expanded in P. armoraciae for further studies. Transcriptome analyses revealed the presence of a pseudogene (PaMFS9-ps) that shares between 43 and 95% nucleotide sequence identity with members of the focal clade. The protein encoded by this pseudogene is predicted to possess only two transmembrane domains due to a premature stop codon caused by frameshift mutations in the coding sequence. To obtain the full-length ORFs of partial transcripts, we synthesized 5′- and 3′-rapid amplification of cDNA (complementary DNA) ends (RACE)–cDNA using the SMARTerRACE cDNA Amplification Kit (Clontech) and performed 5′- and 3′-RACE-PCR according to the manufacturer’s protocols (Clontech). All full-length ORFs were cloned into the pCR™4-TOPO® TA vector (Thermo Fisher Scientific) for sequence verification.Tissue-specific expression of candidate transportersWe used quantitative PCR (qPCR) to analyze the expression of candidate transporter genes in the foregut, midgut, hindgut, Malpighian tubules, and other tissues of 1-day-old adult P. armoraciae beetles (n = 4 biological replicates, each with two technical replicates), respectively. In addition, the expression of candidate transporter genes was analyzed in the proximal, central, and distal Malpighian tubule regions (Supplementary Fig. 2d) dissected from 4-day-old adult P. armoraciae beetles (n = 3 biological replicates, each with two technical replicates). Total RNA was extracted using the InnuPrep RNA Mini Kit (Analytik Jena), treated with TURBO DNase (Thermo Fisher Scientific), and purified using the RNeasy MinElute Cleanup Kit (Qiagen). First-strand cDNA was synthesized using the Verso cDNA Synthesis Kit (Thermo Fisher Scientific) using a 3:1 mixture (v/v) of random hexamer and oligo dT primers according to the manufacturer’s protocol. Quantitative PCR was performed in optical 96-well plates (Bio-Rad) on a Bio-Rad CFX Connect Real-Time System using the Absolute Blue qPCR SYBR Green Kit (Thermo Fisher Scientific). The PCR program was as follows: 95 °C for 15 min, then 40 cycles at 95 °C for 15 s, 57 °C for 30 s, and 72 °C for 30 s, followed by a melt cycle from 55 to 95 °C in 0.5 s increments. Primers (Supplementary Data 4) were designed using Primer3web version 4.1.0. We verified the amplification specificity by sequencing each PCR product after cloning into the pCR™4-TOPO® TA vector (Thermo Fisher Scientific) and by melting-curve analyses. Primer efficiencies were calculated using a cDNA template dilution series. Gene expression was normalized to the expression level of eukaryotic initiation factor 4A (eIF4A), which showed the lowest variability across tissues among four tested reference genes (Supplementary Table 4).Expression of candidate transporters in insect cellsFor protein expression, we cloned each ORF without stop codon into the pIEx-4 expression vector (Novagen) in frame with the vector-encoded carboxy-terminal 6× His-tag and sequenced the resulting constructs. Primer sequences are listed in Supplementary Data 4. One construct of each candidate gene was used for transfection of High Five™ insect cells (Gibco) cultured in Express Five® SFM medium (Gibco) supplemented with 20 mM glutamine (Gibco) and 50 μg/mL gentamicin (Gibco). Confluent insect cells were diluted 1:5, dispensed in 500 μL aliquots into 24-well plates, and incubated at 27 °C. The next day, we transfected the cells using FuGENE HD Transfection Reagent (Promega) according to the manufacturer’s protocol. Cells treated with transfection reagent only were used as a negative control. After 48 h, we harvested the cells for Western blotting and uptake assays.Western blottingTo confirm protein expression, transfected insect cells were washed twice with phosphate-buffered saline (PBS; pH 7.4), collected by centrifugation, and resuspended in hypotonic buffer (20 mM Tris-HCl (pH 7.5), 5 mM EDTA, 1 mM dithiothreitol, 0.1% benzonase nuclease (Merck Millipore) (v/v), and protease inhibitors (cOmplete Mini, ETDA-free, Roche Diagnostics GmbH)). After incubation on ice for 10 min, the samples were frozen in liquid nitrogen, thawed, and centrifuged (16,000 × g for 15 min at 4 °C). The resulting cell pellet was resuspended in hypotonic buffer and used for Western blotting using horseradish peroxidase-conjugated anti-His antibody (1:10,000; Novex, Life technologies).Glucoside uptake assays with transfected insect cellsCells were washed with PBS (pH 5.5) by pipetting and incubated with different glucoside substrates at 200 µM in PBS (pH 5.5) for 1 h at 27 °C. Assays were performed with substrate mixtures containing seven different glucosinolates (2-propenyl glucosinolate (Roth), 4-methylsulfinylbutyl glucosinolate (Phytoplan), 4-methylthiobutyl glucosinolate (Phytoplan), 2-phenylethyl glucosinolate (Phytoplan), benzyl glucosinolate (Phytoplan), 4-hydroxybenzyl glucosinolate (isolated from Sinapis alba seeds61), and I3M glucosinolate (Phytoplan)), or five other plant glucosides (salicin (Sigma-Aldrich), linamarin (BIOZOL), dhurrin (Roth), catalpol (Sigma-Aldrich), and aucubin (Sigma-Aldrich)). After incubation, cells were washed three times with ice-cold PBS (pH 5.5) by pipetting, collected in 300 μL 80% (v/v) methanol, frozen in liquid nitrogen, thawed, and centrifuged at 3220 × g for 10 min at 4 °C. The supernatant was dried by vacuum centrifugation, dissolved in ultrapure water, and analyzed by liquid chromatography coupled with tandem mass spectrometry (LC-MS/MS). All glucoside uptake assays were performed in triplicates. The amount of each substrate in transporter-expressing cells was compared with that detected in control cells. Transporters were considered active when the average glucoside amounts detected in transporter-expressing cells were at least twofold higher than those detected in control cells.Cloning of PaGTR1 into the pNB1u vector and cRNA synthesisWe amplified the ORF of PaGTR1 without stop codon by PCR using uracil-containing primers (Supplementary Data 4). The 3′ primer was designed to encode a human influenza hemagglutinin tag to enable the detection of recombinant protein by Western blotting if necessary. The pNB1u vector was digested overnight at 37 °C with PacI and Nt.BbvCI (New England Biolabs) to generate 8-nt overhangs. One microliter gel-purified PCR product (100 ng/µL) was combined with 1 µL gel-purified vector (50 ng/µL), 1 U USER enzyme, 2 μL 5× PCR reaction buffer and 5 μL H2O, incubated at 37 °C for 25 min, followed by 25 min at room temperature. After the transformation of chemically competent E. coli cells, colonies containing the appropriate insert were identified by Sanger sequencing. The DNA template for complementary RNA (cRNA) synthesis was amplified by PCR from the X. laevis expression construct using pNB1uf/r primers (Supplementary Data 4) and cRNA was synthesized using the mMESSAGE mMACHINE™ T7 Transcription Kit (Invitrogen) according to the manufacturer’s manual.Biochemical characterization of PaGTR1 in Xenopus oocytesThe cRNA concentration was adjusted to 800 ng/μL with RNase-free water for oocyte injection. Xenopus laevis oocytes (Ecocyte Bioscience) were injected with 50 nL containing 40 ng cRNA or with 50 nL pure water as a control using a Drummond NANOJECT II (Drummond Scientific Company) or a Nanoliter 2010 Injector (World Precision Instruments).Injected oocytes were incubated in Kulori buffer (90 mM NaCl, 1 mM KCl, 1 mM MgCl2, 1 mM CaCl2, and 5 mM 4-(2-hydroxyethyl)-1-piperazineethanesulfonic acid, pH 7.4) supplemented with 50 μg/mL gentamicin at 16 °C for 3 days until assaying. All transport assays were performed at room temperature. Injected oocytes were pre-incubated for 5 min in Kulori buffer (90 mM NaCl, 1 mM KCl, 1 mM MgCl2, 1 mM CaCl2, and 5 mM 2-(N-morpholino)ethanesulfonic acid (MES), pH 6.0) before they were transferred into the same buffer containing the substrate(s). To determine the substrate preference of PaGTR1, we incubated oocytes in Kulori buffer (pH 6.0) containing an equimolar mixture of 2-propenyl glucosinolate, 4-methylsulfinylbutyl glucosinolate, 4-methylthiobutyl glucosinolate, 2-phenylethyl glucosinolate, benzyl glucosinolate, 4-hydroxybenzyl glucosinolate, and I3M glucosinolate, each at 200 µM, for 1 h.The pH dependency of glucosinolate transport was determined by incubating the injected oocytes with 100 μM I3M glucosinolate for 10 min in Kulori buffer adjusted to different pH values. The effects of different ionophores on glucosinolate transport was studied by incubating the PaGTR1-expressing oocytes in Kulori buffer at pH 6.0 containing either 20 μM carbonyl cyanide m-chlorophenyl hydrazone (H+ ionophore, Sigma-Aldrich), 20 μM nigericin (K+/H+ exchanger, Abcam), or 20 μM valinomycin (K+ ionophore, Abcam) for 15 min. Afterwards, we incubated the oocytes in Kulori buffer containing 100 μM I3M glucosinolate and the corresponding ionophore for 10 min. Assays performed with oocytes incubated in Kulori buffer without any ionophore served as a control.The time course of I3M glucosinolate uptake was analyzed by incubating oocytes with 100 μM I3M glucosinolate in Kulori buffer (pH 6.0) for 5, 10, 20, 30, 45, 60, 120, 180, and 240 min, respectively. The apparent Km value of PaGTR1 for I3M glucosinolate was determined by incubating injected oocytes in Kulori buffer (pH 6.0) for 10 min with different substrate concentrations. The Km value was calculated by nonlinear regression analysis in SigmaPlot 14.0 (Systat Software Inc.).Each assay consisted of 14–15 oocytes and was stopped by washing oocytes four times with Kulori buffer. Afterwards, 12 of the washed oocytes were distributed into four Eppendorf tubes, with three oocytes per tube and immediately homogenized in 100 µL of 50% (v/v) methanol. After centrifugation (21,380 × g or 16,000 × g for 15 min), the supernatant was incubated at −20 °C for at least 1 h to precipitate proteins, which were pelleted by centrifugation (21,380 × g or 16,000 × g for 15 min). Finally, 60 μL sample was diluted with 120 μL ultrapure water, filtered through a 0.22 μm PVDF-based filter plate (Merck Millipore), and analyzed by LC-MS/MS. The glucosinolate concentration in oocytes was calculated by assuming an oocyte volume of 1 μL62.LC-MS/MSGlucosinolates were quantified by LC-MS/MS using an Agilent 1200 HPLC system connected to an API3200 tandem mass spectrometer (AB SCIEX). Separation was achieved on an EC 250/4.6 NUCLEODUR Sphinx RP column (250 mm × 4.6 mm, 5 μm; Macherey-Nagel) using a binary solvent system consisting of 0.2% (v/v) formic acid in water (A) and acetonitrile (B), with a flow rate of 1 mL/min at 25 °C. The elution gradient was as follows: 0–1 min, 1.5% B; 1–6 min, 1.5–5% B; 6–8 min, 5–7% B; 8–18 min, 7–21% B; 18–23 min, 21–29% B; 23–23.1 min, 29–100% B; 23.1–24 min, 100% B; 24–24.1 min, 100 to 1.5% B; 24.1–28 min, 1.5% B. Glucosinolates were detected in negative ionization mode. The ion spray voltage was set to −4500 V. Gas temperature was set to 700 °C, curtain gas to 20 psi, collision gas to 10, nebulizing gas to 70 psi, and drying gas to 60 psi. Nonhost glucosides were quantified using an Agilent 1200 HPLC system connected to an API5000 tandem mass spectrometer (AB SCIEX). Separation was achieved on an Agilent XDB-C18 column (5 cm × 4.6 mm, 1.8 μm) using a binary solvent system consisting of 0.05% (v/v) formic acid in water (A) and acetonitrile (B) with a flow rate of 1.1 mL/min at 25 °C. The elution gradient was as follows: 0–0.5 min, 5% B; 0.5–2.5 min, 5–31% B; 2.5–2.52 min, 31–100% B; 2.52–3.5 min, 100% B; 3.5–3.51 min, 100 to 5% B; 3.51–6 min, 5% B. Compounds were detected in negative ionization mode with ion spray voltage set to −4500 V. The gas temperature was set to 700 °C, curtain gas to 30 psi, collision gas to 6, and both nebulizing gas and drying gas to 60 psi. Multiple reaction monitoring (MRM) was used to monitor the transitions from precursor ion to product ion for each compound (Supplementary Table 5). Compounds were quantified using external standard curves. Analyst Software 1.6 Build 3773 (AB Sciex) was used for data acquisition and processing.Samples from the pH dependency experiment were analyzed by LC-MS/MS using an Advance UHPLC system (Bruker) connected to an EVOQ Elite TripleQuad mass spectrometer (Bruker) equipped with an electrospray ion source. Separation was achieved on a Kinetex 1.7u XB-C18 column (100 mm × 2.1 mm, 1.7 µm, 100 Å, Phenomenex) using a binary solvent system consisting of 0.05 % (v/v) formic acid in water (A) and acetonitrile with 0.05% (v/v) formic acid (B), with a flow rate of 0.4 mL/min at 40 °C. The elution gradient was as follows: 0–0.2 min, 2% B; 0.2–1.8 min, 2–30% B; 1.8–2.5 min, 30–100% B; 2.5–2.8 min, 100% B; 2.8–2.9 min, 100 to 2% B; and 2.9–4.0 min, 2% B. Glucosinolates were detected in negative ionization mode. The instrument parameters were optimized by infusion experiments with pure standards. The ion spray voltage was set to −4000 V. Cone temperature was set to 350 °C, cone gas to 20 psi, heated probe temperature to 200 °C, and probe gas flow to 50 psi. Nebulizing gas was set to 60 psi, and collision gas to 1.6 mTorr. MRM parameters are provided in Supplementary Table 5. Bruker MS Workstation software (Version 8.2.1, Bruker) was used for data acquisition and processing of glucosinolates. All other samples from experiments using the oocyte expression system were analyzed using the LC-MS/MS method described above for insect cell-based assays. The concentrations of all glucosinolates in the substrate preference assay were determined using external standard curves. Assays performed to characterize I3M glucosinolate transport were quantified using 2-propenyl glucosinolate as an internal standard.Double-stranded RNA synthesisWe synthesized seven different double-stranded RNAs (dsRNAs) between 120 and 298 bp in length, one specific for each PaGTR1/2/3/5/6/7/8, respectively, and a 223-bp fragment of the inducible metalloproteinase inhibitor (IMPI) from the greater wax moth Galleria mellonella (AY330624.1) (dsIMPI) using the T7 RiboMAX™ Express RNAi System (Promega). In silico off-target prediction was done by searches of all possible 21 mers of both RNA strands against the local P. armoraciae transcriptome database allowing for two mismatches. Except for PaGTR5, the prediction did not find any off-target towards putative transporter genes in the transcriptome (Supplementary Data 5).Function of PaGTR1 in vivoTo analyze the function of PaGTR1, we injected newly emerged adult P. armoraciae beetles (reared on B. juncea) with 100 nL ultrapure water containing 80 ng of dsPaGTR1 or 80 ng of dsIMPI, respectively, using a Nanoliter 2010 Injector (World Precision Instruments). Injected beetles were provided with detached leaves of 3–4-week-old B. juncea plants and kept in plastic containers with moistened tissue in the laboratory under ambient conditions. Four days after dsRNA injection, we collected dsIMPI- and dsPaGTR1-injected beetles for gene expression analysis (n = 5 replicates, three beetles per replicate) and glucosinolate analysis (n = 10 replicates, three beetles per replicate), respectively. The remaining beetles were used for a sequestration experiment with Arabidopsis thaliana Col-0 (Arabidopsis) plants that had been cultivated in a growth chamber at 21 °C, 55% relative humidity, and a 10-h photoperiod. To compare the accumulation and excretion of ingested glucosinolates in dsIMPI- and dsPaGTR1-injected beetles, we fed beetles with detached Arabidopsis leaves in Petri dishes (60 mm diameter) that contained 50 μL of ultrapure water and were sealed with parafilm to prevent leaf wilting (n = 10 replicates, five beetles per replicate). Feeding assays were performed in the laboratory under ambient conditions, and leaves were exchanged every day for five consecutive days. To estimate how much the beetles fed, we determined the weight of each leaf before and after feeding. Feces were collected every day using 100 µL of ultrapure water per replicate, combined with 300 μL of pure methanol in a 1.5 mL Eppendorf tube and dried by vacuum centrifugation. Feces samples were then homogenized in 200 μL of 80% (v/v) methanol using metal beads (2.4 mm diameter, Askubal) in a tissue lyzer (Qiagen) for 1 min at 30 Hz. After feeding, adults were starved for one day, weighed, frozen in liquid nitrogen, and stored at −20 °C until glucosinolate extraction. Beetle samples were homogenized using a plastic pestle in 200 μL 80% (v/v) methanol. All samples were then extracted with 1 mL 80% (v/v) methanol containing 25 μM 4-hydroxybenzyl glucosinolate as an internal standard. After centrifugation (16,000 × g for 10 min), glucosinolates were extracted from the supernatant, converted to desulfo-glucosinolates, and analyzed by high-performance liquid chromatography coupled with diode array detection (HPLC-DAD) as described below. The glucosinolate content in adults or feces was calculated in nanomole per adult, respectively.To confirm the specificity of PaGTR1 knockdown, we analyzed the effect of dsPaGTR1 injection on the expression of PaGTR2, PaGTR3, PaGTR9, and PaGTR10. PaGTR9 and PaGTR10 share the highest nucleotide sequence similarity (69% sequence identity) with PaGTR1. PaGTR2 and PaGTR3 expression was analyzed because the recombinant transporters also used I3M glucosinolate as a substrate. RNA extraction, purification, cDNA synthesis, and qPCR were performed as described above.Function of PaGTR5/6/7/8 in vivoTo analyze the function of PaGTR5/6/7/8, we injected newly emerged adult P. armoraciae beetles that had fed for 2 days on B. juncea leaves with 100 nL ultrapure water containing 100 ng of dsIMPI or each 100 ng of dsPaGTR5/6/7/8 using a Nanoliter 2010 Injector (World Precision Instruments). A subset of the dsRNA-injected beetles was fed with detached leaves of Arabidopsis plants in plastic containers with moistened tissue for gene expression analysis (n = 6 replicates, two beetles per replicate). The remaining dsRNA-injected beetles were used for a sequestration experiment with Arabidopsis plants to compare the accumulation and excretion of previously stored and ingested glucosinolates in dsIMPI- and dsPaGTR5/6/7/8-injected beetles. We fed the injected beetles with detached Arabidopsis leaves in Petri dishes (60 mm diameter) that contained 30 μL of ultrapure water and were sealed with parafilm to prevent leaf wilting (n = 10 replicates, six beetles per replicate). Feeding assays were performed in the laboratory under ambient conditions, and leaves were exchanged every day for six consecutive days. To estimate how much the beetles fed, we determined the weight of each leaf before and after feeding. Starting from the second day, feces were collected as above for 5 days. After drying by vacuum centrifugation, feces were homogenized in 1 mL of 80% (v/v) methanol containing 25 μM 4-hydroxybenzyl glucosinolate as an internal standard using metal beads (2.4 mm diameter, Askubal) in a tissue lyzer (Qiagen) for 1 min at 30 Hz. Fed adults were starved for 1 day, weighed, and frozen in liquid nitrogen until glucosinolate extraction. Beetle samples were homogenized using a plastic pestle in 200 μL of 80% (v/v) methanol and then extracted with 1 mL of 80% (v/v) methanol containing 25 μM 4-hydroxybenzyl glucosinolate as an internal standard. After centrifugation (16,000 × g for 10 min), glucosinolates were extracted from the supernatant, converted to desulfo-glucosinolates, and analyzed by HPLC-DAD as described below.To confirm the specificity of PaGTR5/6/7/8 knockdown, we analyzed the effect of dsPaGTR5/6/7/8 injection on the expression of PaGTR9 and PaGTR10, which share the highest nucleotide sequence similarity (67–69% sequence identity) with PaGTR5/6/7/8. RNA extraction, purification, cDNA synthesis, and qPCR were performed as described above.Function of PaGTR2/3 in vivoTo analyze the functions of PaGTR2 and PaGTR3, we injected third instar larvae of P. armoraciae (reared on B. rapa) with 100 nL ultrapure water containing 80 ng of dsPaGTR2 and 80 ng of dsPaGTR3 (dsPaGTR2/3) or 80 ng of dsIMPI, respectively, using a Nanoliter 2010 Injector (World Precision Instruments). Injected larvae were provided with detached B. rapa petioles and kept in plastic tubes with moistened tissue in the laboratory under ambient conditions until pupation. Newly emerged adults were again injected with dsRNAs and provided with detached B. rapa leaves. Three days after the second dsRNA injection, we collected dsIMPI- and dsPaGTR2/3-injected beetles for gene expression analysis (n = 6 replicates, two beetles per replicate) and glucosinolate analysis (n = 10 replicates, three beetles per replicate), respectively. The remaining beetles were used for a feeding experiment with Arabidopsis. Each replicate consisted of one Arabidopsis leaf and three beetles that were placed in a Petri dish (60 mm diameter) with 50 μL of ultrapure water (n = 12–13 replicates). Each leaf was photographed before and after feeding to determine the consumed leaf area. Fed leaves were frozen in liquid nitrogen, freeze-dried, and homogenized using metal beads (2.4 mm diameter, Askubal) in a tissue lyzer (Qiagen) for 2 min at 30 Hz. Fed beetles were starved for 1 day, weighed, frozen in liquid nitrogen, and stored at −20 °C until glucosinolate extraction. Beetle samples were homogenized using a plastic pestle in 200 μL 80% (v/v) methanol. All samples were extracted with 1 mL 80% (v/v) methanol containing 25 μM 1-methylethyl glucosinolate (extracted from Sisymbrium officinale seeds) as an internal standard. Extracts were applied to DEAE-Sephadex A-25 (Sigma-Aldrich) columns in 96-well filter plates (Nunc) that were preconditioned with 1 mL of ultrapure water and 2 × 500 μL of 80% (v/v) methanol. After loading 900 μL of extract, the columns were washed with 500 μL of 80% (v/v) methanol, followed by 2 × 1 mL of ultrapure water. To adjust the pH condition, 500 μL of 0.02 M MES buffer (pH 5.2) was added to each column. After adding 30 μL of Helix pomatia sulfatase solution to each column and overnight incubation at room temperature, desulfo-glucosinolates were eluted using 300 or 500 μL of ultrapure water into 96-deep well plates (Nunc). Samples were analyzed by desulfo-HPLC-DAD as described below. The glucosinolate content in adults and fed leaves was calculated in nanomole per adult and nanomole per cm2 leaf, respectively. The ingested glucosinolate amount was calculated based on the ingested leaf area and the corresponding leaf glucosinolate content. To elucidate which proportion of the ingested glucosinolates was sequestered, we expressed the glucosinolate amount detected in the beetles relative to the ingested glucosinolate amount from the leaves, which was set to 100%.HPLC-DADSamples were analyzed by HPLC on an Agilent Technologies HP1100 Series instrument equipped with a photodiode array detector. After injection of 100 μL of each sample, separation was achieved on an EC 250/4.6 NUCLEODUR Sphinx RP column (250 mm × 4.6 mm, 5 μm; Macherey-Nagel; samples of the RNA interference (RNAi) experiments of PaGTR1, and PaGTR2/3) or an EC 250/4.6 NUCLEODUR 100-5 C18ec column (250 mm × 4.6 mm, 5 μm; Macherey-Nagel; samples of the RNAi experiments of PaGTR5/6/7/8, and PaGTR2/3) using a binary solvent system consisting of ultrapure water (A) and acetonitrile (B), with a flow rate of 1 mL/min. The elution gradient was as follows: 0–1 min, 1.5% B; 1–6 min, 1.5–5% B; 6–8 min, 5–7% B; 8–18 min, 7–21% B; 18–23 min, 21–29% B; 23–23.5 min, 29–100% B; 23.5–26 min, 100% B; 26–26.1 min, 100 to 1.5% B; and 26.1–31 min, 1.5% B. The eluent was monitored by diode array detection between 190 and 360 nm. Desulfo-glucosinolates were identified based on a comparison of retention times and absorption spectra with those of known standards63. The glucosinolate content in each sample was calculated from the peak areas at 229 nm relative to the peak area of the internal standard using relative response factors64.Glucosinolate concentration in the hemolymph of adult P. armoraciae
Hemolymph was collected from 7-day-old adult P. armoraciae reared on B. juncea by cutting off an abdominal leg and collecting the extruding droplet using glass capillaries (0.5 µL, Hirschmann® minicaps®) (n = 6 replicates, 50 beetles per replicate). The capillaries were marked with 1 mm intervals (corresponding to 15.6 nL) to estimate the volume of collected hemolymph. The hemolymph was diluted in 500 µL of 90% (v/v) methanol, homogenized using metal beads (2.4 mm diameter, Askubal) in a tissue lyzer (Qiagen) for 1 min at 30 Hz, and boiled for 5 min at 95 °C. After two centrifugation steps (13,000 × g for 10 min each), the supernatant was dried by vacuum centrifugation, dissolved in 50 µL 50% (v/v) methanol, diluted in ultrapure water, and analyzed by LC-MS/MS.Morphology of the Malpighian tubule system of P. armoraciae
To investigate the structure of the Malpighian tubule system, we dissected the gut and Malpighian tubules of 4-day-old P. armoraciae adults in PBS (pH 6.0) under a stereomicroscope. The tracheae that attach Malpighian tubules to the gut were removed using fine forceps to release the tubules. Pictures were taken with a Canon EOS 600D camera.pH of hemolymph and excretion fluid of isolated Malpighian tubules of P. armoraciae
The pH of the hemolymph and Malpighian tubule excretion fluid of 5-day-old P. armoraciae adults was assessed using the pH indicator bromothymol blue (Alfa Aesar). Hemolymph was collected by cutting off an abdominal leg and collecting the extruding droplet using a pipette (n = 3 replicates, six to ten beetles per replicate). Excretion fluid was collected from dissected Malpighian tubules that were incubated in saline A as described for the Ramsay assay (n = 4 replicates, one tubule per replicate). Hemolymph and excretion fluid were immediately mixed with the same volume of 0.16% (w/v) bromothymol blue dissolved in 10% (v/v) ethanol, respectively, under water-saturated paraffin oil in a Sylgard-coated petri dish. The resulting color of the droplet was compared with those of citric acid–Na2HPO4 buffer solutions ranging from pH 5.2 to 6.6 in 0.2 increments mixed with 0.16% (w/v) bromothymol blue.Fate of plant glucosides injected in P. armoraciae beetlesTo analyze the fate of plant glucosides in vivo, we injected 100 nL of an equimolar mixture of 2-propenyl glucosinolate, 4-hydroxybenzyl glucosinolate, linamarin, salicin, and catalpol, each at 10 mM, and 0.15% (w/v) amaranth into the hemolymph of 2-day-old adult P. armoraciae (reared on B. rapa). One group of beetles was sampled 30 min after injection by freezing beetles in liquid nitrogen (n = 10 replicates, five beetles per replicate). The remaining beetles were fed with detached leaves of Arabidopsis in Petri dishes (60 mm diameter) in the laboratory under ambient conditions (n = 10 replicates, five beetles per replicate). We added 30 μL of ultrapure water to each Petri dish and sealed them with parafilm to prevent leaf wilting. After 1 day, we sampled the beetles as described above (n = 10 replicates, five beetles per replicate). Feces were collected using 100 µL of ultrapure water per replicate and combined with 300 μL of pure methanol in a 1.5 mL Eppendorf tube. All samples were stored at −20 °C until extraction. Beetle and feces samples were homogenized as described in the RNAi experiment. After centrifugation (16,000 × g for 10 min), the supernatant was dried by vacuum centrifugation, dissolved in 100 µL of ultrapure water, and analyzed by LC-MS/MS. The glucoside content in adults or feces was calculated as nanomole per adult, respectively.Ramsay assayTo analyze the excretion of plant glucosides in situ, we performed Ramsay assays65 with dissected Malpighian tubules from 4- to 5-day-old P. armoraciae adults reared on B. rapa (Supplementary Fig. 6). Malpighian tubules were dissected in saline A (100 mM NaCl, 8.6 mM KCl, 2 mM CaCl2, 8.5 mM MgCl2, 4 mM NaH2PO4, 4 mM NaHCO3, 24 mM glucose, 10 mM proline, 25 mM 3-(N-morpholino)propanesulfonic acid (MOPS), pH 6.0)66. Single tubules were transferred into a 10 μL droplet of saline B (60 mM NaCl, 10.3 mM KCl, 2.4 mM CaCl2, 10.2 mM MgCl2, 4.8 mM NaH2PO4, 4.8 mM NaHCO3, 28.8 mM glucose, 12 mM proline, 30 mM MOPS, 1 mM cyclic AMP, pH 6.0) under water-saturated paraffin oil in a Sylgard-coated petri dish. The proximal end of the tubule was pulled out of the droplet, attached to a metal pin, and cut using a glass capillary to allow the collection of excretion fluid. To start the assay, we added 2 μL of an equimolar glucoside mixture consisting of 2-propenyl glucosinolate, 4-methylsulfinylbutyl glucosinolate, 4-hydroxybenzyl glucosinolate, 2-phenylethyl glucosinolate, I3M glucosinolate, linamarin, salicin, and catalpol, each at 40 mM, and 0.6% (w/v) amaranth) to saline B. After 2–3 h, we collected the excretion fluid and 2 μL of the bathing droplet in 300 μL of 80% (v/v) methanol, respectively. The Malpighian tubule was washed three times in ~15 mL of saline A and afterwards transferred into 300 μL of 80% (v/v) methanol. All samples were stored at −20 °C until extraction. Malpighian tubule samples were homogenized using a plastic pestle. After centrifugation (16,000 × g for 10 min), the supernatant was dried by vacuum centrifugation, dissolved in 70 µL ultrapure water, and analyzed by LC-MS/MS. Out of 36 assays, 25 assays were excluded because no excretion fluid was visible within 2–3 h.Statistical analysisNo statistical methods were used to predetermine sample size. Statistical analyses were conducted in R 3.5.167 or in SigmaPlot 14.0 or 14.5 (Systat Software Inc.). Two groups were compared by two-tailed Student’s t test, Mann–Whitney U test, or the method of generalized least squares68, depending on the variance homogeneity and the normality of residuals. Data of the Ramsay assays were compared by paired two-tailed Student’s t test using FDR (false discovery rate)-corrected P values69. Three or more groups were compared by one-way analysis of variance, followed by post hoc multiple comparisons test, or the method of generalized least squares68. If necessary, data were transformed prior to analysis. For comparisons using the method of generalized least squares, we applied the varIdent variance structure, which allows each group to have a different variance. The P value was obtained by comparing models with and without explanatory variables using a likelihood ratio test70. Significant differences between groups were revealed with factor level reductions71. Information about data transformation, statistical methods, and results of the statistical analyses are summarized in Supplementary Tables 2 and 3 and Supplementary Data 6.Reporting summaryFurther information on research design is available in the Nature Research Reporting Summary linked to this article. More

1.McGill, B. J., Enquist, B. J., Weiher, E. & Westoby, M. Rebuilding community ecology from functional traits. Trends Ecol. Evol. 21, 178–185 (2006).PubMed 
Article 
PubMed Central 

Google Scholar 
2.Violle, C. et al. Let the concept of trait be functional! Oikos 116, 882–892 (2007).Article 

Google Scholar 
3.Adler, P. B. et al. Functional traits explain variation in plant life history strategies. Proc. Natl. Acad. Sci. 111, 740–745 (2014).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
4.Pierce, S., Brusa, G., Vagge, I. & Cerabolini, B. E. Allocating C. S. R. plant functional types: the use of leaf economics and size traits to classify woody and herbaceous vascular plants. Funct. Ecol. 27, 1002–1010 (2013).5.Kattge, J. et al. TRY – a global database of plant traits. Glob. Change Biol. 17, 2905–2935 (2011).ADS 
Article 

Google Scholar 
6.Kunstler, G. et al. Plant functional traits have globally consistent effects on competition. Nature 529, 204–207 (2016).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
7.Yang, J., Cao, M. & Swenson, N. G. Why functional traits do not predict tree demographic rates. Trends Ecol. Evol. 33, 326–336 (2018).PubMed 
Article 
PubMed Central 

Google Scholar 
8.Connolly, S. R., Keith, S. A., Colwell, R. K. & Rahbek, C. Process, mechanism, and modeling in macroecology. Trends Ecol. Evol. 32, 835–844 (2017).PubMed 
Article 
PubMed Central 

Google Scholar 
9.Dormann, C. F. et al. Correlation and process in species distribution models: bridging a dichotomy. J. Biogeogr. 39, 2119–2131 (2012).Article 

Google Scholar 
10.Thuiller, W., Pollock, L. J., Gueguen, M. & Münkemüller, T. From species distributions to meta-communities. Ecol. Lett. 18, 1321–1328 (2015).PubMed 
PubMed Central 
Article 

Google Scholar 
11.Zurell, D. et al. Benchmarking novel approaches for modelling species range dynamics. Glob. Change Biol. 22, 2651–2664 (2016).ADS 
Article 

Google Scholar 
12.Alexander, J. et al. Lags in the response of alpine plant communities to climate change. Glob. Change Biol. 24, 563–579 (2018).ADS 
Article 

Google Scholar 
13.Evans, M. E., Merow, C., Record, S., McMahon, S. M. & Enquist, B. J. Towards process-based range modeling of many species. Trends Ecol. Evol. 31, 860–871 (2016).PubMed 
Article 
PubMed Central 

Google Scholar 
14.Urban, M. C. et al. Improving the forecast for biodiversity under climate change. Science 353, aad8466 (2016).PubMed 
Article 
CAS 
PubMed Central 

Google Scholar 
15.Paine, C. E. T., Deasey, Anna, Bradley, DuthieA. & Ken, Thompson Towards the general mechanistic prediction of community dynamics. Funct. Ecol. 32, 1681–1692 (2018).Article 

Google Scholar 
16.Hartig, F. et al. Connecting dynamic vegetation models to data–an inverse perspective. J. Biogeogr. 39, 2240–2252 (2012).Article 

Google Scholar 
17.Levins, R. The strategy of model building in population biology. Am. Sci. 54, 421–431 (1966).
Google Scholar 
18.Kraft, N. J., Godoy, O. & Levine, J. M. Plant functional traits and the multidimensional nature of species coexistence. Proc. Natl. Acad. Sci. 112, 797–802 (2015).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
19.Mayfield, M. M. & Stouffer, D. B. Higher-order interactions capture unexplained complexity in diverse communities. Nat. Ecol. Evol. 1, 0062 (2017).Article 

Google Scholar 
20.Carrara, F., Giometto, A., Seymour, M., Rinaldo, A. & Altermatt, F. Experimental evidence for strong stabilizing forces at high functional diversity of aquatic microbial communities. Ecology 96, 1340–1350 (2015).PubMed 
Article 
PubMed Central 

Google Scholar 
21.Laughlin, D. C. & Messier, J. Fitness of multidimensional phenotypes in dynamic adaptive landscapes. Trends Ecol. Evol. 30, 487–496 (2015).PubMed 
Article 
PubMed Central 

Google Scholar 
22.Curtsdotter, A. et al. Ecosystem function in predator–prey food webs—confronting dynamic models with empirical data. J. Anim. Ecol. 88, 196–210 (2019).PubMed 
Article 
PubMed Central 

Google Scholar 
23.Chesson, P. Mechanisms of maintenance of species diversity. Annu. Rev. Ecol. Syst. 31, 343–366 (2000).24.Terborgh, J. W. Toward a trophic theory of species diversity. Proc. Natl. Acad. Sci. 112, 11415–11422 (2015).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
25.Swenson, N. G. & Enquist, B. J. Ecological and evolutionary determinants of a key plant functional trait: wood density and its community-wide variation across latitude and elevation. Am. J. Bot. 94, 451–459 (2007).PubMed 
Article 
PubMed Central 

Google Scholar 
26.Körner, C. Alpine Plant Life. (Springer, 2003).27.Westoby, M. A leaf-height-seed (LHS) plant ecology strategy scheme. Plant Soil 199, 213–227 (1998).CAS 
Article 

Google Scholar 
28.Adler, P. B. et al. Competition and coexistence in plant communities: intraspecific competition is stronger than interspecific competition. Ecol. Lett. 21, 1319–1329 (2018).PubMed 
Article 
PubMed Central 

Google Scholar 
29.Chalmandrier, L., Albouy, C. & Pellissier, L. Species pool distributions along functional trade-offs shape plant productivity–diversity relationships. Sci. Rep. 7, 15405 (2017).ADS 
PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
30.Nagelkerke, N. J. A note on a general definition of the coefficient of determination. Biometrika 78, 691–692 (1991).MathSciNet 
MATH 
Article 

Google Scholar 
31.Gelman, A., Hwang, J. & Vehtari, A. Understanding predictive information criteria for Bayesian models. Stat. Comput. 24, 997–1016 (2014).MathSciNet 
MATH 
Article 

Google Scholar 
32.Godínez-Alvarez, H., Herrick, J. E., Mattocks, M., Toledo, D. & Van Zee, J. Comparison of three vegetation monitoring methods: Their relative utility for ecological assessment and monitoring. Ecol. Indic. 9, 1001–1008 (2009).Article 

Google Scholar 
33.Grime, J. P. Evidence for the existence of three primary strategies in plants and its relevance to ecological and evolutionary theory. Am. Nat. 111, 1169–1194 (1977).Article 

Google Scholar 
34.Reich, P. B. The world-wide ‘fast–slow’ plant economics spectrum: a traits manifesto. J. Ecol. 102, 275–301 (2014).Article 

Google Scholar 
35.de Bello, Fde et al. Hierarchical effects of environmental filters on the functional structure of plant communities: a case study in the French Alps. Ecography 36, 393–402 (2013).Article 

Google Scholar 
36.Dray, S. et al. Combining the fourth-corner and the RLQ methods for assessing trait responses to environmental variation. Ecology 95, 14–21 (2014).PubMed 
Article 
PubMed Central 

Google Scholar 
37.Violle, C. et al. Competition, traits and resource depletion in plant communities. Oecologia 160, 747–755 (2009).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
38.Michelsen, A., Schmidt, I. K., Jonasson, S., Quarmby, C. & Sleep, D. Leaf 15N abundance of subarctic plants provides field evidence that ericoid, ectomycorrhizal and non-and arbuscular mycorrhizal species access different sources of soil nitrogen. Oecologia 105, 53–63 (1996).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
39.Laughlin, D. C., Joshi, C., van Bodegom, P. M., Bastow, Z. A. & Fulé, P. Z. A predictive model of community assembly that incorporates intraspecific trait variation. Ecol. Lett. 15, 1291–1299 (2012).PubMed 
Article 
PubMed Central 

Google Scholar 
40.Laughlin, D. C. The intrinsic dimensionality of plant traits and its relevance to community assembly. J. Ecol. 102, 186–193 (2014).Article 

Google Scholar 
41.Díaz, S. et al. The global spectrum of plant form and function. Nature 529, 167–171 (2016).ADS 
Article 
CAS 

Google Scholar 
42.Letten, A. D., Ke, P.-J. & Fukami, T. Linking modern coexistence theory and contemporary niche theory. Ecol. Monogr. 87, 161–177 (2017).Article 

Google Scholar 
43.O’Dwyer, J. P. Whence Lotka-Volterra? Theor. Ecol. 11, 441–452 (2019).Article 

Google Scholar 
44.Pichler, M., Boreux, V., Klein, A.-M., Schleuning, M. & Hartig, F. Machine learning algorithms to infer trait-matching and predict species interactions in ecological networks. Methods Ecol. Evol. 11, 281–293 (2020).Article 

Google Scholar 
45.Maestre, F. T., Callaway, R. M., Valladares, F. & Lortie, C. J. Refining the stress-gradient hypothesis for competition and facilitation in plant communities. J. Ecol. 97, 199–205 (2009).Article 

Google Scholar 
46.Damgaard, C. A critique of the space-for-time substitution practice in community ecology. Trends Ecol. Evol. 34, 416–421 (2019).PubMed 
Article 
PubMed Central 

Google Scholar 
47.Arnoldi, J.-F., Loreau, M. & Haegeman, B. The inherent multidimensionality of temporal variability: how common and rare species shape stability patterns. Ecol. Lett. 22, 1557–1567 (2019).PubMed 
PubMed Central 
Article 

Google Scholar 
48.May, R. M. Patterns of species abundance and diversity. in Ecology and Evolution of Communities (eds Cody, M. L. & Diamond, J. M.) 81–120 (Harvard University Press, 1975).49.Falster, D. S., Brännström, Å., Westoby, M. & Dieckmann, U. Multitrait successional forest dynamics enable diverse competitive coexistence. Proc. Natl. Acad. Sci. 114, E2719–E2728 (2017).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
50.Boulangeat, I., Georges, D. & Thuiller, W. FATE-HD: a spatially and temporally explicit integrated model for predicting vegetation structure and diversity at regional scale. Glob. Change Biol. 20, 2368–2378 (2014).Article 

Google Scholar 
51.Lischke, H. & Löffler, T. J. Finding all multiple stable fixpoints of n-species Lotka–Volterra competition models. Theor. Popul. Biol. 115, 24–34 (2017).PubMed 
MATH 
Article 
PubMed Central 

Google Scholar 
52.Jamil, T., Ozinga, W. A., Kleyer, M. & ter Braak, C. J. Selecting traits that explain species–environment relationships: a generalized linear mixed model approach. J. Veg. Sci. 24, 988–1000 (2013).Article 

Google Scholar 
53.ter Braak, C. J. New robust weighted averaging-and model-based methods for assessing trait–environment relationships. Methods Ecol. Evol. 10, 1962–1971 (2019).Article 

Google Scholar 
54.Pielou, E. C. Ecological Diversity. (Wiley, New York, 1975).55.Gibert, A., Gray, E. F., Westoby, M., Wright, I. J. & Falster, D. S. On the link between functional traits and growth rate: meta-analysis shows effects change with plant size, as predicted. J. Ecol. 104, 1488–1503 (2016).Article 

Google Scholar 
56.Durand, Y. et al. Reanalysis of 44 Yr of climate in the French Alps (1958–2002): methodology, model validation, climatology, and trends for air temperature and precipitation. J. Appl. Meteorol. Climatol. 48, 429–449 (2009).ADS 
Article 

Google Scholar 
57.Chalmandrier, L. et al. Spatial scale and intraspecific trait variability mediate assembly rules in alpine grasslands. J. Ecol. 105, 277–287 (2017).Article 

Google Scholar 
58.Cornelissen, J. H. C. et al. A handbook of protocols for standardised and easy measurement of plant functional traits worldwide. Aust. J. Bot. 51, 335–380 (2003).Article 

Google Scholar 
59.Reich, P. B. et al. Generality of leaf trait relationships: a test across six biomes. Ecology 80, 1955–1969 (1999).Article 

Google Scholar 
60.Poorter, H. & Bergkotte, M. Chemical composition of 24 wild species differing in relative growth rate. Plant Cell Environ. 15, 221–229 (1992).CAS 
Article 

Google Scholar 
61.Farquhar, G. D., O’leary, M. H. & Berry, J. A. On the relationship between carbon isotope discrimination and the intercellular carbon dioxide concentration in leaves. Funct. Plant Biol. 9, 121–137 (1982).CAS 
Article 

Google Scholar 
62.Thioulouse, J. et al. Multivariate Analysis of Ecological Data with ade4. (Springer, 2018).63.Rapisarda, F., Brigo, D. & Mercurio, F. Parameterizing correlations: a geometric interpretation. IMA J. Manag. Math. 18, 55–73 (2007).MathSciNet 
MATH 
Article 

Google Scholar 
64.Blumenson, L. E. A derivation of n-dimensional spherical coordinates. Am. Math. Mon. 67, 63–66 (1960).MathSciNet 

Google Scholar 
65.Banner, K. M., Irvine, K. M. & Rodhouse, T. The use of bayesian priors in ecology: the good, the bad, and the not great. Methods Ecol. Evol. 00, 1–8 (2020).
Google Scholar 
66.Hartig, F., Minunno, F. & Paul, S. BayesianTools: General-Purpose MCMC and SMC samplers and tools for bayesian statistics. R package (2017).67.Warton, D. I. et al. So Many Variables: Joint Modeling in Community Ecology. Trends Ecol. Evol. 30, 766–779 (2015). 68.Pichler, M. & Hartig, F. A new method for faster and more accurate inference of species associations from novel community data. Preprint at https://arxiv.org/abs/2003.05331 (2020).69.Advanced Research Computing Center (ARCC). Teton Computing Environment. https://doi.org/10.15786/m2fy47 (2018). More

Study areaChina is one of the countries with the poorest per capita water resources in the world while also having the largest water consumption in the world. In 2018, China’s total water consumption was 601.55 billion m3, with 369.31 billion m3 of water used in agriculture, accounting for 61.4% of the total water2. Agriculture is the most important industrial sector in water resource consumption. However, due to regional and climate differences, the distribution of agricultural water resources is uneven, and the shortage of water resources seriously affects agricultural development in water-deficient areas.Figure 1 shows the agricultural water consumption in China by province for 2007 and 2018. The agricultural water consumption includes farmland irrigation water consumption (classified as paddy field, irrigated land, vegetable field, groundwater exploitation), forest, animal husbandry, fishery, and livestock (classified as forest and fruit, grassland, fish pond, animal husbandry, groundwater exploitation), domestic water consumption of rural residents and rural ecological environment water consumption. Previous studies have mainly considered the irrigation water consumption of the planting industry as the research object at the provincial or regional levels (e.g., eastern, central, and western regions). Few were able to consider all 31 provinces in China and have comprehensively assessed water consumption and water use efficiency in the various types of agricultural production3,4,5,6,10,16,17,22,23,24,25,30. In this study, the agricultural water use efficiency and its influencing factors are assessed based on the agricultural water consumption of agriculture, forestry, animal husbandry, and fishery in China.Figure 1Agricultural water consumption in China by province for (a) 2007 and (b) 2018. Note: Map created using ArcGIS [10.2], (http://www.esri.com/software/arcgis).Full size imageResearch methodIn this study, the agricultural water use efficiency (under the common frontier and the group frontier) is calculated using the super-efficiency slacks-based measure (Super-SBM) model. The significant factors affecting water-use efficiency are then analyzed through the threshold regression model.Super-efficiency SBM modelData envelopment analysis (DEA) is an efficiency evaluation method proposed by Charnes31, a famous American operational research scientist. While traditional radial and angular DEA models do not require the specific form of the estimation function, they ignore the relaxation of variables and result in efficiency values in the range of 0 to 1. If there are multiple efficiency value of decision making units(DMUs) with an efficiency value of 1, these values cannot be compared. The efficiency of the super efficiency DEA model can be greater than 1, which means that the efficiency level of all decision-making units can be compared.To avoid the problem of slack variables, Tone (2001) proposed the SBM model, which is a non-radial and non-angular DEA analysis method based on the relaxation variable measure16,17,18,19,20,32. The SBM model of unexpected output solves the slack problem of input and output variables, minimizing deviations in the efficiency measurement. The super-efficiency SBM model combines the super-efficiency DEA model and the SBM model. It is also one of the methods based on data envelopment analysis, which can measure the efficiency of all decision-making units and the slack of input and output variables.Assume n to be the decision-making units, each of which has m inputs, expected output r1, and unexpected output r2. Let X (X ∈ Rm), Yd (Yd ∈ Rs1), and Yu (Yu ∈ Rs2) be matrices, such that (X=[{x}_{1},dots ,{x}_{n}]in {R}^{m*n}) and (Y=[{y}_{1}^{d}, dots ,{ y}_{n}^{d}in {R}^{{r}_{1}*n}). The form of the super-efficiency SBM model is as follows1,17,19,54:$$min=frac{frac{1}{m}sum_{i=1}^{m}(overline{x}/{x}_{ik})}{1/left({r}_{1}+{r}_{2}right)*(sum_{r=1}^{{r}_{1}}overline{{y}^{d}}/{y}_{rk}^{d}+sum_{q=1}^{{r}_{2}}overline{{y}^{u}}/{y}_{qk}^{u}}.$$
(1)
Among them,$$overline{x}ge sum_{j=1ne k}^{n}{x}_{ij}{lambda }_{j}, i=1,dots ,m;$$
(2)
$$overline{{y}^{d}}le sum_{j=1,ne k}^{n}{y}_{rj}^{d}{lambda }_{j}, r=1,dots ,{s}_{1};$$
(3)
$$overline{{y}^{d}}ge sum_{j=1,ne k}^{n}{y}_{qj}^{u}{lambda }_{j}, q=1,dots ,{s}_{2};$$
(4)
$${lambda }_{y}ge 0,j=1,dots ,n;jne 0;$$
(5)
$$overline{x}ge {x}_{k},k=1,dots ,m;$$
(6)
$$overline{{y}^{d}}le {y}_{k}^{d},d=1,dots ,{r}_{1};$$
(7)
$$overline{{y}^{u}}ge {y}_{k}^{u},b=1,dots ,{r}_{2}.$$
(8)
Based on the Super-SBM model (Eq. 1) and its constraint formula, the agricultural water use efficiency for the different provinces was calculated for the period 2007–2018 using Maxdea 8 ultra software.Threshold effectConsidering the differences in economic development and technical levels, the agricultural water use in different regions of China shows characteristics of time-series evolution, spatial heterogeneity, and unbalanced spatial distribution. There is a non-linear relationship between the influencing factors of agricultural water use efficiency, which suggests the existence of certain threshold characteristics33,34. This means that for a particular determinant, agricultural water use efficiency would be affected differently depending on whether the parameter has crossed the threshold. In this study, the threshold panel model proposed by Hansen is used. The threshold value of the threshold variable is taken as the critical point, and the regression equation is divided into different stage intervals to analyze the influence of threshold variables on the explained variables at different stages . Therefore, according to the relationship between agricultural water use efficiency and its influencing factors in different regions, the following single threshold regression model is set:$${Y}_{it}=alpha {X}_{it}+{beta }_{1}{T}_{it}Ileft({T}_{it}le {gamma }_{1}right)+{beta }_{2}{T}_{it}Ileft({T}_{it} >{gamma }_{1}right)+C+{varepsilon }_{it},$$
(9)
such that i is the province; t is the year; Yit and Tit are the explanatory variables and explained variables, respectively; Xit is the control variable that has a significant impact on the explained variables; Tit is threshold variable, which changes with the different explanatory variables; γ is a specific threshold value; α is the corresponding coefficient vector; β1 and β2 represent the influence coefficients of the threshold variable Tit on the explained variable Yit in the case of ({T}_{it}le {gamma }_{1}) and ({T}_{it} >{gamma }_{1}) , respectively; C is a constant; ε is random disturbance term, ({varepsilon }_{it}sim i.i.d.N(0,{sigma }^{2})); and, I (·) is an indicative function. After obtaining the estimated value of each parameter, two tests need to be carried out: (1) establish whether the threshold effect is significant; and (2) determine whether the estimated threshold value is equal to the true value. In addition, the above equation assumes that only one threshold exists. For two or more thresholds, the model would have to be adjusted according to the data.Based on the panel data of 31 provinces in China from 2007 to 201844,45,46, Stata15.0 software was used to perform threshold regression on seven variables: per capita water resources, rural labor force, disposable income, government’s attention, foreign trade dependence, industrial structure, and gross domestic product (GDP). The threshold effect of each factor can be analyzed, and the impact on agricultural water consumption can be assessed using the threshold value.Variable selection and data sourceThe super-efficiency SBM model was used in calculating the agricultural water use efficiency for the 31 provinces in China from 2007 to 2018. The input–output indicators were defined before the calculations, as shown in Extended Data Table 1.The selection of input–output factors to measure the utilization efficiency of agricultural water resources follows the principles of availability and operability. The input variables included: (1) agricultural water consumption, (2) the number of employees in agriculture, forestry, animal husbandry, and fishery, (3) the total power of agricultural machinery, and (4) the expenditure of local finance on agriculture, forestry, and water affairs. In terms of output, the added value in agriculture, forestry, animal husbandry, and fishery (based on 2007) was used as the expected output, while ammonia nitrogen emission, agricultural chemical oxygen demand emission, and agricultural carbon emission comprised the unexpected output.This study considered the scale of carbon emissions released by the agricultural system. According to existing research, agricultural carbon emissions are associated with rural environmental pollution35. The main consequence of agricultural pollutant emissions is soil pollution, which leads to rural groundwater pollution36,37,39,40,41,41. The deterioration of groundwater quality adversely affects the development of the agricultural economy and threatens the safety of the drinking water supply for rural residents.The threshold regression model was used to investigate the convergence of agricultural water use efficiency and observe the changes in agricultural water consumption under different influencing factors. The control variables include the following: water resource endowment, the number of agricultural labor, the income level of rural residents, industrial structure, the degree of government’s attention, the degree of dependence on foreign trade, and the level of economic development, as shown in Extended Data Table 2. For water resource endowment (WR), WR is expressed in per capita water resource (m3 / person). Zhang Lixiao45,46 and previous studies have shown a negative correlation between water resource endowment and water resource utilization. For agricultural labor (ah), the variable is expressed by the number of people engaged in agriculture, forestry, animal husbandry, and fishery (10,000 people). Past studies suggest rural population affects the consumption of agricultural water resources47,50,51,52,53,52. For income levels, rural residents’ income level is indicated by the per capita disposable income of rural households. Wang Xueyuan et al.3 and Han Qing et al.53 argue that the increase in the rural residents’ income would limit agricultural water consumption. For industrial structure (× 2), which is expressed by the proportion of industrial added value in GDP, research has shown water resource efficiency would vary under different industrial structures54,57,56. For the government’s attention degree (GA), the variable is expressed by the proportion of agriculture, water affairs, and forestry spending in the total financial expenditure. The government’s support for comprehensive agricultural development and infrastructure and technology upgrading for agricultural, forestry, and water conservation significantly affects water resource utilization efficiency16,56,59,58. For the degree of dependence on foreign trade (open), the parameter is indicated by the proportion of the total import and export of agricultural and sideline products in the GDP. Changes in import demand can reduce or increase the consumption and pollution of water resources. Likewise, export demand changes, especially in high water-consuming and high polluting products, can significantly improve or degrade water resource efficiency. And for the level of economic development, expressed in terms of GDP, the level of regional economic development plays a positive role in promoting the efficiency of water resource utilization59,62,61. More

1.Turner, I. M. Species loss in fragments of tropical rain forest: a review of the evidence. J. Appl. Ecol. 33, 200–209 (1996).Article 

Google Scholar 
2.Pimm, S. L. & Raven, P. Biodiversity: extinction by numbers. Nature 403, 843 (2000).ADS 
CAS 
PubMed 
Article 

Google Scholar 
3.Robinson, W. D. et al. Distribution of bird diversity in a vulnerable Neotropical landscape. Conserv. Biol. 18, 510–518 (2004).Article 

Google Scholar 
4.Rompré, G., Robinson, W. D. & Desrochers, A. Causes of habitat loss in a Neotropical landscape: The Panama Canal corridor. Landsc. Urban Plan. 87, 129–139 (2008).Article 

Google Scholar 
5.Diamond, J. Dammed experiments. Science 294, 1847 (2001).CAS 
PubMed 
Article 

Google Scholar 
6.Şekercioḡlu, Ç. H. et al. Disappearance of insectivorous birds from tropical forest fragments. Proc. Natl. Acad. Sci. USA 99, 263 (2002).ADS 
PubMed 
Article 
CAS 

Google Scholar 
7.Henle, K., Davies, K. F., Kleyer, M., Margules, C. & Settele, J. Predictors of species sensitivity to fragmentation. Biodivers. Conserv. 13, 207–251 (2004).Article 

Google Scholar 
8.Stratford, J. A. & Robinson, W. D. Gulliver travels to the fragmented tropics: geographic variation in mechanisms of avian extinction. Front. Ecol. Environ. 3, 85–92 (2005).Article 

Google Scholar 
9.Robinson, W. D. & Sherry, T. W. Mechanisms of avian population decline and species loss in tropical forest fragments. J. Ornithol. 153, 141–152 (2012).Article 

Google Scholar 
10.Terborgh, J. Preservation of natural diversity: the problem of extinction prone species. Bioscience 24, 715–722 (1974).Article 

Google Scholar 
11.Karr, J. R. Population variability and extinction in the avifauna of a tropical land bridge island. Ecology 63, 1975–1978 (1982).Article 

Google Scholar 
12.Sieving, K. E. Nest predation and differential insular extinction among selected forest birds of central Panama. Ecology 73, 2310–2328 (1992).Article 

Google Scholar 
13.Bierregaard, R. O., Lovejoy, T. E., Kapos, V., dos Santos, A. A. & Hutchings, R. W. The biological dynamics of tropical rainforest fragments. Bioscience 42, 859–866 (1992).Article 

Google Scholar 
14.Laurance, W. F. Forest-climate interactions in fragmented tropical landscapes. Philos. Trans. Roy. Soc. Lond. B: Biol. Sci. 359, 345–352 (2004).Article 

Google Scholar 
15.Laurance, W. F. & Curran, T. J. Impacts of wind disturbance on fragmented tropical forests: a review and synthesis. Austral. Ecol. 33, 399–408 (2008).Article 

Google Scholar 
16.Stratford, J. A. & Stouffer, P. C. Forest fragmentation alters microhabitat availability for Neotropical terrestrial insectivorous birds. Biol. Conserv. 188, 109–115 (2015).Article 

Google Scholar 
17.Patten, M. A. & Smith-Patten, B. D. Testing the microclimate hypothesis: light environment and population trends of Neotropical birds. Biol. Conserv. 155, 85–93 (2012).Article 

Google Scholar 
18.Ausprey, I. J., Newell, F. L. & Robinson, S. K. Adaptations to light predict the foraging niche and disassembly of avian communities in tropical countrysides. Ecology 102, e03213 (2021).PubMed 
Article 

Google Scholar 
19.Busch, D. S., Robinson, W. D., Robinson, T. R. & Wingfield, J. C. Influence of proximity to a geographical range limit on the physiology of a tropical bird. J. Anim. Ecol. 80, 640–649 (2011).PubMed 
Article 

Google Scholar 
20.Stouffer, P. C. & Bierregaard, R. O. Use of Amazonian forest fragments by understory insectivorous birds. Ecology 76, 2429–2445 (1995).Article 

Google Scholar 
21.Ferraz, G. et al. Rates of species loss from Amazonian forest fragments. Proc. Natl. Acad. Sci. 100, 14069–14073 (2003).ADS 
CAS 
PubMed 
Article 

Google Scholar 
22.Brooks, T. M., Pimm, S. L. & Oyugi, J. O. Time lag between deforestation and bird extinction in tropical forest fragments. Conserv. Biol. 13, 1140–1150 (1999).Article 

Google Scholar 
23.Kuussaari, M. et al. Extinction debt: a challenge for biodiversity conservation. Trends Ecol. Evol. 24, 564–571 (2009).PubMed 
Article 

Google Scholar 
24.Ewers, R. M. & Didham, R. K. Confounding factors in the detection of species responses to habitat fragmentation. Biol. Rev. 81, 117–142 (2006).PubMed 
Article 

Google Scholar 
25.Kattan, G. H., Alvarez-López, H. & Giraldo, M. Forest fragmentation and bird extinctions: San Antonio eighty years later. Conserv. Biol. 8, 138–146 (1994).Article 

Google Scholar 
26.Christiansen, M. B. & Pitter, E. Species loss in a forest bird community near Lagoa Santa in southeastern Brazil. Biol. Conserv. 80, 23–32 (1997).Article 

Google Scholar 
27.Laurance, W. F. et al. Ecosystem decay of Amazonian forest fragments: a 22-year investigation. Conserv. Biol. 16, 605–618 (2002).Article 

Google Scholar 
28.Sigel, B. J., Sherry, T. W. & Young, B. E. Avian community response to lowland tropical rainforest isolation: 40 years of change at La Selva biological station, Costa Rica. Conserv. Biol. 20, 111–121 (2006).PubMed 
Article 

Google Scholar 
29.Stouffer, P. C., Bierregaard, R. O., Strong, C. & Lovejoy, T. E. Long-term landscape change and bird abundance in amazonian rainforest fragments. Conserv. Biol. 20, 1212–1223 (2006).PubMed 
Article 

Google Scholar 
30.Moura, N. G. et al. Two hundred years of local avian extinctions in Eastern Amazonia. Conserv. Biol. 28, 1271–1281 (2014).PubMed 
Article 

Google Scholar 
31.Haddad, N. M. et al. Habitat fragmentation and its lasting impact on Earth’s ecosystems. Sci. Adv. 1, e1500052 (2015).ADS 
PubMed 
PubMed Central 
Article 

Google Scholar 
32.Foster, R. B. & Brokaw, N. V. Structure and History of the Vegetation of Barro Colorado Island (1982).33.Leigh, E. G. Tropical Forest Ecology: A View from Barro Colorado Island (Oxford University Press, 1999).
Google Scholar 
34.Panama Canal Authority (ACP). Meteorology and Hydrology Branch. http://www.pancanal.com (2016).35.ANAM. Informe Final de Resultados de la Cobertura Boscosa y uso del Suelo de la Republica de Panamá 1992–2000 (La Autoridad Nacional para el Ambiente (ANAM) y The International Tropical Timber Organization Panamá, 2003).36.Paton, S. 2017 Meterological and Hydrological Summary for Barro Colorado Island (2018).37.Rompré, G., Robinson, W. D., Desrochers, A. & Angehr, G. Environmental correlates of avian diversity in lowland Panama rain forests. J. Biogeogr. 34, 802–815 (2007).Article 

Google Scholar 
38.Karr, J. R. Avian extinction on Barro Colorado island, Panama: a reassessment. Am. Nat. 119, 220–239 (1982).Article 

Google Scholar 
39.Willis, E. O. Populations and local extinctions of birds on Barro Colorado Island, Panama. Ecol. Monogr. 44, 153–169 (1974).Article 

Google Scholar 
40.Robinson, W. D. Long-term changes in the avifauna of Barro Colorado Island, Panama, a tropical forest isolate. Conserv. Biol. 13, 85–97 (1999).Article 

Google Scholar 
41.Robinson, W. D., Robinson, T. R., Robinson, S. K. & Brawn, J. D. Nesting success of understory forest birds in central Panama. J. Avian Biol. 31, 151–164 (2000).Article 

Google Scholar 
42.Robinson, W. D. & Robinson, T. R. Observations of predation events at bird nests in central Panama. J. Field Ornithol. 72, 43–48 (2001).Article 

Google Scholar 
43.Robinson, W. D., Rompré, G. & Robinson, T. R. Videography of Panama bird nests shows snakes are principal predators. Ornitol. Neotrop. 16, 187–195 (2005).
Google Scholar 
44.Chapman, F. M. My Tropical Air Castle (D. Appleton and Co., 1929).
Google Scholar 
45.Chapman, F. M. Life in an Air Castle: Nature Studies in the Tropics (D. Appleton-Century Company, Incorporated, 1938).
Google Scholar 
46.Eisenmann, E. Annotated List of Birds of Barro Colorado Island, Panama Canal Zone Vol. 117 (Smithsonian Institution, 1952).
Google Scholar 
47.Willis, E. O. & Eisenmann, E. A revised list of birds of Barro Colorado Island, Panamá. Smithson. Contrib. Zool. https://doi.org/10.5479/si.00810282.291 (1979).Article 

Google Scholar 
48.Robinson, W. D. Changes in abundance of birds in a Neotropical forest fragment over 25 years: a review. Anim. Biodivers. Conserv. 24, 51–65 (2001).
Google Scholar 
49.Robinson, W. D., Brawn, J. D. & Robinson, S. K. Forest bird community structure in central Panama: influence of spatial scale and biogeography. Ecol. Monogr. 70, 209–235 (2000).Article 

Google Scholar 
50.Sodhi, N. S., Liow, L. H. & Bazzaz, F. A. Avian extinctions from tropical and subtropical forests. Annu. Rev. Ecol. Evol. Syst. 35, 323–345 (2004).Article 

Google Scholar 
51.Dunning, J. B. Jr. CRC Handbook of Avian Body Masses (CRC Press, 2007).Book 

Google Scholar 
52.R Core Team. R: A Language and Environment for Statistical Computing (R Foundation for Statistical Computing, 2018).
Google Scholar 
53.McCune, B. & Mefford, M. J. PC-ORD. Multivariate Analysis of Ecological Data (MjM Software, 2011).
Google Scholar 
54.Kursa, M. B. & Rudnicki, W. R. Feature selection with the Boruta package. J. Stat. Softw. 36, 1–13 (2010).Article 

Google Scholar 
55.Chevan, A. & Sutherland, M. Hierarchical partitioning. Am. Stat. 45, 90–96 (1991).
Google Scholar 
56.Navarrete, C. B. & Soares, F. C. dominanceanalysis: Dominance Analysis. R package version 1.0.0. (2019).57.McFadden, D. Conditional Logit Analysis of Qualitative Choice Behavior (1973).58.Menard, S. Coefficients of determination for multiple logistic regression analysis. Am. Stat. 54, 17–24 (2000).
Google Scholar 
59.McFadden, D. Quantitative methods for analyzing travel behaviour of individuals: some recent developments, Cowles Foundation Discussion Papers No. 474 (Cowles Foundation for Research in Economics, Yale University, 1977).60.Clark, W. A. & Hosking, P. L. Statistical Methods for Geographers. (1986).61.Walsh, C. & MacNally, R. Hier.Part: Hierarchical Partitioning. R package version 1.0-4. (2013).62.Harrell Jr, F. E. RMS: Regression Modeling Strategies. R package version 5.1-3. City (2019).63.Le Cessie, S. & Van Houwelingen, J. C. A goodness-of-fit test for binary regression models, based on smoothing methods. Biometrics 47, 1267–1282 (1991).MATH 
Article 

Google Scholar 
64.Guisan, A. & Zimmermann, N. E. Predictive habitat distribution models in ecology. Ecol. Model. 135, 147–186 (2000).Article 

Google Scholar 
65.Barbosa, A. M., Brown, J. A., Jimenez-Valverde, A. & Real, R. modEvA: Model Evaluation and Analysis. R package version 1.3.2. (2016).66.Suzuki, R., Shimodaira, H., Suzuki, M. R. & Suggests, M. Package ‘pvclust’. R Top. Doc. 14, 1540–1542 (2015).
Google Scholar 
67.Suzuki, R. & Shimodaira, H. Pvclust: an R package for assessing the uncertainty in hierarchical clustering. Bioinformatics 22, 1540–1542 (2006).CAS 
PubMed 
Article 

Google Scholar 
68.Anderson, M. J. A new method for non-parametric multivariate analysis of variance. Austral. Ecol. 26, 32–46 (2001).
Google Scholar 
69.Anderson, M. J. Distance-based tests for homogeneity of multivariate dispersions. Biometrics 62, 245–253 (2006).MathSciNet 
PubMed 
MATH 
Article 

Google Scholar 
70.Oksanen, J. et al. Vegan: Community Ecology Package (2013).71.Moore, R. P. Biogeographic and Experimental Evidence for Local Scale Dispersal Limitation in Central Panamanian Forest Birds (Oregon State University, 2005).
Google Scholar 
72.Lande, R. Risks of population extinction from demographic and environmental stochasticity and random catastrophes. Am. Nat. 142, 911–927 (1993).PubMed 
Article 

Google Scholar 
73.Moore, R. P., Robinson, W. D., Lovette, I. J. & Robinson, T. R. Experimental evidence for extreme dispersal limitation in tropical forest birds. Ecol. Lett. 11, 960–968 (2008).CAS 
PubMed 
Article 

Google Scholar 
74.Asquith, N. M. & Mejía-Chang, M. Mammals, edge effects, and the loss of tropical forest diversity. Ecology 86, 379–390 (2005).Article 

Google Scholar 
75.Wolda, H. Trends in abundance of tropical forest insects. Oecologia 89, 47–52 (1992).ADS 
PubMed 
Article 

Google Scholar 
76.Franks, N. R. A new method for censusing animal populations: the number of Eciton burchelli army ant colonies on Barro Colorado Island, Panama. Oecologia 52, 266–268 (1982).ADS 
PubMed 
Article 

Google Scholar 
77.Socolar, J. B. & Wilcove, D. S. Forest-type specialization strongly predicts avian responses to tropical agriculture. Proc. R. Soc. B 286, 20191724 (2019).PubMed 
Article 

Google Scholar 
78.Şekercioğlu, Ç. H., Primack, R. B. & Wormworth, J. The effects of climate change on tropical birds. Biol. Conserv. 148, 1–18 (2012).Article 

Google Scholar 
79.Karr, J. R. & Freemark, K. E. Habitat selection and environmental gradients: dynamics in the” stable” tropics. Ecology 64, 1481–1494 (1983).Article 

Google Scholar 
80.Ibarra-Macias, A., Robinson, W. D. & Gaines, M. S. Experimental evaluation of bird movements in a fragmented Neotropical landscape. Biol. Conserv. 144, 703–712 (2011).Article 

Google Scholar 
81.Stouffer, P. C. et al. Long-term change in the avifauna of undisturbed Amazonian rainforest: ground-foraging birds disappear and the baseline shifts. Ecol. Lett. 24, 186–195 (2021).PubMed 
Article 

Google Scholar 
82.Blake, J. G. & Loiselle, B. A. Enigmatic declines in bird numbers in lowland forest of eastern Ecuador may be a consequence of climate change. PeerJ 3, e1177 (2015).PubMed 
PubMed Central 
Article 

Google Scholar 
83.Legendre, P. & Condit, R. Spatial and temporal analysis of beta diversity in the Barro Colorado Island forest dynamics plot, Panama. For. Ecosyst. 6, 7 (2019).Article 

Google Scholar 
84.Condit, R., Pérez, R., Lao, S., Aguilar, S. & Hubbell, S. P. Demographic trends and climate over 35 years in the Barro Colorado 50 ha plot. For. Ecosyst. 4, 17 (2017).Article 

Google Scholar 
85.Aguilar, E. et al. Changes in precipitation and temperature extremes in Central America and northern South America, 1961–2003. J. Geophys. Res.: Atmos. 110, 2064–2082 (2005).Article 

Google Scholar 
86.Gonzalez, A. & Loreau, M. The causes and consequences of compensatory dynamics in ecological communities. Annu. Rev. Ecol. Evol. Syst. 40, 393–414 (2009).Article 

Google Scholar 
87.Kaspari, M. & Weiser, M. D. Ant activity along moisture gradients in a neotropical forest 1. Biotropica 32, 703–711 (2000).Article 

Google Scholar 
88.Wall, D. H. et al. Global decomposition experiment shows soil animal impacts on decomposition are climate-dependent. Glob. Change Biol. 14, 2661–2677 (2008).ADS 
Article 

Google Scholar 
89.Levings, S. C. & Windsor, D. M. Litter moisture content as a determinant of litter arthropod distribution and abundance during the dry season on Barro Colorado Island, Panama. Biotropica 16, 125–131 (1984).Article 

Google Scholar 
90.Brawn, J. D., Benson, T. J., Stager, M., Sly, N. D. & Tarwater, C. E. Impacts of changing rainfall regime on the demography of tropical birds. Nat. Clim. Chang. 7, 133 (2017).ADS 
Article 

Google Scholar 
91.Karp, D. S. et al. Agriculture erases climate-driven β-diversity in Neotropical bird communities. Glob. Change Biol. 24, 338–349 (2018).ADS 
Article 

Google Scholar 
92.Frishkoff, L. O. et al. Loss of avian phylogenetic diversity in neotropical agricultural systems. Science 345, 1343 (2014).ADS 
CAS 
PubMed 
Article 

Google Scholar 
93.Wright, S. J. How isolation affects rates of turnover of species on islands. Oikos 44, 331–340 (1985).Article 

Google Scholar 
94.Chadwick, R., Good, P., Martin, G. & Rowell, D. P. Large rainfall changes consistently projected over substantial areas of tropical land. Nat. Clim. Chang. 6, 177–181 (2016).ADS 
Article 

Google Scholar 
95.Esquivel-Muelbert, A. et al. Compositional response of Amazon forests to climate change. Glob. Change Biol. 25, 39–56 (2019).ADS 
Article 

Google Scholar  More