Reducing herbivory in mixed planting by genomic prediction of neighbor effects in the field

Genetic diversity is a key facet of biodiversity that can alter ecosystem services, including plant resistance to herbivores. Mixed stands of plant genotypes can lead to either associational resistance or associational susceptibility, depending on neighbor interactions mediated by plant–herbivore, plant–plant, and plant–carnivore processes. While intraspecific diversity can enhance resistance, outcomes are inconsistent and beneficial genotype combinations are difficult to identify at scale. Prior work has mainly focused on limited pairwise combinations or growth-related traits, constrained by the combinatorial effort needed to control neighbor composition. This study addresses whether genome-wide polymorphism data combined with a Neighbor GWAS framework can predict genotype pairs that reduce herbivory in the field. The authors test this by planting 199 A. thaliana genotypes in randomized mixtures across two field sites (Zurich, Switzerland; Otsu, Japan), quantifying herbivore damage and insect communities, modeling neighbor genotype effects, and validating predicted beneficial pairs in the field.

The paper builds on literature documenting associational resistance and susceptibility in plant–herbivore systems and the multifaceted ecological mechanisms (volatile cues, apparency, physical and chemical defenses, and carnivore attraction) underlying these effects. Previous GWAS in A. thaliana have dissected stand-level growth and neighbor effects in controlled settings but not large field-scale anti-herbivore resistance due to logistical constraints. Known single-gene effects such as GLABRA1 (GL1) influencing flea beetle resistance and glucosinolate loci (e.g., GS-OH, MAM1) have been reported, though glucosinolate effects may be attenuated with specialist herbivores. The study leverages and extends Neighbor GWAS, inspired by the Ising model, to analyze spatially structured traits and infer positive or negative neighbor effects at the SNP level.

- Study system and design: 199 selfed A. thaliana accessions genotyped via RegMap and 1001 Genomes were planted in randomized block designs at two field sites: Zurich (native range) and Otsu (exotic), over 2017–2018 for GWAS and 2019 for validation. Each site had eight 200-plant blocks (199 accessions + Col-0) in 2017–2018 and three blocks in 2019, with blocks >2 m apart. Plants were pre-grown under controlled conditions and transplanted to pots arranged in checkered layouts. - Phenotyping: Herbivore damage was measured as number of flea beetle leaf holes (Zurich) or categorical leaf area loss score (Otsu). Insects (18 species) were visually counted every 2–3 days; extended phenotypes included numbers of external chewers, other herbivores, and total insect species per plant. Initial plant size and presence/absence of inflorescences were recorded as covariates. Data were log-transformed where appropriate. - Standard GWAS: A linear mixed model framework assessed focal genotype effects and confirmed heritability. Known GL1 glabrous mutants were included; GL1-associated SNPs significantly affected flea beetle damage in Zurich after Bonferroni correction. - Neighbor GWAS: A mixed model including focal genotype (β1) and neighbor genotype effects (β2) was fitted. Phenotypic variation explained (PVE) by neighbor genotypes was tested via stepwise likelihood ratio tests. Manhattan plots mapped β2 association across the genome, correcting for relatedness with weighted kinship matrices and Bonferroni thresholds; permutation tests (neighbor shuffling and genome rotation) assessed significance. - Genomic prediction (LASSO): To capture polygenic neighbor effects, the Neighbor GWAS was recast as a multiple regression (y = X0β0 + X1β1 + X2β2 + e). SNPs (MAF>0.05) were LD-pruned to 1,242,128 markers, and the same number of neighbor-related features were included. Models were trained on 2017–2018 and validated on 2019 using Spearman correlations between predicted and observed residual phenotypes (controlling covariates). The neighbor-including LASSO that maximized predictability provided averaged β2 estimates across regularization values, yielding 756 SNPs with non-zero neighbor effects for Zurich herbivory. - In silico pairwise effect sizes: Using the β2 vector, pairwise mixture vs monoculture effect sizes for all 19,701 genotype pairs were calculated as [xi∘xj]β2 - [xi∘xi]β2 to predict associational resistance (>0) or susceptibility (<0). Simulations on lattices examined how genotype richness affects predicted damage under random or positive-pair selection. - GO enrichment: Candidate genes within 10 kb of LASSO-selected SNPs were tested for GO enrichment (Fisher’s exact, FDR correction), summarized with REVIGO, and evaluated separately for positive and negative β2 sets. - Field validation of predicted pairs: Three pairs covering a gradient of positive predicted effect sizes were planted in Zurich in 2019 and 2021 under monoculture and mixture, and in large (8×8) and small (4×4) patches: Bg-2/Uod-1 (large positive), Vastervik/Jm-0 (moderate positive), Bla-1/Bro1-6 (slight positive). Leaf holes were counted after three weeks. Linear mixed models with log-transformed counts assessed effects; marginal means compared mixture vs monoculture per accession (Sidak correction). - Laboratory choice experiments: Black flea beetle Phyllotreta astrachanica adults were identified via COI barcoding and used in choice arenas with the same pairs. Negative binomial GLMs tested accession preferences. - Software: rNeighborGWAS, gaston, glmnet (Python), vegan, lme4, lmerTest, emmeans, and custom scripts (GitHub/Zenodo).

