A total of 225 AtLRR-RLKs protein sequences were download from the TAIR database¹ (Sun J. et al., 2017). Protein sequences of the semi-winter rapeseed cultivar ‘‘Zhongshuang11’’ (ZS11) were download from the National Genomics Data Center (NGDC),² and the accession number was PRJCA002883. The longest transcript of each gene was used for subsequent identification. Hidden Markov Model (HMM) and BLASTP search were used to classify the BnLRR-RLKs. The HMM profiles of LRR domain LRR1 (PF00560), LRR3 (PF07725), LRR4 (PF12799), LRR5 (PF13306), LRR6 (PF13516), LRR8 (PF13855), LRR9 (PF14580), LRR10 (PF18805), LRR11 (PF18831), LRR12 (PF18837), LRRNT-2 (PF08263), and Pkinase domain (PF00069) as well asPkinase_Tyr (PF07714) were obtained from the Pfam database.³

HUMMER3.1 software of Linux version was used to take the intersection of LRR and Pkinase domain, and the threshold was E-value < e^–10. The 225 AtLRR-RLK protein sequences were used as queries to perform a BLASTP search in the local protein database of ZS11 with an E-value < e^–10 and identity ≥ 50% (Camacho et al., 2009), and the putative BnLRR-RLKs were obtained by taking the intersection of the HUMMER and Blastap methods. Finally, these BnLRR-RLKs were submitted to the NCBI-CDD server⁴ and the SMART database⁵ to screen for proteins with the LRR and Pkinase domains. In addition, the theoretical isoelectric points (pI) and molecular weight (MW) of the LRR-RLKs were calculated by using the ExPASy database.⁶

Multiple sequence alignment of the full-length BnLRR-RLK proteins were performed with the FFT-NS-I method of the MAFFT software (Katoh et al., 2019). Fasttree software was used as the maximum likelihood method for constructing phylogenetic trees (Price et al., 2009). The EvolView software was used to visualize the phylogenetic tree (Zhang et al., 2012). Conserved motifs (the maximum is set to 10) were predicted using the MEME v5.2.0 with the maximum number of motifs and the optimum width of each motif falls between 10 and 200 residues (Bailey et al., 2009). The structures of the BnLRR-RLK genes were displayed based on GFF annotation files by TBtools software (Chen et al., 2020).

Chromosome location information of the BnLRR-RLK genes was extracted from the GFF file, and plotted by the MapChart version 3.0 software (Voorrips, 2002). In order to accurately analyze BnLRR-RLK duplication events, the e-value of e^–100 was used to align all of the protein sequences in ZS11 with the BLASTP program (Camacho et al., 2009). The duplication pattern of these genes were detected with the default parameters of the MCScanX software, and divided into tandem duplication and segmental duplication (Krzywinski et al., 2009). Similarly, the e-value of e^–100 was used to align all of the protein sequences between ZS11 and Arabidopsis with the BLASTP program, and all syntenic blocks were mapped with JCVI software (Tang et al., 2008).

RNA-seq data from 12 tissues (root, stem, leaf, flower, silique, sepal, pistil, stamen, ovule, pericarp, wilting pistil, and blossomy pistil) of ZS11 were downloaded from the NCBI database (project ID: PRJNA394926) (Sun F. M. et al., 2017). RNA-seq data of ZS11 leaf under abiotic stress, including dehydration, salt, ABA and cold stress treatments, were downloaded from the NGDC database (project ID: CRA001775) (Zhang et al., 2019). These transcriptome data were mapped to the ZS11 reference genome using HISAT2 software (Kim et al., 2015). The TPM values (Transcripts Per Million) were calculated by FeatureCounts R package and the heatmaps were presented using TBtools software (Liao et al., 2014; Chen et al., 2018).

To construct the CRISPR/Cas9 recombinant plasmid, a common target sequence of six BnBri1 genes was designed (Supplementary Table 1), and this target dsDNA was generated via direct annealing of two oligonucleotides primers, and assembled into the pHSE401 vector (Xing et al., 2014). The resulting pHSE401-BRI1 constructs were transformed into the agrobacterium tumefaciens strain GV3101 by electroporation and subsequently transformed into the selfing line K407 of B. napus using the method described previously (Bhalla and Singh, 2008).

