Skip to main content

Genetic analysis reveals three novel QTLs underpinning a butterfly egg-induced hypersensitive response-like cell death in Brassica rapa



Cabbage white butterflies (Pieris spp.) can be severe pests of Brassica crops such as Chinese cabbage, Pak choi (Brassica rapa) or cabbages (B. oleracea). Eggs of Pieris spp. can induce a hypersensitive response-like (HR-like) cell death which reduces egg survival in the wild black mustard (B. nigra). Unravelling the genetic basis of this egg-killing trait in Brassica crops could improve crop resistance to herbivory, reducing major crop losses and pesticides use. Here we investigated the genetic architecture of a HR-like cell death induced by P. brassicae eggs in B. rapa.


A germplasm screening of 56 B. rapa accessions, representing the genetic and geographical diversity of a B. rapa core collection, showed phenotypic variation for cell death. An image-based phenotyping protocol was developed to accurately measure size of HR-like cell death and was then used to identify two accessions that consistently showed weak (R-o-18) or strong cell death response (L58). Screening of 160 RILs derived from these two accessions resulted in three novel QTLs for Pieris brassicae-induced cell death on chromosomes A02 (Pbc1), A03 (Pbc2), and A06 (Pbc3). The three QTLs Pbc1–3 contain cell surface receptors, intracellular receptors and other genes involved in plant immunity processes, such as ROS accumulation and cell death formation. Synteny analysis with A. thaliana suggested that Pbc1 and Pbc2 are novel QTLs associated with this trait, while Pbc3 also contains an ortholog of LecRK-I.1, a gene of A. thaliana previously associated with cell death induced by a P. brassicae egg extract.


This study provides the first genomic regions associated with the Pieris egg-induced HR-like cell death in a Brassica crop species. It is a step closer towards unravelling the genetic basis of an egg-killing crop resistance trait, paving the way for breeders to further fine-map and validate candidate genes.

Peer Review reports


Plant-insect interactions often start with herbivore egg deposition on plant tissues. Through millions of years of co-evolution with insects, plants have evolved mechanisms to recognize insect eggs as non-self to induce defence responses [1]. Different egg-killing traits have been described, such as neoplasm formation [2,3,4], secretion of toxic chemicals [5], tissue crushing [6], and hypersensitive response (HR)-like cell death [7,8,9]. Such defenses represent an additional component in the plant-insect arms race, but their potential for sustainable crop protection has so far largely been overlooked and under-utilized [10].

Eggs of Pieris spp. butterflies (Lepidoptera: Pieridae) induce a HR-like cell death on their natural host plants belonging to the Brassicaceae family [11]. The large cabbage white (P. brassicae L.) and the small cabbage white (P. rapae L.) represent major pests of Brassica crops worldwide [12, 13]. Despite available knowledge on plant defences against them, Pieris spp. are specialists well equipped for feeding on Brassicaceae [14,15,16,17]. Their caterpillars effectively detoxify the secondary metabolites produced by their host plants, so-called glucosinolates or “mustard oils” [18, 19]. Considering this, Pieris egg-induced HR-like cell death may represent a genuine and unexplored defense response of Brassicaceae against adapted specialist herbivores.

The HR-like cell death response was initially observed underneath eggs of P. rapae and P. napi deposited on leaves of wild populations of black mustard (Brassica nigra L.), on which it caused egg-killing by desiccating or dropping off [7]. Later, this egg-killing trait was also observed underneath eggs of P. brassicae and, interestingly, it was found to work in concert with the attraction of egg parasitoid wasps through the release of oviposition-induced plant volatiles [20, 21]. Under field conditions, the synergistic effect of HR-like cell death and egg parasitism reduced up to 80% of Pieris egg survival on B. nigra [21]. The direct egg-killing effect of the cell death seems to work mainly against singly laid eggs, irrespective of whether they are from solitary species such as P. rapae and P. napi, or from the gregarious species P. brassicae [11, 22]. These studies suggest that HR-like cell death may be an effective egg-killing trait, for which the plant molecular and genetic mechanisms are still poorly understood.

Most of the initial knowledge on the molecular aspects of the Pieris egg-plant interaction was obtained using the model plant A. thaliana [23]. Upon P. brassicae oviposition, A. thaliana responds with reactive oxygen species (ROS) and cell death, further accompanied by the accumulation of phytohormone salicylic acid (SA) and the induction of the SA-responsive gene PATHOGENESIS-RELATED PROTEIN 1 (PR1) [24, 25]. These responses were shown to be dependent on ENHANCED DISEASE SUSCEPTIBILITY 1 (EDS1), ISOCHORISMATE SYNTHASE 1/SALICYLIC ACID INDUCTION DEFICIENT 2 (ICS1/SID2) and, partially, NONEXPRESSER OF PR GENES 1 (NPR1), which are known signaling components of plant defense responses against biotrophic pathogens [26]. Transcriptomic studies in different plants species have confirmed that insect oviposition induces genes associated with SA- and ROS-mediated immune responses and PR1 gene expression [15, 24, 27,28,29,30,31,32,33], including in B. nigra and B. rapa [21, 34, 35]. Further, it  has been suggested that there is a conserved transcriptional response amongst different plant-insect egg interactions [36].

The similarities between the plant defenses induced against insect eggs and biotrophic pathogens suggest that insect eggs are also recognized by the plant immune system [23, 37], but it is yet not known how. The induction of plant defenses partly relies on the specific recognition of non-self molecules released by biotic attackers that are detected by plasma membrane pattern recognition receptors (PRRs) [38, 39] or intracellular nucleotide-binding leucin-rich repeat receptors (NLRs) [40]. Feeding of herbivorous insects induce plant immunity through the release of herbivore-associated molecular patterns (HAMPs) contained in oral secretions of insect larvae and/or damage-associated molecular patterns (DAMPs) resulting from damaged plant tissues [41]. Both signals have been associated with the perception by different types of PRRs [16, 42].

Contrary to cues of larval feeding, only a few insect egg-associated molecular patterns (EAMPs) have been identified [2, 27, 43, 44]. In A. thaliana, some candidate PRRs involved in perception of P. brassicae eggs were recently discovered. Several L-type lectin receptor-like kinases (LecRKs), a class of PRRs, were upregulated upon P. brassicae oviposition [24]. One of them, LecRK-I.8, was found to be required for the induction of downstream ROS production, cell death and PR1 expression [26, 45]. More recently, a genome-wide association study (GWAS) in A. thaliana identified LecRK-I.1 as a candidate gene underlying one of two loci involved in the induction of cell death upon treatment with P. brassicae egg extract [46].

To date, only a few studies have attempted to map genetic loci associated with insect oviposition-induced responses [46,47,48,49]. A strong HR-like cell death that eventually leads to egg-killing has been mainly shown for plant species of the tribe Brassiceae (Lineage II), which includes wild species such as Brassica nigra, Sinapis spp., Crambe spp., as well as diverse Brassica crops such as B. napus, B. oleracea and B. rapa but not A. thaliana [11, 21, 50]. Interestingly, species belonging to the tribe Brassiceae are known host plants for Pieris spp. while A. thaliana is not [51].

Next to interspecific variation between Brassicaceae species we also identified intraspecific variation in HR-like cell death among accessions of several species [11, 21], suggesting that genetic analysis to identify casual loci should be feasible. Up to now, classical forward genetics, such a linkage mapping and/or GWAS, helped to identify quantitative trait loci (QTLs) involved in both upstream (perception) and downstream mechanisms associated with plant resistance to insect feeding [15, 52,53,54,55,56]. Currently, genetic mapping efforts in Brassica crop species are increasingly made possible given a growing availability of high quality genomic and genetic resources [57,58,59,60]. It is thus timely to use genetic approaches to unravel the genetics underlying plant-insect egg interactions in non-model species. Genetic mapping of insect egg-induced defenses in Brassica crops can help both the fundamental understanding of HR-like cell death and its applied use as novel defense trait in plant breeding.

Here we present the genetic analysis of P. brassicae butterfly egg-induced HR-like cell death in Brassica rapa by QTL mapping. First, we investigated the phenotypic variation for HR-like cell death within B. rapa germplasm using a core collection previously assembled and curated [61, 62]. Then, we assessed the robustness of the phenotype and we quantitatively measured cell death size with a novel image-based phenotyping protocol. We identified two accessions with a significant difference in size of HR-like cell death and we screened a recombinant inbred line (RIL) population resulting in the identification of three novel QTLs. This study provides the first QTLs and candidate genes associated with butterfly egg-induced cell death in B. rapa, an important crop species and natural host plant of Pieris spp.


Screening of a B. rapa core collection

As a first objective, we investigated whether there was intraspecific variation for HR-like cell death in our B. rapa core collection. Out of the whole collection, we screened a subset of 56 accessions representing geographical and morphological diversity within the collection (Additional file 2: Supplementary Table S1). Plants were treated with P. brassicae egg wash, which was previously reported as a reliable egg-mimicking treatment in B. nigra [35]. The main phenotypic diversity in HR-like response among the B. rapa accessions was limited to variation in cell death size (Additional file 1: Supplementary Fig. S1, Additional file 2: Supplementary Table S2). Egg wash induced a cell death on most of the accessions which appeared as necrotic black/dark spots of varying size on the leaf abaxial side (score 1–2). However, the spots never developed into the fully expanded and brown necrotic tissue, also visible on the adaxial side (score 3) (Additional file 1: Supplementary Fig. S1a). Such a strong cell death was only observed in the B. nigra accession included as positive control (Additional file 1: Supplementary Fig. S1b). Nevertheless, we found differences in HR-like cell death between B. rapa accessions (Kruskal-Wallis: χ257 = 130.59, P < 0.001). Six accessions, i.e. CC-106, FT-086, MIZ-019, R500, R-o-18 and VT-089, showed no cell death (score 0) in all biological replicates (Additional file 1: Supplementary Fig. S1b). Most of the accessions developed only a weak response, with a within-accession variation between individual plants ranging from no cell death (score 0) to small dark necrotic spots (score 1). At the other end of the phenotypic distribution, eleven accessions developed an HR-like cell death of score 2 in most of the biological replicates (i.e. ZCT, PC-184, IMB211, L58, PC-078, CC-114, CC-048, CC-168, CC-050, CC-Z16, CC-058). A specific morphotype was not associated with HR-like cell death as most of the major crop types (Pak choi, turnip, oil types) were found at both extremes of the phenotypic distribution (Additional file 1: Supplementary Fig. S1b). The only exception were the Chinese cabbage (CC) accessions, of which 8 out of 14 developed an HR-like cell death with large black/dark spots (score 2) on most of the biological replicates. Genetic heterogeneity of accessions appeared to be not associated with cell death variation as heterogeneous accessions and homogenous inbred lines and DH lines were found on both side of the phenotypic distribution.