- Neighbor genotype effects contribute to herbivory and insect community variation: Including neighbor genotypes significantly increased PVE for herbivore damage (likelihood ratio χ² > 7.4, d.f. 1, p < 0.01) at both sites; permutation tests supported non-random neighbor effects in Zurich (p < 0.05). - Mobility-dependent neighbor influence: In Zurich, external chewers (mostly jumping flea beetles) were significantly influenced by neighbor genotypes (χ² = 6.0, d.f. 1, p < 0.05). In Otsu, other herbivores (including mobile thrips) showed significant neighbor effects (χ² = 5.6, d.f. 1, p < 0.05). Total insect species was neighbor-influenced at both sites (χ² > 4.8, p < 0.05). - Polygenic basis: GWAS of β2 detected no genome-wide significant SNPs after Bonferroni correction, consistent with polygenic neighbor effects. - Asymmetry in β2 sign and allele frequency: Among top 0.1% SNPs, both positive (associational resistance) and negative (associational susceptibility) β2 were present; positive β2 SNPs tended to have higher minor allele frequencies and showed more signatures of balancing selection than negative β2 SNPs. - Improved prediction with neighbor-inclusive models: For Zurich herbivory, neighbor-including LASSO slightly outperformed neighbor-excluding LASSO on the 2019 test set (max Spearman ρ = 0.416 vs 0.391). This yielded 756 SNPs with non-zero neighbor β2. - Mechanistic insight: GO enrichment for SNPs with positive β2 highlighted jasmonic acid biosynthetic process, including LOX2/LOX6, consistent with volatile-mediated defenses and potential neighbor repellency/communication; such enrichment was absent for negative β2. - Predicted pair effects and simulations: Distribution of pairwise mixture effect sizes skewed negative, indicating prevalent associational susceptibility. Pairwise effect size was negatively correlated with monoculture damage (Spearman ρ = -0.38, p = 3.4×10⁻⁸). Simulations suggested increasing genotype richness generally increased predicted herbivore damage, unless carefully selecting positive pairs. - Field validation: Of 19,701 pairs, 823 were predicted to reduce herbivory. Three selected pairs were tested in Zurich: Bg-2/Uod-1 (large positive), Vastervik/Jm-0 (moderate), Bla-1/Bro1-6 (slight). Mixtures of Bg-2/Uod-1 and Vastervik/Jm-0 showed significant reductions in flea beetle damage relative to monocultures (average reductions within pairs ~24.8% and ~22.7%, respectively; marginal means t > 2.02, d.f. 125, adjusted p < 0.05). Bla-1/Bro1-6 showed no significant reduction. Effects were generally stronger in larger patches. - Lab choice assays: Flea beetles showed significant preferences consistent with field results (Bg-2 vs Uod-1: χ² = 13.35, p < 0.001; Vastervik vs Jm-0: χ² = 5.71, p < 0.05; Bla-1 vs Bro1-6: χ² = 0.87, p = 0.35).

The study demonstrates that genome-wide neighbor effects shape herbivore damage and certain aspects of insect community composition, particularly for mobile herbivores. Despite the prevalence of associational susceptibility among genotype mixtures, integrating Neighbor GWAS with genomic prediction allowed identification of beneficial genotype pairs that confer associational resistance. The lack of significant β2 SNPs in GWAS coupled with improved predictive performance using neighbor-inclusive LASSO supports a polygenic architecture for neighbor interactions. Functional enrichment near positive β2 loci implicates jasmonate pathway genes (e.g., LOX2/LOX6), suggesting volatile-mediated repellency or plant–plant signaling as plausible mechanisms behind associational resistance in this system, whereas carnivore-mediated effects appear less likely. Field validations translated predictions into measurable reductions of herbivore damage in selected pairs, indicating practical potential for designing resistant genotype mixtures. The approach extends genomic selection to include neighbor effects, enabling prediction and deployment of genotype combinations with beneficial biodiversity effects in ecologically realistic field conditions.

By combining field experiments on 199 Arabidopsis genotypes with Neighbor GWAS and LASSO-based genomic prediction, the study identified and validated genotype pairs that reduce herbivory in mixed plantings, despite widespread associational susceptibility. The work advances understanding of polygenic neighbor effects on herbivory, highlights potential roles of jasmonate-mediated volatile defenses, and provides a predictive framework to assemble beneficial genotype mixtures. This strategy can inform pest management and biodiversity-based agriculture by targeting specific genotype pairs that enhance resistance without complicating agronomy (e.g., compatible phenology). Future research should dissect causal genes and interactions through multiplex genome editing and expand predictions to longer-lived plant systems and more complex communities to generalize neighbor effect mechanisms and optimize mixture designs.

- No genome-wide significant SNPs for neighbor effects were detected, indicating a highly polygenic architecture and limiting locus-specific mechanistic inference. - Predictive gains from neighbor-inclusive LASSO were modest; community composition phenotypes lacked improved predictability, likely due to inter-annual shifts in dominant species. - Validation tested only three genotype pairs and at one site for field reduction effects; broader validation across environments, herbivore communities, and management conditions is needed. - Effects of mixtures diminished with increasing genotype richness due to confounding negative pairwise interactions; scaling from pairs to multi-genotype stands remains challenging. - Mechanistic attribution (e.g., volatile signaling vs other pathways) remains inferential; direct functional validation of implicated genes and pathways is required. - The Otsu site showed weaker neighbor effects for sedentary external chewers, suggesting context dependence tied to herbivore mobility and life history.