Genomic DNA extracted from T₃ plants was used for PCR amplification with genotyping primers (Supplementary Table 1), the PCR products were conducted by sanger sequencing, and cloned into a pEASY-Blunt cloning vector (TransGen Biotech, Beijing, China), positive single clones were sequenced to identify mutations (Aoke Biotech, Beijing, China).

To analyze the hypocotyl variation of the gene-edited mutants, the seeds were sown in humus soil under 16 h light/8 h dark photoperiod at 22^°C, and the hypocotyl length was measured after 6 days. To verify the response to Brassinozole (BRZ), the seeds of mutants were sterilized (75% alcohol for 30 s, 15% NaClO for 8 min, rinse 5 times with sterile water) and bunch planting at 1/2 strength Murashige and Skoog solid medium (1/2 MS), which contained 0 and 1,000 nM BRZ, respectively (Tokyo Chemical Industry Co., Ltd., Tokyo, Japan). After 5 days of darkness cultivation at 22^°C, the hypocotyl length was measured. Each treatment was repeated three times.

To detect the expression level of BnBRI1, total RNA was extracted from seedlings of the L18, L24 and WT by RNAiso Plus kit (TaKaRa, Dalian, China). One microgram of RNA was reverse-transcribed to cDNA using One-Step gDNA Removal and cDNA Synthesis SuperMix (TransGen Biotech, Beijing, China). qRT-PCR was performed using a SYBR Green Master Mix kit (TaKaRa, Dalian, China). The BnUBC21 gene was used as an internal reference, the primer pairs are listed in Supplementary Table 1. Each sample contains 3 biological replicates and 3 technical replicates.

All data were analyzed with One-way analysis of variance (ANOVA) of SPSS 26 software. GraphPad prism (version 8.0.2) was used for drawing graphs. Each experiment was repeated at least three times, and significant differences are indicated at the significance levels of P ≤ 0.05.

We detected a total of 444 BnLRR-RLK genes in the genome of B. napus cultivar ZS11, which is nearly double the number of LRR-RLK genes in Arabidopsis (Table 1). Specifically, we uncovered 215 LRR-RLK genes in the A subgenome and 225 LRR-RLK genes in the C subgenome. We renamed all of the BnLRR-RLK genes on the basis of chromosomal localization. Gene position and other details are presented in Supplementary Table 2. We found that BnLRR-RLK protein lengths were between 459 and 1,303 amino acids (aa). In addition, the protein molecular weights (MV) ranged from 51.31 to 143.19 kDa, and the isoelectric point (pI) from 4.77 to 10.27.

To evaluate the evolutionary relationships between LRR-RLK genes in B. napus and A. thaliana, we constructed a phylogenetic tree using the 444 BnLRR-RLK full-length proteins together with 225 AtLRR-RLK proteins. The 444 BnLRR-RLKs were subdivided into 22 groups (I–XV) (Figure 1). The largest group, LRR-RLK-III, contained 99 genes. In contrast with A. thaliana, no genes were identified in group I-2. There were at least twice as many LRR-RLK genes in each subgroup of B. napus vs. A. thaliana, except for groups I-1, VI-2, VIII-1, and XV, which had a similar number of LRR-RLK genes in B. napus and A. thaliana. These data suggest that allopolyploidization contributed to the expansion of LRR-RLKs in B. napus, but that not every subfamily has an increased the number of genes.

In addition, we analyzed conserved motifs and gene structures to define BnLRR-RLK family characteristics. Using MEME tools, we identified 10 putative motifs (Supplementary Table 3). We found that the type, number, and arrangement of the N-terminus determined the classification of different subfamilies, and expanded the variety and quantity of LRR-RLKs (Supplementary Figure 1). The N-termini of these proteins were enriched for many different tandem repeat motifs, such as motif 1, motif 3 and motif 4. Intriguingly, motif 2 and motif 5 were often found at the C-terminus.