Overall, we found statistical differences in HR-like cell death (Dunn’s test, P < 0.01) between the accessions that showed no cell death (score 0) and the accessions that developed large dark necrotic spots (score 2) upon egg wash treatment (Additional file 2: Supplementary Table S2). We then selected ten accessions either showing no response (CC-106, R-o-18, R500, SO-040), little cell death (BRO-030) or a strong cell death (score 2) in at least few replicates (BRO-127, CC-AO3, IMB211, CC-168, L58) for a further evaluation of their cell death phenotype. These accessions were selected based on specific criteria (see Material and Methods), and also because they were available as homozygous lines; being either inbred due to repeated selfing (self-compatible accessions) or previously used to generate homozygous DH lines (self-incompatible accessions).

Image-based phenotyping of HR-like cell death size on selected B. rapa homozygous lines

The selected B. rapa homozygous lines (inbred and DH lines) were re-evaluated to assess the robustness of their HR-like cell death phenotype with the aim to identify ideal parental lines to generate biparental mapping populations. Plants were treated with both P. brassicae egg clutches (10–20 eggs) and egg wash droplets to compare to which extent the egg wash could mimic the HR-like cell death induced by eggs. Further, we measured HR-like cell death size as quantitative trait using an image-based phenotyping protocol (Additional file 1: Supplementary Fig. S2-S3). Image analysis confirmed the differences in HR-like cell death between the selected accessions (Fig. 1, Additional file 2: Supplementary Table S3). Overall, we found differences in mean cell death sizes in response to both eggs (ANOVA: F74 = 8.55, P < 0.001) and egg wash (ANOVA: F74 = 15.88, P < 0.001). The two accessions that developed the smallest HR-like response (CC-106, R-o-18) were statistically different in cell death size from the ones with the largest HR-like response (IMB211, L58) for both eggs and egg wash (Tukey’s HSD, P < 0.01). Overall, accessions IMB211 and L58 showed the largest cell death size for both treatments (Fig. 1). In fact, mean cell death size induced by either eggs or egg wash were similar for IMB211 (1.20 and 1.24 mm2, respectively) and L58 (1.23 and 1.33 mm2, respectively). In contrast, accessions CC-AO3, CC-168 and SO-040 showed a cell death induced by eggs that was two to three times larger than the response induced by egg wash. To a lesser extent, R500, BRO-030 and BRO-127 also showed a higher cell death induced by eggs compared to egg wash. CC-106 and R-o-18 showed the smallest mean cell death underneath the eggs (0.07 and 0.24 mm2, respectively), in contrast to the total absence of cell death upon egg wash treatment observed in the germplasm screening. Overall, the two treatments showed limited correlation across the ten accessions (Pearson’s r = 0.55, P < 0.001), mainly because they resulted in comparable responses only for the accessions at the extremes of the distribution. The broad-sense heritability (H2) was slightly lower for cell death size induced by eggs (0.47) than by egg wash (0.64) (Additional file 2: Supplementary Table S3).

Fig. 1
figure 1

Phenotypic variation in hypersensitive response (HR)-like cell death size between ten B. rapa homozygous lines. A Cell death induced by 10–15 P. brassicae eggs. B Differential response on DH lines R-o-18 and L58 leaves underneath P. brassicae eggs. C Cell death upon spot-inoculation with 5 μl droplets of P. brassicae egg wash. D Differential response on DH lines R-o-18 and L58 leaves at egg washed-treated spots. For each accession, N = 6–10 plants were used for both experiments, each plant was treated with eggs or egg wash on two leaves. Cell death size was quantified using a custom image-based phenotyping protocol. Each plant was assigned the cell death of the most severe spot. Boxplots represents the interquartile range (1st and 3rd quantile) and the median, each dot represents a single plant. Letters report differences in mean size of HR-like cell death between accessions (Tukey’s HSD test, P < 0.01). Broad-sense heritability (H2) is indicated at top right corner of each graph. Magnification bars inside photos = 1 mm

In summary, accessions with a small/intermediate HR-like response showed a larger cell death size under eggs compared to egg wash, while the overall ranking was similar. Thus, we concluded that using egg deposition worked better than egg wash to screen for HR-like cell death in order to not underestimate the cell death induced by low responsive B. rapa accessions. Overall, IMB211 and L58 were confirmed as lines with a strong HR-like cell death while CC-106, R-o-18, R500 confirmed to be lines with a weak cell death, validating the results of the germplasm screening. Out of these accessions, L58 and R-o-18 represented ideal candidates for crossings because of their self-compatibility, similar flowering time, and leaf size/shape. Thus, we used the L58 (♀) x R-o-18 (♂) RIL population that was previously generated by Bagheri et al. [63] to pursue QTL linkage mapping.

Phenotypic analysis and QTL mapping on a RIL population

The RIL population L58 x R-o-18 consisting of 160 lines (F10) was used to identify QTLs underlying P. brassicae egg-induced HR-like cell death. We generated a new linkage map combining markers from previous studies (Additional file 2: Supplementary Tables S4-S5) [63, 64]. The final genetic map consisted of 485 loci and covered a total of 1154.4 cM, with a mean density of 2.38 cM (Additional file 1: Supplementary Fig. S4, Additional file 2: Supplementary Tables S6-S7). Image-based phenotyping of egg-induced cell death from three experiments was used to estimate best linear unbiased estimators (BLUEs) of cell death size for each parental and RIL genotype. Overall, the parents R-o-18 and L58 showed BLUE values of 0.49 (SD = 0.4) and 1.53 (SD = 0.42) mm2, respectively (Fig. 2, Table 1). Their within-accession variation in HR-like cell death size, i.e. their phenotypic range, was larger than what we observed in previous germplasm evaluations, thus resulting in a smaller difference in mean cell death size between the two parents. The RILs showed an approximate normal distribution of cell death size with a mean BLUE value of 0.77 (SD = 0.51) mm2 (Fig. 2, Table 1). The RILs phenotypic distribution was skewed towards the R-o-18 phenotypic value and only seven RILs developed a cell death size larger than L58. The broad-sense heritability across the three experiments was similar to what was previously observed for egg-induced cell death size (H2 = 0.49).

Fig. 2
figure 2

Phenotypic distribution of P. brassicae egg-induced cell death in the B. rapa RIL population L58 x R-o-18. Blue (R-o-18) and red (L58) dots indicate single plants used across three experiments (N = 7) and that were used to estimate single parental BLUE values. Green dots indicate a single BLUE value for each RIL (N = 3). Each plant was oviposited with two egg clutches and cell death size was quantified using a custom image-based phenotyping protocol. The largest cell death out of the two clutches was assigned to each plant. Boxplots represents the interquartile range (1st and 3rd quantile) and the median. Black diamonds represent mean BLUE value of the two parents and the whole RIL population. Broad-sense heritability (H2) is indicated at top right corner of the graph

Table 1 Summary statistics of cell death phenotypic data (BLUEs) of the L58 x R-o-18 RIL population

A total of three QTLs associated with HR-like cell death size were identified on three B. rapa chromosomes using an interval mapping method (Haley-Knott regression). The QTLs were named P. brassicae egg-induced cell death (Pbc) (Fig. 3a, Table 2, Additional file 1: Supplementary Fig. S5). First, phenotypic data (BLUEs) were analyzed using single-QTL models, resulting in the identification of two QTLs, i.e. Pbc1 on chromosome A02 (LOD 5.63) and Pbc3 on chromosome A06 (LOD 4.15). Additionally, multi-QTL model (MQM) mapping detected another QTL, Pbc2, on chromosome A03 (LOD 3.33). Two-QTL models revealed absence of epistatic interactions from any pairwise comparison among Pbc1–3, and weak additive interactions between Pbc1:Pbc2 and Pbc1:Pbc3 (Additional file 1: Supplementary Fig. S6). Pbc1 explained 17.9% of the additive phenotypic variance, with BrID11121 as top marker (85.4 cM) and a 1.5-LOD confidence interval spanning about 27 cM between markers 899,118|9,904,922 and BrID11907 (Table 2). The minor QTLs Pbc2 and Pbc3 explained a smaller proportion of the additive phenotypic variance, 6.35 and 7.32% respectively, with BrID90099 (129.2 cM) and BrID90095 (63.9 cM) as respective top markers (Table 2). Pbc1 was the only QTL with a 1.5-LOD confidence interval lying entirely above the LOD significance threshold (Fig. 3a). As the RIL phenotypic distribution was skewed toward the R-o-18 values, we expected L58 alleles contributing to a larger cell death size for all QTLs. Interestingly, this was true only for Pbc1 which showed opposite effect size compared to Pbc2 and Pbc3 (Fig. 3b). In fact, the allele of L58 contributed to an increase in HR-like size of 0.45 mm2 for Pbc1, while the allele of R-o-18 determined an increase in HR-like of 0.27 mm2 for Pbc2 and 0.28 mm2 Pbc3 (Table 2).

Fig. 3
figure 3

Chromosomal locations and allelic effects of three QTLs for Pbrassicae egg-induced cell death (Pbc) in B. rapa. A LOD score of chromosomes A02, A03 and A06 from MQM mapping of HR-like cell death size induced by Pieris eggs on 160 RILs. Labels indicate the closest marker to the peak LOD score. LOD threshold is indicated with a dashed horizontal line (2.59 after 1000 permutations at 5% error rate). Marks on the x-axis indicate the position of makers on the genetic map. Coloured boxes above markers indicate the 1.5-LOD confidence interval of each QTL. B Effect plots of each QTL. Cell death size across 160 RILs grouped by the parental allele (L58, red; R-o-18, blue). Black diamonds represent mean cell death size of all RILs within each allelic group

Table 2 Quantitative trait loci associated with HR-like cell death size in the L58 x R-o-18 RIL population

Validation of QTL effects on selected RIL lines

The validation of QTL effects was carried out on twelve selected RIL lines which showed contrasting genotypes at the peak markers of the three QTLs Pbc1–3 (Additional file 1: Fig. S7). Overall, we observed differences in egg-induced cell death between RILs (ANOVA: F11,24 = 5.06, P < 0.001), mostly due to allelic differences at Pbc1 (BrID11121), as the RILs with the L58 allele showed larger cell death size. Analysis of QTL effects with a three-way ANOVA showed a significant main effect for Pbc1 (F1, 32 = 84.02, P < 0.001) and Pbc3 (F1, 32 = 10.91, P = 0.002). The effect of Pbc3 was only significant upon inclusion of Pbc1 in the model (Fig. 4b), while no effect was detected for Pbc2 (F1, 32 = 3.05, P = 0.09) (Fig. 4a). Analysis of the Pbc1-Pbc3 haplotypes highlighted the large effect of the Pbc1-L58 allele and a marginal effect of both Pbc2-R-o-18 and Pbc3-R-o-18 alleles.