The gene structure characteristics revealed that different subfamilies shared different intron-exon distribution, whereas members within the same subfamily had a similar gene structure, suggesting a similar evolution within subgroups but not between subgroups, based on GFF annotation files (Supplementary Figure 1). The large number of intron-exons in subgroups I-1, II, XIII, V, and VI-2, indicate a complex gene structure, whereas other subfamilies contained few intron-exons and showed simpler gene structures.

Knowledge of the chromosome distribution of gene family members is essential for studying family member duplication and collinearity. We found that the 444 BnLRR-RLK genes were widely distributed on 19 chromosomes, and the distribution was extremely uneven, with about 40 genes on chromosome 3C and 13 genes on chromosomes 4A and 8A, respectively (Supplementary Figure 2). We did not observe any preference for the distribution of LRR-RLK genes on each chromosome, with the exception that there were no LRR-RLK genes on the upper arm of chromosome 4A.

Gene duplication is a major mechanism underlying the expansion of gene family members. Based on BLAST and MCScanX software, we identified a total of 6 tandem duplication events and 395 segmental duplications in the ZS11 genome (Figure 2 and Supplementary Tables 4, 5). Only 9 events were detected within the same chromosome, whereas 386 segmental duplications occurred across chromosomes, which suggests that segmental duplication events across chromosomes plays a key role in LRR-RLKs expansion. Further analysis showed that 51 duplication events occurred on the AA subgenome, 55 events on the CC subgenome, whereas 289 events occurred across AA/CC subgenomes. Considering that B. napus (AACC genome) derived from a spontaneous hybridization between B. rapa (AA genome) and B. oleracea (CC genome), and contains the complete diploid chromosome sets of both parental genomes. Therefore, we infer that allopolyploidization plays a significant role in the expansion of the LRR-RLKs in B. napus.

To investigate the evolution of the LRR-RLKs, we analyzed their synteny between B. napus and A. Thaliana at the whole genome level. A total of 374 gene pairs were detected between the two genomes (Supplementary Figure 3 and Supplementary Table 6). Most AtLRR-RLK genes have multiple orthologous genes in B. napus. For instance, AtBRI1 (AT4G39400) corresponds to 6 BnBRI1 genes, indicating that BnBRI1 genes expanded in rapeseed. However, the 29 AtLRR-RLKs have only one collinear BnLRR-RLK gene.

The specific expression of genes in different tissues sheds light on gene function. In order to explore the function of the BnLRR-RLKs, we investigated their expression profiles in 12 different tissues, using public available RNA-seq data. We detected the expression of 426 BnLRR-RLKs in at least one tissue, whereas 18 genes were not expressed in any of the tissues analyzed (Figure 3 and Supplementary Table 7). We divided the 426 genes into two groups: group 1 displayed high expression levels in different tissues, whereas group 2 displayed low expression in tissues. For instance, BnLRR-RLK-6A-15 was highly expressed in ovule, BnLRR-RLK-8A-12 and BnLRR-RLK-8C-6 were highly expressed in root. These results reflected that the expression of LRR-RLKs vary in different tissues, suggesting that BnLRR-RLKs play functionally diverse roles in tissue development.

LRR-RLK genes play a major role in abiotic stress responses. To explore the expression of LRR-RLK genes in response to abiotic stress, we examined their expression patterns under four abiotic stresses (dehydration, salt, abscisic acid (ABA), cold) using the published transcriptome data of B. napus. A total of 132 genes were extracted from the expression matrix based on STDEV value > 2 among different treatments (Figure 4 and Supplementary Table 8). These genes clustered into two groups.

Group 1 contained 41 genes with significantly different expression patterns under different stress treatments. For example, BnLRR-RLK-4A-8, BnLRR-RLK-3C-3, and BnLRR-RLK-3A-8 were all up-regulated after 4 h-salt or 4 h-cold treatment, suggesting they may be involved in the related stress response. Interestingly, some genes displayed opposing expression profiles under different treatments. For instance, BnLRR-RLK-9A-6 was sharply down-regulated after 1 h-dehydration, but up-regulated after 4 h-ABA treatment. These data suggest that the responses to dehydration and ABA stress have opposing molecular mechanisms.