Fig. 4
figure 4

Validation of QTL effects and additive interactions for QTLs Pbc1-3 on selected RILs. Twelve RILs (N = 3) with contrasting genotypes at the peak markers of QTLs Pbc1–3 were selected randomly for a second phenotypic evaluation with P. brassicae egg clutches. RILs are grouped by genotype at the peak markers to show pairwise effects and additive interactions between QTLs. A) interaction between Pbc1 (BrID11121) and Pbc2 (BrID90099). B) interaction between Pbc1 and Pbc3 (BrID90095). Blue box with “L” = L58 allele, red box with “R” = R-o-18 allele. Boxplots represents the interquartile range (1st and 3rd quantile) and the median. White diamonds represent mean cell death of each QTL genotype. Letters report differences in mean size of HR-like cell death between haplotypes (Tukey’s HSD test, P < 0.01)

Identification of candidate genes underlying the QTLs

We investigated the genomic locations of the three QTLs for potential candidate genes associated with HR-like cell death using the B. rapa reference genome cv. Chiifu (v3.0). We searched for annotated genes that encode for cell surface receptors (PRRs), intracellular receptors (NLRs), or that are involved in general plant defense mechanisms, such as ROS production and cell death (Additional file 2: Supplementary Tables S8–10). The QTL Pbc1 showed the largest effect with the allele of the L58 parent contributing to a large cell death and it was located at the interval 8.65–25.29 Mb (± 1.5 LOD) on chromosome A02. This region contains 2012 annotated genes, of which 69 are related to the plant immunity and defense (Additional file 2: Supplementary Table S8). Among them, we found 14 cell surface receptors (of both the RLK and RLP type) and 19 intracellular TIR-NBS-LRR (TNL) receptors. Sixteen of the TNLs are closely located in three clusters, at the intervals 12.47–12.55 Mb, 21.64–21.73 Mb and 22.68–22.99 Mb. Moreover, Pbc1 also includes B. rapa homologs to three RLKs previously found to be upregulated upon oviposition in A. thaliana [24]: i.e., BraA02g017190.3C homolog of a LRKL domain-kinase protein (AT1G66880), BraA02g022170.3C and BraA02g022180.3C homologs of NEMATODE-INDUCED LRR-RLK 1 (NILR1, AT1G74360) and BraA02g033550.3C and BraA02g033570.3C, homologs of PBS1-LIKE 19 (AT5G47070). Further, Pbc1 region includes genes involved in mediating cell death processes, such as BraA02g032910.3C and BraA02g032940.3C that are both homologs of A. thaliana ACCELERATED CELL DEATH 11 (ACD11), and BraA02g033670.3C, a homolog of BAX-INHIBITOR 1 (BI-1, AT5G47120).

The minor QTL Pbc2 was located in the interval 25.34–38.15 Mb on A03 (± 1.5 LOD) and included 1594 genes in total, of which 49 being plant immunity-related genes (Additional file 2: Supplementary Table S9). Within Pbc2 we found different types of PRRs such as a cluster of 15 cysteine-rich RLKs (CRKs) at the interval 25.91–26.32 Mb and BraA03g059300.3C, homolog of L-TYPE LECTIN RECEPTOR KINASE I.9 (LecRK-I.9). In this region we found also NLRs, specifically a cluster of four TIR-NBS-LRR at the interval 25.52–25.56 Mb. Further, this region included also BraA03g053480.3C and BraA03g057870.3C, homologs of two known regulators of plant immunity, i.e. SUPPRESSOR-OF-NPR1 CONSTITUTIVE 4 (SNC4) and BRI1-ASSOCIATED RECEPTOR KINASE (BAK1), respectively, and two homologs of putative RESPIRATORY BURST OXIDASE HOMOLOGUE G (RbohG) genes. The third QTL, Pbc3, was located between 6.77 and 16.13 Mb on A06 (± 1.5 LOD). This region included a total of 2292 genes, of which 28 plant defense-related genes (Additional file 2: Supplementary Table S10). Within Pbc3 we found homologs of RbohD and RbohJ, different types of RLKs, i.e. homologs to two WALL-ASSOCIATED RECEPTOR KINASES 1 (WAK1, AT1G21250) and 2 (WAK2, AT1G21270), and, interestingly, a cluster of four L-type LecRKs including LecRK-I.1, that was recently associated to P. brassicae egg extract-induced cell death in A. thaliana [46].

Given that Pbc3 appeared to overlap with one of the two loci identified in A. thaliana by Groux et al. [46], we investigated the syntenic relationship between Pbc1–3 regions and A. thaliana genome. Pbc1 was syntenic to regions on A. thaliana chromosomes 1, on the top of chromosome 4 and on middle of chromosome 5 (Additional file 1: Supplementary Fig. S8). Pbc2, which was located to the bottom of chromosome A03, showed synteny to the bottom half of A. thaliana chromosome 4. Pbc3, which is located on the center of chromosome A06, was syntenic to regions on both A. thaliana chromosomes 3 and 5. Indeed, Pbc3 was syntenic to the region of A. thaliana chromosome 3 that included LecRK-I.1, candidate gene associated with egg extract-induced cell death [46]. Overall, Pbc1 and Pbc2 represented novel loci mediating HR-like cell death as they did not show synteny to the two loci previously identified in A. thaliana.


Here, we report the first genetic analysis of a butterfly egg-induced defence trait in an economically important crop B. rapa. We showed intraspecific variation for P. brassicae egg-induced HR-like cell death within a B. rapa germplasm collection. By developing an automated image-based phenotyping protocol, we could accurately measure HR-like cell death size in a RIL population and identified three new QTLs associated with this trait. The three QTLs include many candidate genes that are involved in plant immunity processes such as extra- and intra-cellular receptors, ROS production, and cell death.

Genetic mapping of egg-induced cell death identified three QTLs, Pbc1 (on chromosome A02), Pbc2 (A03) and Pbc3 (A06) that together explained about 31.5% of the phenotypic variance. Thus, in B. rapa, HR-like cell death size appears to be a polygenic trait similar to whitebacked planthopper egg-induced lesions in rice [47] or other pathogen-induced leaf necrotic symptoms [65, 66]. None of the three QTLs reported here have been validated yet, for example by using alternative segregating populations. Nevertheless, Pbc1 may represent a stable QTL as it explained the larger proportion of variance (17.9%), its confidence interval was entirely above the LOD threshold and it contributed to larger cell death with the allele of L58, the parent showing a stronger HR-like cell death. On the contrary, Pbc2 and Pbc3 represented minor QTLs, their LOD peaks were just above the LOD threshold, and their positive effect was due to alleles of R-o-18, the parent showing a smaller HR-like cell death. The unexplained phenotypic variance may be due to other undetected minor QTLs for which we expect that L58 alleles contribute to a larger cell death. In fact, only few RILs showed transgressive segregation beyond the mean value of L58 while the phenotypic distribution of the whole population was skewed towards the value of R-o-18, the parent showing a small cell death. Future research should validate the stability of the QTLs identified in this study, their (epistatic) interactions, and the effect of the plant genetic background by using other genetic populations/association panels and/or testing different environments.

In this study, we implemented the first image-based phenotyping method to assess insect egg-induced cell death on plant tissues and to perform QTL mapping. So far, image-based methods were used for genetics studies of plant disease symptoms [65,66,67,68] or insect feeding damage [54, 69, 70], but never for egg-induced responses. Nevertheless, our experiments showed an intermediate broad-sense heritability (H2 = 0.47–0.49) which points at a considerable environmental effect on the trait that our current bioassays and/or phenotyping protocol could not disentangle. Insect genetic variation is known to contribute to low heritability of genetic studies of plant defense traits to insects [71]. Unfortunately, at the moment it is difficult to control for the genetic variation of the P. brassicae butterflies used in our bioassay, beyond using a large number of adults and refreshing them regularly during experiments. Our imaged-based phenotyping protocol allowed us to measure cell death size with increased precision (highly repeatable measurements) but it can possibly be improved in accuracy (measuring the true size of cell death spots). Alternatively, the measurement of cell death size could be combined with other cell death-related traits, e.g. severity (variation in lesion colour), to possibly determine traits with higher heritability and increase the power of QTL detection.

The three QTLs Pbc1–3 provide a new source of candidate genes that will help to understand the molecular mechanisms underlying the interaction between P. brassicae eggs and Brassica host plants. Within these QTL regions, we identified different type of genes that are commonly involved in plant immunity processes, such as signalling/stress perception, ROS accumulation and cell death formation. Perception of pathogen-associated molecular patterns (PAMPs) by PRRs surface receptors triggers plant immunity [38]. Many PRRs belonging to different classes were found within the QTLs Pbc1–3, although involvement in plant immunity has been experimentally validated only for a few. For example, the Pbc1 region includes LecRK-V.5 which was previously found to negatively modulate plant immunity against biotrophic pathogens [72,73,74]. Further, Pbc1 includes also orthologs of three (predicted) RLKs (AT1G66880, AT1G74360, AT5G47070), which were previously found to be upregulated upon P. brassicae oviposition in A. thaliana [24]. One of these three, AT1G74360, encodes the LRR-RLK NILR1 that was found to be required for plant immunity to nematodes in a BAK1-dependent manner [75]. BAK1 is a known central regulator of pathogen-triggered immunity which works as a co-receptor of PRRs mediating the perception of many different biotic stresses [76], including feeding by insects, such as aphids and caterpillars [77, 78]. Interestingly, one B. rapa ortholog of BAK1, BraA03g057870.3C, is present within our QTL Pbc2. Whether BAK1 is also involved in the regulation of defences against Pieris eggs is an intriguing question that awaits future experimental validation. Pbc2 includes also many cysteine-rich RLKs (CRKs), of which CRK5 was experimentally shown to mediate pathogen-induced cell death [79] and CRK11 was upregulated upon P. brassicae egg deposition [24]. Finally, QTL Pbc3 includes a cluster of B. rapa orthologs of LecRK-I genes, among which LecRK-I.1, which was recently identified as one of the two A. thaliana loci associated with Pieris egg extract-induced cell death [46]. Considering the large confidence intervals of the QTL reported here, we cannot point yet at specific PRRs and/or RLKs as the casual genes responsible for the variation in HR-like cell death observed within our plant material.