Group 2 contained 91 genes with low expression levels under each stress treatment. Unusually, part of the genes exhibited sharply changes in expression levels after a specific stress treatment, for instance, BnLRR-RLK-4C-17 and BnLRR-RLK-5A-10 were significantly up-regulated after 4 h-cold treatment, BnLRR-RLK-7A-5 and BnLRR-RLK-8A-2 were up-regulated after 4 h-ABA treatment, whereas BnLRR-RLK-1A-4 and BnLRR-RLK-3C-11 were strikingly down-regulated under dehydration treatment. The results provided useful information for studying the function of LRR-RLKs in response to abiotic stress.

Based on previous research in Arabidopsis and rice, we identified six BRI1 orthologs (renamed as BnBRI1.a∼BnBRI1.f), and six BRI1-Like (BRL) genes in B. napus (Figure 5A and Supplementary Table 9). The BnBRI1 and BnBRL genes both belonged to the Xb-1 subfamily, while the BnBRI1.b and the BnBRI1.f had higher homology with AtBRI1. Clustering analysis based on expression between 12 tissues showed that the six BnBRI1 genes clustered into one subclass and six BnBRL genes into another subclass, and the expression levels of BnBRI1 genes were generally higher than the BnBRL genes (Figure 5B). Furthermore, BnBRI1.c, BnBRI1.e, and BnBRI1.f had higher expression in various tissues compared to other BnBRI1 genes. It is noteworthy BnBRI1.f gene was highly expressed in all tissues.

To generate new germplasm resources with a dwarf and optimized plant structure of B. napus, we designed single guide RNAs targeted to a common conserved region of BnBRI1 genes (Figure 5C), inserted these sgRNA into the pHSE401 expression vector and transformed into hypocotyls of rapeseed cultivar K407 (Figure 5D). We obtained a total of 130 transgenic positive plants (Figure 5E). Interestingly, more than 34 T₀ plants displayed retarded growth with dark green rolling leaves at the stage of seedling growth (Supplementary Figure 4), consistent with phenotypes of previously reported classical bri1 mutants (Huang et al., 2021). About 42 plants (44%) exhibited an extreme dwarf phenotype (plant height 30–60 cm), and 45 plants (47%) showed dwarf phenotypes (61–80 cm), and 9 plants (9%) showed semi-dwarf phenotypes (81–120 cm) (Figure 5F).

To determine whether the transgenic lines had gene editing events, we selected the T₃ generation plants of the semi-dwarf lines L18 and L24 for subsequent experiments (Figure 6A). The sequencing results of the single clones and PCR products near the sgRNA-targeted sites indicated that all target sites were edited in L24, except the BnBRI1.b and BnBRI1.f genes, similarly, all target sites were edited in L18, except the BnBRI1.f gene. The variations between the two lines were single base insertions of A, T or C (Figures 6B,C and Supplementary Figures 5, 6).

Furthermore, qRT-PCR was conducted (Supplementary Figure 7), and showed that the BnBRI1 expression of the L18 and L24 were both reduced significantly. As a result, the hypocotyl length of both mutants shortened obviously (Supplementary Figures 8A,C). The results of the BRZ treatment showed that the hypocotyl length of the two mutants were significantly decreased under 1,000 nM, so the L18 and L24 were BRZ sensitive mutants (Supplementary Figures 8B,D). These verified that the both gene edited lines are indeed related to brassinostreoid signaling transduction pathway.

Agronomic traits were closely related to agricultural production. We planted the L18 and L24 lines in a greenhouse, and investigated the major agronomic traits, plant height, effective branch height, secondary branch number, silique number per plant, silique number of main inflorescence, length of silique and yield per plant (Figures 6C–I). We found that almost all traits exhibited significant differences between the control WT and the two gene edited lines. Specifically, compared with the control, L24 had fewer secondary branches, fewer silique number of main inflorescences, fewer siliques in the whole plant, and shorter siliques, resulting in a significantly lower yield than the control. In contrast, when the plant height and the branch length of L18 were significantly diminished compared to the control, the yield per plant was not reduced. Among yield-related traits, although L18 had shorter silique and a lower silique number of main inflorescence, the increased number of secondary branches increased the total number of silique and prevented a decline in yield.