Despite the many different RLKs present within the QTL Pbc1–3, it is still possible that variation in cell death size is dependent on other genetic mechanisms. Within our QTLs, we identified clusters of intracellular TIR-NBS-LRRs (TNLs) and a B. rapa homolog N REQUIREMENT GENE 1 (NRG1). NRG1 was shown to interact with ENHANCED DISEASE SUSCEPTIBILITY 1 (EDS1) and SENESCE-ASSOCIATED GENE 101 (SAG101) to form a protein complex that is required for TNL cell death signalling [80, 81]. It is interesting to recall that EDS1 and PHYTOALEXIN DEFICIENT 4 (PAD4), which forms with EDS1 another protein complex that is also required for TNL-mediated cell death [81, 82], were upregulated in A. thaliana upon P. brassicae egg deposition [24]. Despite all this, it is still largely speculative to state that insect eggs may be perceived by intracellular TNLs in absence of evidences on whether and how egg-associated molecular patterns can get into contact with the inside of plant cells.

Some plant species have been shown to accumulate ROS underneath insect eggs, such as A. thaliana, Pinus sylvestris or Brassica spp. [24, 35, 83], while others use ROS accumulation to directly kill eggs such as Solanum dulcamara [4]. Hence, it is suggested that ROS signalling in response to insect eggs is conserved in different plant species [36]. We identified different NADPH oxidases (Rboh proteins) which are known to be involved in production of ROS [84]. Further, we found genes involved in cell death regulation such as the B. rapa orthologs to BI-1 and ACD11 within QTL Pbc1. While BI-1 is a known suppressor of H2O2-dependent cell death [85] and has been associated with the cell death regulation in plant-powdery mildew interaction [86], ACD11 is involved in autoimmunity [87] and activation of cell death and defence responses [88]. Finally, we also found a homolog of ICS1/SID2 within Pbc1. This enzyme produces the precursor of SA [89], the main phytohormone so far associated with plant defences to insect eggs [24,25,26, 31]. Overall, the QTLs here presented include relatively large regions containing thousands of genes. Thus, future fine-mapping efforts are necessary to increase the resolution on the QTLs here reported and determine the exact genetic mechanisms of P. brassicae egg-induced cell death in B. rapa.

The identification of QTLs involved in HR-like cell death size in B. rapa also allows comparisons with the model plant A. thaliana. QTLs Pbc1 and Pbc2 showed no synteny to the two loci recently identified in A. thaliana [46], thus representing novel genomic regions associated with P. brassicae egg-induced cell death. Interestingly, the QTL Pbc3 is syntenic to the A. thaliana locus on chromosome 3 containing a cluster of clade I L-type LecRKs. The comparison in genetic architectures of HR-like cell death between the B. rapa and A. thaliana is intriguing but too premature to point at communalities and specificities. In fact, there are crucial differences in the plants’ responses, experimental setup and plant genetic material employed by the different studies. HR-like cell death in A. thaliana was visually scored in discrete categories which accounted also for light/strong chlorosis and light/strong cell death [46] On the contrary, in our study no chlorosis was observed and only cell death size was measured. Different phenotyping methods/parameters can affect genetic mapping results and indeed it was previously shown that shape and size of leaf disease symptoms can be genetically unliked in different pathosystems [66, 90]. Moreover, bioassays on A. thaliana were conducted by treating plants with P. brassicae egg extract, which was shown to also induce cell death in B. nigra [91]. Nonetheless, egg extract likely contains lipidic compounds from the inner egg tissues [44] and it is still unclear whether and how they are able to diffuse through the egg and reach the leaf surface. Clearly, more research is needed to further disentangle the genetic architectures of egg-induced cell death in A. thaliana and B. rapa.

It is remarkable that B. rapa showed a phenotypic variation in HR-like cell death that was limited to necrotic black spots varying in size. A similar mild cell death appearing as black spots was also observed underneath P. brassicae eggs on a limited number of B. oleracea accessions [11, 23]. This mild cell death contrasts sharply with the strong and severe cell death that we regularly observe on leaves of wild species of the tribe Brassiceae (Brassicaceae Lineage II), such as in Brassica nigra, Sinapis spp., and Crambe spp., which leads to reduced egg survival [11, 22]. The differences in HR-like cell death between wild Brassicaceae, e.g. Brassica nigra, and other Brassica crops raise questions on the role of crop domestication on this defence trait. Brassica crops are characterized by an extraordinary intraspecific diversity in morphotypes which differ significantly from the progenitor wild relatives as result of domestication [92,93,94]. The selection for specific crop morphotypes, but also for quality traits, such as flavour, taste and storage, mostly targeted leaf morphological and/or biochemical traits, which often show trade-offs with overall plant defense traits [95, 96]. Whether similar trade-offs also impacted the HR-like cell death expressed by current Brassica crop types should be tested. Certainly, we cannot yet conclude to have captured the full extent of intraspecific phenotypic variation as we screened only 56 accessions from one representative B. rapa core collection [61, 62]. Our choice aimed at encompassing the overall genetic diversity of the core collection while representing all species morphotypes. Nevertheless, we may have missed accessions with a stronger HR-like phenotype and/or alternative variants at the identified QTLs. In summary, future germplasm screenings should not only include more accessions, but also include Brassica wild material. What is the genetic basis for the differences in HR-like cell death severity between wild Brassicaceae, e.g. Brassica nigra, and other Brassica crops is a fascinating question that deserves to be addressed in future research.

The mild HR-like cell death observed in our B. rapa germplasm was shown to not affect egg survival [11]. This is expected to have certain implications regarding its deployment as egg-killing defense trait in plant breeding. In future, screening of more germplasms within Brassica crops, including crop wild relatives, could still identify strong HR-like cell death. An alternative possibility for deploying it as crop defense trait could be via introgression from B. nigra. Interspecific introgression of other disease resistance traits is a viable option as it is already being pursued within the Brassica genus by using interspecific crosses, embryo rescue, and marker-assisted selection [97, 98].


We report the identification of the first QTLs associated with a HR-like cell death induced by Pieris butterfly eggs in the economically important crop B. rapa. Our study confirms that plant genetic factors are involved in the elicitation of a HR-like cell death, a plant defense response against insect eggs. This work provides the basis for further identification of genes mediating the interaction between butterfly eggs and plants. Future studies should validate the QTLs by screening other genetic populations and/or association panels. Fine-mapping of the identified QTLs would then help to increase the resolution of the loci and further elucidate the genetic regulation of the egg-induced HR-like cell death.

Material and methods

Plant material

For germplasm screening, 56 Brassica rapa L. (Brassicaceae) accessions were selected representing all major B. rapa crop types (e.g., Chinese cabbage, Pak choi, turnip, oil types), and to include different levels of genetic heterogeneity, such as feral populations, landraces, breeding material, and doubled haploid (DH) lines (Additional file 1: Supplementary Fig. S1, Table S1). Most B. rapa accessions and all DH lines were obtained from the core collection of Dr. Bonnema at Plant Breeding, Wageningen University and Research [61, 62] with a few additional accessions obtained from the Centre for Genetic Resources (CGN, The Netherlands). A B. nigra accession previously reported to induce a strong HR-like cell death was used as positive control [22]. After the germplasm screening, ten accessions, considered suitable as potential parents of biparental mapping population, were selected for a second HR-like cell death evaluation. Criteria for the selection were: the accession i) displays an HR-like cell death score at the extremes of the phenotypic distribution and is consistent across individual plants; ii) is fast flowering (< 1 year); iii) is self-compatible; iv) was multiplied by selfings or was used to generate a DH line (in order to have homozygous material to repeat experiments); and, v) it was preferably showing comparable leaf phenotypes to minimize the segregation of leaf morphological traits after crossing. Finally, a mapping population of 160 recombinant inbred lines (RILs), previously generated from a cross between the B. rapa DH lines L58 (caixin type, ssp. parachinensis) and R-o-18 (yellow sarson type, spp. tricolaris) was used for QTL mapping of HR-like cell death size [63].

Plant growing conditions

Plants were grown in a greenhouse compartment under standardized conditions (21/18 °C day/night minimum temperature, 16/8 h light/dark photoperiod; and 50–70% relative humidity). The daily maximum temperature was not controlled and subjected to some variation (max + 5 °C). Seeds were vernalized at 4 °C for 2 days and then sown in small trays with sowing soil (Lensli, Bleiswijk, The Netherlands). Seedlings were transplanted 1 week after germination to 17 cm diameter pots with potting soil (Lensli, Bleiswijk, The Netherlands). Plants were grown for 5 weeks before being subjected to P. brassicae oviposition or treatment with egg wash.

Insect rearing

Pieris brassicae L. butterflies were obtained from a rearing facility of the Laboratory of Entomology, Wageningen University. Insects were kept in a greenhouse compartment (21 °C, 16/8 h dark photoperiod, 60–80% relative humidity). Larvae were reared on Brussel sprouts (Brassica oleracea var. gemmifera cv. Cyrus), while the adults were fed with a 10% honey solution and allowed to oviposit on the same plant.

Oviposition and egg wash treatment

Freshly eclosed P. brassicae female butterflies were allowed to mate, subsequently separated from the males, kept without plants for 2 days, and then used for no-choice oviposition experiments. Butterflies were allowed to freely oviposit on one B. rapa plant at the time in small cages each containing 2–4 female butterflies. A maximum of two egg clutches consisting of 10–20 eggs were laid on the two youngest fully developed leaves of each plant, while the other leaves were covered by a net. Then a new plant was placed into the cage. After every 4–5 egg-laden plants, mated female butterflies in each of the small cages were replaced to randomize the effect of insect genetic diversity.

For the egg wash treatment, P. brassicae egg clutches were collected from Brussel sprouts leaves within 24 h after oviposition. Eggs were carefully removed with a stainless-steel lab spatula and collected in an Eppendorf tube in a ratio of ~ 1000 eggs per 1 ml demineralized water. Eggs were incubated overnight at room temperature, after which the liquid phase was directly used or stored at − 20 °C.

Experimental design

Germplasm screening was carried out in September 2017 by application of egg wash. Two 5 μl droplets of egg wash were applied on the two youngest fully developed leaves of each plant. Droplets of an equivalent amount of demineralised water were applied as negative control. Each genotype/accession was represented by 3–5 replicates (individual plants). Plants were arranged in a randomized complete block design with five blocks and one plant per accession within each block.

Re-evaluation of ten homozygous lines, that were either inbred (CC-106, SO-040, BRO-127, IMB211, CC-168 BRO-030, CC-AO3, L58) or DH lines (R-o-18, R500, BRO-030, CC-AO3, L58), was carried out in February 2018 using both no-choice oviposition and egg wash treatment. Plants were arranged in a randomized complete block design with two blocks and five plants per accession within each block.

Three QTL mapping experiments were carried out in August/September 2018 using eggs deposited by P. brassicae females. The whole RIL population was grown three times over three consecutive weeks, each time with one replicate per RIL and three replicates for the two parents L58 and R-o-18. Validation of QTL effects and additive interactions was carried in September 2019 using twelve RILs (RIL_19, RIL_22, RIL_32, RIL_45, RIL_73, RIL_77, RIL_93, RIL_97, RIL_100, RIL_106, RIL_130, RIL_137) which were selected randomly for their contrasting genotypes at the peak markers of QTLs Pbc1–3. For all experiments with RILs, plants were subjected to no-choice oviposition as described above.

Assessment of HR-like cell death

Egg wash-induced HR-like cell death was scored in the germplasm screening on a scale from 0 to 3 with: 0, no visible symptoms; 1, a grey/dark spot smaller than droplet size; 2, a black necrotic spot covering the whole treated area; 3, strong cell death visible also on the adaxial side (Additional file 1: Supplementary Fig. S1a). For the re-evaluation of homozygous (inbred and DH) lines from ten selected accessions and for the QTL experiment, egg and egg wash-induced HR-like cell death size (area) was measured with an image-based phenotyping protocol (see next section). For all experiments, individual plants were assigned the highest HR-like cell death score or size out of all egg- and egg wash-treated spots. To account for variability in egg clutches size, cell death size measured underneath egg clutches was divided by the number of eggs in the clutch, and normalized to 10 eggs:

$${\mathrm{cell}\ \mathrm{death}}_{clutch}=\frac{\mathrm{cell}\ \mathrm{death}\ }{n\ \mathrm{eggs}} \ast 10$$

Development of image-based phenotyping protocol

To obtain a reliable and reproducible quantification of egg- and egg wash-induced HR-like cell death size (area), we developed a custom image-based phenotyping protocol (Additional file 1: Supplementary Fig. S2). Leaf discs containing egg clutches- or egg wash-treated spots were sampled with a cork borer of 6 mm diameter and placed in Petri dishes with 1% phytoagar (Duchefa Biochemie, Haarlem, The Netherlands) or wet filter paper (Additional file 1: Supplementary Fig. S2a). Just prior to the sampling of the leaf disks with deposited egg clutches, eggs were counted and then gently removed with adhesive tape to prevent leaf damaging. Spots treated with egg wash were sampled directly. Leaf discs were imaged with a Dino-Lite Edge Digital microscope (AnMo Electronics Corporation, Hsinchu, Taiwan) connected to a laptop (Additional file 1: Supplementary Fig. S2b). Each leaf disc was imaged with the light polarizer filter “fully open” using the following settings: LED zone 2 and 4: ON; LED zone 1 and 3: OFF; autoexposure: ON; white balance: STANDARD; output file format: PNG; resolution: 2592 × 1944 pixel.

Image analysis was performed on Fiji with ImageJ v1.52 software [99] using the image segmentation plugin Trainable WEKA Segmentation v3.2.28 [100]. Image analysis was executed through a custom Fiji macro script. In WEKA, the image segmentation was performed using the training features Minimum, Maximum, Mean, Variance, Median and the classifier algorithm FastRandomForest. For each leaf disc, first a classifier model was trained using a training set composed by representative image pixels that were labelled either as “healthy leaf tissue” or “HR-like cell death”. The trained classifier was then applied to the whole image to generate an 8-bit segmentation of HR-like cell death spots (Additional file 1: Supplementary Fig. S2d). The 8-bit segmented HR-like area was finally measured in Fiji using the command Analyze particles with Area as measurement (Additional file 1: Supplementary Fig. S2e). The use of the WEKA automated segmentation resulted in more reproducible measurements of cell death size (Additional file 1: Supplementary Fig. S3).

Phenotypic data analysis

All data analyses were performed in R 3.5.3 [101]. Raw data were firstly checked for assumptions of normality (Shapiro-Wilk normality test) and homogeneity of variances (Fligner-Killeen test). Non-normal data were analyzed after data transformation (root square on cell death size) or with a non-parametric test (cell death score). Phenotypic data obtained from germplasm screening were not normally distributed and thus analyzed with the non-parametric Kruskall-Wallis test. Post-hoc analysis was conducted with the Dunn test using Benjamini–Hochberg correction as implemented by the dunnTest function from FSA package [102]. Phenotypic data from re-evaluation of ten B. rapa accessions were analyzed on square root-transformed data. HR-like cell death size was analyzed by using the following model:

$$\mathrm{y}=\upmu +\mathrm{Block}+\mathrm{Block}:\mathrm{Row}+\mathrm{Block}:\mathrm{Col}+\mathrm{G}+\mathrm{r_{error}}$$

in which μ is the overall trait mean, Block is the blocking factor of the experimental design, Row and Col represents the spatial location of plants within a Block, G represents the B. rapa accession/genotype. Parsimonious models were explored by stepwise removal of each factor and comparison of the full versus reduced model with a Likelihood Ratio Test. The most parsimonious model resulted with only G as factor. Finally, data were analyzed using one-way ANOVA, followed by Tukey’s HSD post hoc test with Benjamini–Hochberg correction with α < 0.01 as implemented in multcomp package [103]. Phenotypic data from the QTL experiments and the re-evaluation of RILs were also analyzed on square root transformed data with the model in Eq. (2). Genotypic means of RILs across the three QTL experiments were calculated as the Best Linear Unbiased Estimator (BLUEs), using the mixed model:

$$\mathrm{BLUE}=\upmu +\mathrm{G}+\mathrm{Exp}+\mathrm{r_{error}}$$

in which μ is the overall trait mean, G represents the RIL genotype, Exp the QTL experiment (1–3). BLUEs were calculated by fitting G as fixed effect and Exp as random effect. The model was analyzed by REML procedure using the function lmer from lme4 package [104].

Estimation of variance components

Variance components for genetic and experimental residual error were estimated with Eq. (2) fitted as mixed model with all factors included as random [105]. The model was analyzed by R EML procedure using the function lmer from lme4 package [104]. Classic broad-sense heritability (H2) was calculated by using the estimated variance components with the formula σG / (σG + σE) as previo usly described (Hollander 2003). Genetic and environmental coefficient of variation (CV) was calculated according to the equation:


in which n is the grand mean of the population, and σ2X is a variance component (σ2G or σ2E).

Linkage map construction and QTL analysis

A combined and denser genetic map for the RIL population L58 x R-o-18 was created using marker data previously generated in two separate studies, that is AFLP, SSR and SNP markers [63], and InDel PCR markers [64] (Additional file 2: Supplementary Tables S4-S5). We used Haldane’s mapping function with default setting as implemented in JoinMap 4.0 [106], to convert recombination frequencies to centiMorgan (cM). The final linkage map constituted of 485 markers and 10 linkage groups corresponding to the 10 chromosomes of the B. rapa A genome for a total of 1154.44 cM (Additional file 1: Supplementary Fig. S4, Additional file 2: Supplementary Tables S6-S7).

QTL analysis was performed using the R/qtl package in R 3.5.3 [107]. Genotype probabilities at positions not covered by the linkage map were estimated every 1 cM with the calc.genoprob function (step size = 1). First, single QTL models were searched with the scanone function using an interval mapping method (Haley-Knott regression). Subsequently, multiple-QTL model (MQM) interval mapping using Haley-Knott regression was performed to investigate multiple-QTL models which included the previously identified QTLs and additional (potential) QTLs with the mqmscan function. As separate analysis, MQM was also implemented with the stepwiseqtl function which gave similar results. Finally, epistatic additive and interactions were investigated with pairwise two-QTLs models as implemented in the scantwo command. For all analysis, LOD score significance threshold at 5% error rate was estimated with a 1000 permutations test.

Identification of candidate genes withing the QTLs regions

In order to investigate the gene content underlying the identified QTLs, the linkage map was anchored to the B. rapa reference genome Chiifu v3.0 (downloaded from BRAD Brassica database, accessed in August 2021) [108]. Sequences of the QTL-flanking markers were aligned to the reference genome using Geneious Prime v8 [109] to extract their physical location and all B. rapa genes contained within them. B. rapa gene functional annotation was assigned as the best match to A. thaliana protein database (genome TAIR 10) using BLAST+ v2.12.0 (E-value = 1e− 5). Candidate genes associated with plant defense, biotic stress and cell death, including cell surface and intracellular receptors, were manually searched within the description of the A. thaliana orthologs to each B. rapa gene.

Analysis of syntenic relationships between B. rapa and Arabidopsis thaliana genomes was performed using the comparative genomic tool SynMap on the CoGe web platform [110,111,112]. SynMap legacy version was used with the following settings: DAGChainer Options “Maximum distance between two matches (-D): 20 genes”; “Minimum number of aligned pairs (-A): 5 genes”; Merge Syntenic Blocks “Algorithm: Quota align merge”; Syntenic Depth “Algorithm: Quota Align”, “Ratio of coverage depth (A. thaliana) 1 -to- 3 (B. rapa)” , “Overlap distance 40”; Fractionation Bias “Run OFF”; CodeML “Calculate syntenic CDS pairs: Synonymous (Ks) substitution rate; “Color scheme: Rainbow 2”, “Max Value: 2”, “Log10 Transform: OFF”; Advance Options “Tandem duplication distance: 10”.

Availability of data and materials

All data generated or analysed during this study are included in this published article (Additional file 3). Datasets and scripts used for data analysis are also available in a Zenodo repository (



Quantitative trait locus


Hypersensitive response


Reactive oxygen species


Salicylic acid


Pattern recognition receptor


Receptor-like protein


Receptor-like kinase


L-type lectin receptor kinase


Cysteine-rich receptor kinase


Wall-associated receptor kinase


Nucleotide-binding leucin-rich repeat receptor


Toll-interleukin 1 receptor domain NLR


  1. Hilker M, Fatouros NE. Plant responses to insect egg deposition. Annu Rev Entomol. 2015;60:493–515.

    Article  CAS  PubMed  Google Scholar 

  2. Doss RP, Oliver JE, Proebsting WM, Potter SW, Kuy S, Clement SL, et al. Bruchins: insect-derived plant regulators that stimulate neoplasm formation. Proc Natl Acad Sci. 2000;97:6218–23.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Petzold-Maxwell J, Wong S, Arellano C, Gould F. Host plant direct defence against eggs of its specialist herbivore, Heliothis subflexa. Ecol Entomol. 2011;36:700–8.

    Article  Google Scholar 

  4. Geuss D, Stelzer S, Lortzing T, Steppuhn A. Solanum dulcamara’s response to eggs of an insect herbivore comprises ovicidal hydrogen peroxide production. Plant Cell Environ. 2017;40:2663–77.

    Article  CAS  PubMed  Google Scholar 

  5. Seino Y, Suzuki Y, Sogawa K. An Ovicidal substance produced by rice plants in response to oviposition by the whitebacked planthopper, Sogatella furcifera (Horváth) (Homoptera: Delphacidae). Appl Entomol Zool. 1996;31:467–73.

    Article  CAS  Google Scholar 

  6. Desurmont GA, Weston PA. Aggregative oviposition of a phytophagous beetle overcomes egg-crushing plant defences. Ecol Entomol. 2011;36:335–43.

    Article  Google Scholar 

  7. Shapiro AM, DeVay JE. Hypersensitivity reaction of Brassica nigra L. (Cruciferae) kills eggs of Pieris butterflies (Lepidoptera: Pieridae). Oecologia. 1987;71:631–2.

    Article  CAS  PubMed  Google Scholar 

  8. Balbyshev NF, Lorenzen JH. Hypersensitivity and egg drop: a novel mechanism of host plant resistance to Colorado potato beetle (Coleoptera: Chrysomelidae). J Econ Entomol. 1997;90:652–7.

    Article  Google Scholar 

  9. Garza R, Vera J, Cardona C, Barcenas N, Singh SP. Hypersensitive response of beans to Apion godmani (Coleoptera: Curculionidae). J Econ Entomol. 2001;94:958–62.

    Article  CAS  PubMed  Google Scholar 

  10. Fatouros NE, Cusumano A, Danchin EGJ, Colazza S. Prospects of herbivore egg-killing plant defenses for sustainable crop protection. Ecol Evol. 2016;6:6906–18.

    Article  PubMed  PubMed Central  Google Scholar 

  11. Griese E, Caarls L, Bassetti N, Mohammadin S, Verbaarschot P, Bukovinszkine’Kiss G et al. Insect egg-killing: a new front on the evolutionary arms-race between brassicaceous plants and pierid butterflies. New Phytol 2021; 230: 341–353.

  12. Kumar S, Sarwan KC. Assessment of avoidable yield losses in crop brassicas by insect-pests. J Entomol Zool Stud. 2017;5:1814–8.

    Google Scholar 

  13. Ryan SF, Lombaert E, Espeset A, Vila R, Talavera G, Dinca V, et al. Global invasion history of the agricultural pest butterfly Pieris rapae revealed with genomics and citizen science. Proc Natl Acad Sci. 2019;116:20015–24.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Kourelis J, Van Der Hoorn RAL. Defended to the nines: 25 years of resistance gene cloning identifies nine mechanisms for R protein function. Plant Cell. 2018;30:285–99.

    Article  CAS  PubMed  Google Scholar 

  15. Nallu S, Hill JA, Don K, Sahagun C, Zhang W, Meslin C, et al. The molecular genetic basis of herbivory between butterflies and their host plants. Nat Ecol Evol. 2018;2:1418–27.

    Article  PubMed  PubMed Central  Google Scholar 

  16. Erb M, Reymond P. Molecular interactions between plants and insect herbivores. Annu Rev Plant Biol. 2019;70:527–57.

    Article  CAS  PubMed  Google Scholar 

  17. De-la-Cruz IM, Velázquez-Márquez S, Núñez-Farfán J. What do we know about the genetic basis of plant defensive responses to herbivores? A Minireview. In: Núñez-Farfán J, Valverde PL, editors. Evolutionary ecology of plant-herbivore interaction. Cham: Springer; 2020:295–314.

  18. Wheat CW, Vogel H, Wittstock U, Braby MF, Underwood D, Mitchell-Olds T. The genetic basis of a plant-insect coevolutionary key innovation. Proc Natl Acad Sci. 2007;104:20427–31.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Kliebenstein DJ, Kroymann J, Mitchell-Olds T. The glucosinolate-myrosinase system in an ecological and evolutionary context. Curr Opin Plant Biol. 2005;8:264–71.

    Article  CAS  PubMed  Google Scholar 

  20. Fatouros NE, Lucas-Barbosa D, Weldegergis BT, Pashalidou FG, van Loon JJA, Dicke M, et al. Plant volatiles induced by herbivore egg deposition affect insects of different trophic levels. PLoS One. 2012;7:e43607.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Fatouros NE, Pineda A, Huigens ME, Broekgaarden C, Shimwela MM, Figueroa Candia IA, et al. Synergistic effects of direct and indirect defences on herbivore egg survival in a wild crucifer. Proc R Soc. 2014;281:20141254.

    Google Scholar 

  22. Griese E, Dicke M, Hilker M, Fatouros NE. Plant response to butterfly eggs: inducibility, severity and success of egg-killing leaf necrosis depends on plant genotype and egg clustering. Sci Rep. 2017;7:7316.

    Article  PubMed  PubMed Central  Google Scholar 

  23. Reymond P. Perception, signaling and molecular basis of oviposition-mediated plant responses. Planta. 2013;238:247–58.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Little D, Gouhier-Darimont C, Bruessow F, Reymond P. Oviposition by pierid butterflies triggers defense responses in Arabidopsis. Plant Physiol. 2007;143:784–800.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Bruessow F, Gouhier-Darimont C, Buchala A, Metraux JP, Reymond P. Insect eggs suppress plant defence against chewing herbivores. Plant J. 2010;62:876–85.

    Article  CAS  PubMed  Google Scholar 

  26. Gouhier-Darimont C, Schmiesing A, Bonnet C, Lassueur S, Reymond P. Signalling of Arabidopsis thaliana response to Pieris brassicae eggs shares similarities with PAMP-triggered immunity. J Exp Bot. 2013;64:665–74.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Fatouros NE, Broekgaarden C, Bukovinszkine’Kiss G, Van Loon JJA, Mumm R, Huigens ME, et al. Male-derived butterfly anti-aphrodisiac mediates induced indirect plant defense. Proc Natl Acad Sci. 2008;105:10033–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Firtzlaff V, Oberländer J, Geiselhardt S, Hilker M, Kunze R. Pre-exposure of Arabidopsis to the abiotic or biotic environmental stimuli ‘chilling’ or ‘insect eggs’ exhibits different transcriptomic responses to herbivory. Sci Rep. 2016;6:1–13.

    Article  Google Scholar 

  29. Baruah IK, Panda D, Jagadale MV, Das DJ, Acharjee S, Sen P, et al. Bruchid egg induced transcript dynamics in developing seeds of black gram (Vigna mungo). PLoS One. 2017;12:e0176337.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Drok S, Bandoly M, Stelzer S, Lortzing T, Steppuhn A. Moth oviposition shapes the species-specific transcriptional and phytohormonal response of Nicotiana attenuata to larval feeding. Sci Rep. 2018;8:1–14.

    Article  CAS  Google Scholar 

  31. Lortzing V, Oberländer J, Lortzing T, Tohge T, Steppuhn A, Kunze R, et al. Insect egg deposition renders plant defence against hatching larvae more effective in a salicylic acid-dependent manner. Plant Cell Environ. 2019;42:1019–32.

    Article  CAS  PubMed  Google Scholar 

  32. Valsamakis G, Bittner N, Fatouros NE, Kunze R, Hilker M, Lortzing V. Priming by timing: Arabidopsis thaliana adjusts its priming response to Lepidoptera eggs to the time of larval hatching. Front Plant Sci. 2020;11:1969.

    Article  Google Scholar 

  33. Das D, Baruah IK, Panda D, Paswan RR, Acharjee S, Sarmah BK. Bruchid beetle ovipositioning mediated defense responses in black gram pods. BMC Plant Biol. 2021;21:1–22.

    Article  Google Scholar 

  34. Fatouros NE, Voirol LRP, Drizou F, Doan QT, Pineda A, Frago E, et al. Role of large cabbage white butterfly male-derived compounds in elicitation of direct and indirect egg-killing defenses in the black mustard. Front Plant Sci. 2015;6:794.

    Article  PubMed  PubMed Central  Google Scholar 

  35. Caarls L, Bassetti N, van Doesburg F, Verbaarschot P, van Loon JJA, Schranz ME, et al. Deciphering Brassica plant defence responses to cabbage white butterfly egg-associated molecular patterns. bioRxiv. 2021; Preprint at

  36. Lortzing T, Kunze R, Steppuhn A, Hilker M, Lortzing V. Arabidopsis, tobacco, nightshade and elm take insect eggs as herbivore alarm and show similar transcriptomic alarm responses. Sci Rep. 2020;10:1–16.

    Article  Google Scholar 

  37. Erb M, Meldau S, Howe GA. Role of phytohormones in insect-specific plant reactions. Trends Plant Sci. 2012;17:250–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Couto D, Zipfel C. Regulation of pattern recognition receptor signalling in plants. Nat Rev Immunol. 2016;16:537–52.

    Article  CAS  PubMed  Google Scholar 

  39. Tang D, Wang G, Zhou JM. Receptor kinases in plant-pathogen interactions: more than pattern recognition. Plant Cell. 2017;29:618–37.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. van der Burgh AM, Joosten MHAJ. Plant immunity: thinking outside and inside the box. Trends Plant Sci. 2019;24:587–601.

    Article  PubMed  Google Scholar 

  41. Acevedo FE, Peiffer M, Tan CW, Stanley BA, Stanley A, Wang J, et al. Fall armyworm-associated gut bacteria modulate plant defense responses. Mol Plant-Microbe Interact. 2017;30:127–37.

    Article  CAS  PubMed  Google Scholar 

  42. Gust AA, Pruitt R, Nürnberger T. Sensing danger: key to activating plant immunity. Trends Plant Sci. 2017;22:779–91.

    Article  CAS  PubMed  Google Scholar 

  43. Hilker M, Stein C, Schröder R, Varama M, Mumm R. Insect egg deposition induces defence responses in Pinus sylvestris: characterisation of the elicitor. J Exp Biol. 2005;208:1849–54.

    Article  PubMed  Google Scholar 

  44. Stahl E, Brillatz T, Queiroz EF, Marcourt L, Schmiesing A, Hilfiker O, et al. Phosphatidylcholines from Pieris brassicae eggs activate an immune response in Arabidopsis. eLife. 2020;9:1–21.

    Article  Google Scholar 

  45. Gouhier-Darimont C, Stahl E, Glauser G, Reymond P. The Arabidopsis lectin receptor kinase Lecrk-I.8 is involved in insect egg perception. Front Plant Sci. 2019;10:623.

  46. Groux R, Stahl E, Gouhier-Darimont C, Kerdaffrec E, Jimenez-Sandoval P, Santiago J, et al. Arabidopsis natural variation in insect egg-induced cell death reveals a role for LECTIN RECEPTOR KINASE-I.1. Plant Physiol. 2021;185:240–55.

    CAS  PubMed  Google Scholar 

  47. Yang Y, Xu J, Leng Y, Xiong G, Hu J, Zhang G, et al. Quantitative trait loci identification, fine mapping and gene expression profiling for ovicidal response to whitebacked planthopper (Sogatella furcifera Horváth) in rice (Oryza sativa L.). BMC Plant Biol. 2014;14:145.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Mariyammal I, Seram D, Samyuktha SM, Karthikeyan A, Dhasarathan M, Murukarthick J, et al. QTL mapping in Vigna radiata × Vigna umbellata population uncovers major genomic regions associated with bruchid resistance. Mol Breed. 2019;39:1–13.

    Article  CAS  Google Scholar 

  49. Tamiru A, Paliwal R, Manthi SJ, Odeny DA, Midega CAO, Khan ZR, et al. Genome wide association analysis of a stemborer egg induced “call-for-help” defence trait in maize. Sci Rep. 2020;10:1–12.

    Article  Google Scholar 

  50. Pashalidou FG, Fatouros NE, van Loon JJA, Dicke M, Gols R. Plant-mediated effects of butterfly egg deposition on subsequent caterpillar and pupal development, across different species of wild Brassicaceae. Ecol Entomol. 2015;40:444–50.

    Article  Google Scholar 

  51. Harvey JA, Witjes LMA, Benkirane M, Duyts H, Wagenaar R. Nutritional suitability and ecological relevance of Arabidopsis thaliana and Brassica oleracea as foodplants for the cabbage butterfly, Pieris rapae. Plant Ecol. 2007;189:117–26.

    Article  Google Scholar 

  52. Dogimont C, Chovelon V, Pauquet J, Boualem A, Bendahmane A. The vat locus encodes for a CC-NBS-LRR protein that confers resistance to Aphis gossypii infestation and A. gossypii -mediated virus resistance. Plant J. 2014;80:993–1004.

    Article  CAS  PubMed  Google Scholar 

  53. Liu Y, Wu H, Chen H, Liu Y, He J, Kang H, et al. A gene cluster encoding lectin receptor kinases confers broad-spectrum and durable insect resistance in rice. Nat Biotechnol. 2015;33:301–7.

    Article  CAS  PubMed  Google Scholar 

  54. Thoen MPM, Davila Olivas NH, Kloth KJ, Coolen S, Huang PP, Aarts MGM, et al. Genetic architecture of plant stress resistance: multi-trait genome-wide association mapping. New Phytol. 2017;213:1346–62.

    Article  CAS  PubMed  Google Scholar 

  55. Gust AA, Nürnberger T. A plant surface receptor for sensing insect herbivory. Proc Natl Acad Sci. 2020;117:32839–41.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Sun M, Voorrips RE, van’t Westende W, van Kaauwen M, Visser RGF, Vosman B. Aphid resistance in Capsicum maps to a locus containing LRR-RLK gene analogues. Theor Appl Genet. 2020;133:227–37.

    Article  CAS  PubMed  Google Scholar 

  57. Zhang L, Cai X, Wu J, Liu M, Grob S, Cheng F, et al. Improved Brassica rapa reference genome by single-molecule sequencing and chromosome conformation capture technologies. Hortic Res 2018 51. 2018;5:1–11.

    Article  Google Scholar 

  58. Belser C, Istace B, Denis E, Dubarry M, Baurens FC, Falentin C, et al. Chromosome-scale assemblies of plant genomes using nanopore long reads and optical maps. Nat Plants. 2018;4:879–87.

    Article  CAS  PubMed  Google Scholar 

  59. Lou P, Woody S, Greenham K, VanBuren R, Colle M, Edger PP, et al. Genetic and genomic resources to study natural variation in Brassica rapa. Plant Direct. 2020;4:e00285.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  60. Li Y, Liu GF, Ma LM, Liu TK, Zhang CW, Xiao D, et al. A chromosome-level reference genome of non-heading Chinese cabbage [Brassica campestris (syn. Brassica rapa) ssp. chinensis]. Hortic Res. 2020;7:1–13.

    Article  PubMed  PubMed Central  Google Scholar 

  61. Zhao J, Artemyeva A, Del Carpio DP, Basnet RK, Zhang N, Gao J, et al. Design of a Brassica rapa core collection for association mapping studies. Genome. 2010;53:884–98.

    Article  CAS  PubMed  Google Scholar 

  62. Del Carpio DP, Basnet RK, de Vos RCH, Maliepaard C, Visser R, Bonnema G. The patterns of population differentiation in a Brassica rapa core collection. Theor Appl Genet. 2011;122:1105–18.

    Article  Google Scholar 

  63. Bagheri H, El-Soda M, van Oorschot I, Hanhart C, Bonnema G, van den Bosch TJ, et al. Genetic analysis of morphological traits in a new, versatile, rapid-cycling Brassica rapa recombinant inbred line population. Front Plant Sci. 2012;3:183.

  64. Zhang X, Meng L, Liu B, Hu Y, Cheng F, Liang J, et al. A transposon insertion in FLOWERING LOCUS T is associated with delayed flowering in Brassica rapa. Plant Sci. 2015;241:211–20.

    Article  CAS  PubMed  Google Scholar 

  65. Corwin JA, Copeland D, Feusier J, Subedy A, Eshbaugh R, Palmer C, et al. The quantitative basis of the Arabidopsis innate immune system to endemic pathogens depends on pathogen genetics. PLoS Genet. 2016;12:e1005789.

    Article  PubMed  PubMed Central  Google Scholar 

  66. Yates S, Mikaberidze A, Krattinger SG, Abrouk M, Hund A, Yu K, et al. Precision phenotyping reveals novel loci for quantitative resistance to septoria tritici blotch. Plant Phenomics. 2019;2019.

  67. Stewart EL, Croll D, Lendenmann MH, Sanchez-Vallet A, Hartmann FE, Palma-Guerrero J, et al. Quantitative trait locus mapping reveals complex genetic architecture of quantitative virulence in the wheat pathogen Zymoseptoria tritici. Mol. Plant Pathol. 2017;19:201–16.

    Google Scholar 

  68. Fordyce RF, Soltis NE, Caseys C, Gwinner R, Corwin JA, Atwell S, et al. Digital imaging combined with genome-wide association mapping links loci to plant-pathogen interaction traits. Plant Physiol. 2018;178:1406–22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  69. Kloth KJ, Busscher-Lange J, Wiegers GL, Kruijer W, Buijs G, Meyer RC, et al. SIEVE ELEMENT-LINING CHAPERONE1 restricts aphid feeding on Arabidopsis during heat stress. Plant Cell. 2017;29:2450–64.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  70. Visschers IGS, Peters JL, Timmermans LLH, Edwards E, Ferrater JB, Balatero CH, et al. Resistance to three thrips species in Capsicum spp. depends on site conditions and geographic regions. J Appl Entomol. 2019;143:929–41.

    Article  Google Scholar 

  71. Kliebenstein DJ. Quantitative genetics and genomics of plant resistance to insects. Annu Plant Rev. 2017;47:235–62.

    Article  Google Scholar 

  72. Arnaud D, Desclos-Theveniau M, Zimmerli L. Disease resistance to Pectobacterium carotovorum is negatively modulated by the Arabidopsis lectin receptor kinase LecRK-V.5. Plant Signal Behav. 2012;7:1070–2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Desclos-Theveniau M, Arnaud D, Huang TY, Lin GJC, Chen WY, Lin YC, et al. The Arabidopsis lectin receptor kinase LecRK-V.5 represses stomatal immunity induced by pseudomonas syringae pv. Tomato DC3000. PLoS Pathog. 2012;8:e1002513.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Wang Y, Bouwmeester K, Beseh P, Shan W, Govers F. Phenotypic analyses of Arabidopsis T-DNA insertion lines and expression profiling reveal that multiple L-type lectin receptor kinases are involved in plant immunity. Mol Plant-Microbe Interact. 2014;27:1390–402.

    Article  PubMed  Google Scholar 

  75. Mendy B, Wang’ombe MW, Radakovic ZS, Holbein J, Ilyas M, Chopra D, et al. Arabidopsis leucine-rich repeat receptor–like kinase NILR1 is required for induction of innate immunity to parasitic nematodes. PLoS Pathog. 2017;13:e1006284.

    Article  PubMed  PubMed Central  Google Scholar 

  76. Roux M, Schwessinger B, Albrecht C, Chinchilla D, Jones A, Holton N, et al. The Arabidopsis leucine-rich repeat receptor-like kinases BAK1/SERK3 and BKK1/SERK4 are required for innate immunity to hemibiotrophic and biotrophic pathogens. Plant Cell. 2011;23:2440–55.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  77. Yang D-H, Hettenhausen C, Baldwin IT, Wu J. BAK1 regulates the accumulation of jasmonic acid and the levels of trypsin proteinase inhibitors in Nicotiana attenuata’s responses to herbivory. J Exp Bot. 2011;62:641–52.

    Article  CAS  PubMed  Google Scholar 

  78. Vincent TR, Avramova M, Canham J, Higgins P, Bilkey N, Mugford ST, et al. Interplay of plasma membrane and vacuolar ion channels, together with BAK1, elicits rapid cytosolic calcium elevations in Arabidopsis during aphid feeding. Plant Cell. 2017;29:1460–79.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. Chen K, Fan B, Du L, Chen Z. Activation of hypersensitive cell death by pathogen-induced receptor-like protein kinases from Arabidopsis. Plant Mol Biol. 2004;56:271–83.

    Article  CAS  PubMed  Google Scholar 

  80. Lapin D, Kovacova V, Sun X, Dongus JA, Bhandari D, Von Born P, et al. A coevolved EDS1-SAG101-NRG1 module mediates cell death signaling by TIR-domain immune receptors. Plant Cell. 2019;31:2430–55.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  81. Sun X, Lapin D, Feehan JM, Stolze SC, Kramer K, Dongus JA, et al. Pathogen effector recognition-dependent association of NRG1 with EDS1 and SAG101 in TNL receptor immunity. Nat Commun. 2021;12:1–15.

    Google Scholar 

  82. Wagner S, Stuttmann J, Rietz S, Guerois R, Brunstein E, Bautor J, et al. Structural basis for signaling by exclusive EDS1 heteromeric complexes with SAG101 or PAD4 in plant innate immunity. Cell Host Microbe. 2013;14:619–30.

    Article  CAS  PubMed  Google Scholar 

  83. Bittner N, Trauer-Kizilelma U, Hilker M. Early plant defence against insect attack: involvement of reactive oxygen species in plant responses to insect egg deposition. Planta. 2017;245:993–1007.

    Article  CAS  PubMed  Google Scholar 

  84. Torres MA. ROS in biotic interactions. Physiol Plant. 2010;138:414–29.

    Article  CAS  PubMed  Google Scholar 

  85. Ishikawa T, Uchimiya H, Kawai-Yamada M. The role of plant Bax inhibitor-1 in suppressing H2O2-induced cell death. In: Cadenas E, Packer L, editors. Methods in enzymology. Cambridge: Academic Press; 2013:239–256.

  86. Weis C, Pfeilmeier S, Glawischnig E, Isono E, Pachl F, Hahne H, et al. Co-immunoprecipitation-based identification of putative BAX INHIBITOR-1-interacting proteins involved in cell death regulation and plant–powdery mildew interactions. Mol Plant Pathol. 2013;14:791–802.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  87. Brodersen P, Petersen M, Pike HM, Olszak B, Skov S, Ødum N, et al. Knockout of Arabidopsis accelerated-cell-death11 encoding a sphingosine transfer protein causes activation of programmed cell death and defense. Genes Dev. 2002;16:490–502.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  88. Li Q, Ai G, Shen D, Zou F, Wang J, Bai T, et al. A Phytophthora capsici effector targets ACD11 binding partners that regulate ROS-mediated defense response in Arabidopsis. Mol Plant. 2019;12:565–81.

    Article  CAS  PubMed  Google Scholar 

  89. Lefevere H, Bauters L, Gheysen G. Salicylic acid biosynthesis in plants. Front Plant Sci. 2020;11:338.

    Article  PubMed  PubMed Central  Google Scholar 

  90. Méline V, Caldwell DL, Kim B-S, Baireddy S, Yang C, Sparks EE, et al. Image-based assessment of plant disease progression identifies new genetic loci for resistance. bioRxiv. 2021; Preprint at

  91. Groux R, Fouillen L, Mongrand S, Reymond P. Sphingolipids are involved in Pieris brassicae egg-induced cell death in Arabidopsis thaliana. bioRxiv. 2021; Preprint at

  92. Cheng F, Sun R, Hou X, Zheng H, Zhang F, Zhang Y, et al. Subgenome parallel selection is associated with morphotype diversification and convergent crop domestication in Brassica rapa and Brassica oleracea. Nat Genet. 2016;48:1218–24.

    Article  CAS  PubMed  Google Scholar 

  93. McAlvay AC, Ragsdale AP, Mabry ME, Qi X, Bird KA, Velasco P, et al. Brassica rapa domestication: untangling wild and feral forms and convergence of crop morphotypes. Mol Biol Evol. 2021;38:3358–72.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  94. Mabry ME, Turner-Hissong SD, Gallagher EY, McAlvay AC, An H, Edger PP, et al. The evolutionary history of wild, domesticated, and feral Brassica oleracea (Brassicaceae). Mol Biol Evol. 2021;38:4419–34.

    Article  PubMed  PubMed Central  Google Scholar 

  95. Turcotte MM, Turley NE, Johnson MTJ. The impact of domestication on resistance to two generalist herbivores across 29 independent domestication events. New Phytol. 2014;204:671–81.

    Article  PubMed  Google Scholar 

  96. Whitehead SR, Turcotte MM, Poveda K. Domestication impacts on plant-herbivore interactions: a meta-analysis. Philos Trans R Soc B Biol Sci. 2017;372:20160034.

    Article  Google Scholar 

  97. Diederichsen E, Frauen M, Linders EGA, Hatakeyama K, Hirai M. Status and perspectives of clubroot resistance breeding in crucifer crops. J Plant Growth Regul. 2009;28:265–81.

    Article  CAS  Google Scholar 

  98. Lv H, Fang Z, Yang L, Zhang Y, Wang Y. An update on the arsenal: mining resistance genes for disease management of Brassica crops in the genomic era. Hortic Res. 2020;7:34.

    Article  PubMed  PubMed Central  Google Scholar 

  99. Schindelin J, Arganda-Carreras I, Frise E, Kaynig V, Longair M, Pietzsch T, et al. Fiji: an open-source platform for biological-image analysis. Nat Methods. 2012;9:676–82.

    Article  CAS  PubMed  Google Scholar 

  100. Arganda-Carreras I, Kaynig V, Rueden C, Eliceiri KW, Schindelin J, Cardona A, et al. Trainable Weka segmentation: a machine learning tool for microscopy pixel classification. Bioinformatics. 2017;33:2424–6.

    Article  CAS  PubMed  Google Scholar 

  101. R Core Team. R: a language and environment for statistical computing. 2021.

    Google Scholar 

  102. Ogle D, Doll J, Wheeler P, Dinno A. FSA: simple fisheries stock assessment methods. 2021.

    Google Scholar 

  103. Hothorn T, Bretz F, Westfall P. multcomp: simultaneous inference in general parametric models. 2021.

    Google Scholar 

  104. Bates D, Maechler M, Bolker B, Walker S. lme4: linear mixed-effects models using eigen and S4. 2021.

    Google Scholar 

  105. Khan MS, Yin X, van der Putten PEL, Jansen HJ, van Eck HJ, van Eeuwijk FA et al. A model-based approach to analyse genetic variation in potato using standard cultivars and a segregating population. II. Tuber bulking and resource use efficiency. F Crop Res. 2019;242:107582.

  106. Kyazma. JoinMap.

  107. Broman KW, Wu H, Sen S, Churchill G. R/qtl: QTL mapping in experimental crosses. Bioinformatics, 19: 889–890.

  108. BRAD Brassica Database. Accessed 30 Aug 2021.

  109. Biomatters. Geneious Prime.

  110. Lyons E, Pedersen B, Kane J, Alam M, Ming R, Tang H, et al. Finding and comparing syntenic regions among Arabidopsis and the outgroups papaya, poplar, and grape: CoGe with rosids. Plant Physiol. 2008;148:1772–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  111. Haug-Baltzell A, Stephens SA, Davey S, Scheidegger CE, Lyons E. SynMap2 and SynMap3D: web-based whole-genome synteny browsers. Bioinformatics. 2017;33:2197–8.

    Article  CAS  PubMed  Google Scholar 

  112. CoGe: Comparative Genomics. Accessed 30 Aug 2021.

Download references


We are thankful to: Johan Bucher (WUR Laboratory of Plant Breeding) for providing us seeds of the B. rapa core collection and the Centre for Genetic Resources the Netherlands (CGN) for providing additional accessions; the employees of Unifarm (WUR) for taking care of the plants used in our experiments; Pieter Rouweler and André Gidding for rearing of Pieris brassicae; Eline Ampt for the introduction to ImageJ. We are grateful to Dr. Ben Vosman for his input on the experimental designs of both germplasm and RILs screening. Oviposition on the RIL population was possible thanks the help of Dr. Thijmen Breeschoten, Frank Becker, Marco Busscher, Nynke Groendijk-Wilders, Dr. Dêêdi Sogbohossou and Xiuji Wang. We finally thank Dr. Marcos Malosetti and Roel van Bezouw for help with statistical analysis and QTL mapping.

Research involving plants

All methods were performed in accordance with the relevant regulations and legislation. All B. rapa plant material used in this study has been acquired via mutual transfer agreements (MTAs and/or SMTAs) and it is compliant with the Nagoya Protocol (EU regulation 511/2014).


This research was made possible by support of the Dutch Technology Foundation TTW, which is part of the Netherlands Organisation for Scientific Research (NWO), and which is partly funded by Ministry of Economic Affairs (NWO/TTW VIDI grant 14854 to N.E.F.).

Author information

Authors and Affiliations



NB, LC, KB, BJZ, MES, GB and NEF designed the experiments. NB, LC and GBK performed the experiments. NB and JV developed the image-based phenotyping method. ME-S prepared the genetic map. NB analysed the data. NB wrote a first draft of manuscript. LC, ME-S, KB, BJZ, MES, GB and NEF critically commented and revised it. All authors read and approved the manuscript.

Corresponding author

Correspondence to Nina E. Fatouros.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1. Supplementary Figure S1.

Germplasm screening of 56 B. rapa accessions reveals variation in HR-like cell death. Supplementary Figure S2. Image-based phenotyping protocol developed to quantify HR-like cell death size on leaves. Supplementary Figure S3. Replicability of two image segmentation methods to quantify HR-like cell death size. Supplementary Figure S4. Genetic linkage map for the L58 x R-o-18 RIL population. Supplementary Figure S5. Quantitative trait loci for HR-like cell death size in the L58 x R-o-18 RIL population. Supplementary Figure S6. Heatmap of genome-wide LOD scores of two-QTL models to investigate epistatic interaction and additive effects. Supplementary Figure S7. Phenotypic distribution of P. brassicae egg-induced cell death (Pbc) in twelve selected RILs of B. rapa to validate QTL effects. Supplementary Figure S8. Synteny analysis between B. rapa quantitative trait loci for cell death size and A. thaliana.

Additional file 2. Supplementary Table S1.

List of B. rapa accessions used for the germplasm screening. Supplementary Table S2. Evaluation of 56 B. rapa accessions for P. brassicae egg wash-induced cell death. Supplementary Table S3. Summary statistics of B. rapa homozygous lines re-evaluated for egg- and egg wash-induced cell death. Supplementary Table S4. List of SNPs that were used to design PCR markers and construct a genetic map. Supplementary Table S5. List of InDels that were used to design PCR markers and construct a genetic map. Supplementary Table S6. Genetic linkage map constructed with the L58 x R-o-18 RIL population. Supplementary Table S7. Summary of the L58 x R-o-18 RIL population genetic map. Supplementary Table S8. Candidate genes related to plant immunity found within the region of QTLs Pbc1 (A02). Supplementary Table S9. Candidate genes related to plant immunity found within the region of QTLs Pbc2 (A03). Supplementary Table S10. Candidate genes related to plant immunity found within the region of QTLs Pbc3 (A06).

Additional file 3. Supplementary Table S11.

Phenotypic data of the 56 B. rapa accessions tested for the germplasm screening. Supplementary Table S12. Phenotypic data of B. rapa homozygous lines re-evaluated for egg- and egg wash-induced cell death. Supplementary Table S13. Phenotypic data of parents and L58 x R-o-18 RIL population across three repeated QTL experiments. Supplementary Table S14. Marker data used to construct the L58 x R-o-18 RIL population genetic map. Supplementary Table S15. Phenotypic data of twelve selected RIL lines that were re-evaluated for egg-induced cell death.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Bassetti, N., Caarls, L., Bukovinszkine’Kiss, G. et al. Genetic analysis reveals three novel QTLs underpinning a butterfly egg-induced hypersensitive response-like cell death in Brassica rapa. BMC Plant Biol 22, 140 (2022).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: