The B. rapa reference genome line Chiifu, characterized by a CI value < 1.0 at the anthesis stage, served as the receptor parent. Its S-haplotype is S⁶⁰ (class II), as detailed in a recent study (Yamamoto et al., 2023). A Yellow Sarson line, R-o-18, exhibiting a CI value > 10, was employed as the donor parent. The initial cross involved mating the donor parent R-o-18 with the recurrent parent Chiifu, resulting in the generation of F₁ plants. Subsequent backcrossing was systematically performed with Chiifu, leading to the development of backcross populations denoted as BC₁, BC₂, and BC₃. From the BC₃ generation, the BC₃-1, BC₃-6, and BC₃-7 plants with the Mm genotype and high proportion of recurrent parent genome (PR) were self-pollinated to yield BC₃S₁ populations. Individuals with both phenotypic and genetic background similarity to Chiifu while demonstrating self-compatibility were selected from the BC₃S₁ generation to be self-pollinated to produce self-compatible Chiifu lines ( Figure 1 ). Additionally, this study included 841 accessions of B. rapa from our laboratory to identify individuals with mm genotypes. These accessions encompassed 207 Chinese cabbage, 189 Pak choi (ssp. chinensis), 150 Caixin (ssp. parachinensis), 119 European turnip (ssp. rapa), 76 Oil seeds (ssp. oleifera), 41 Yellow Sarson (ssp. tricolaris), 24 Japanese komatsuna (ssp. perviridis), 23 Japanese mizuna (ssp. nipposinica), and 12 Taicai (ssp. chinensis var. tai-tsai) (Cheng et al., 2016).

Genomic DNA was extracted from all samples by the modified CTAB method (Mavrodiev et al., 2021). Plant materials selected for DNA library construction included R-o-18, BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14. These samples were chosen to represent key stages in the breeding program. For WGR, DNA libraries were prepared and subsequently sequenced on a BGISEQ-T7 platform employing paired-end 150 bp (PE150) reads. To ensure the production of high-quality data, the raw sequencing reads were filtered using the fastp software, adhering to default parameters (Chen et al., 2018).

High-quality reads were aligned to the B. rapa Chiifu reference genome version 3.0 using Burrows-Wheeler-Alignment (BWA) software with the “mem” method (Li and Durbin, 2009; Zhang et al., 2018). PCR-duplicated reads were filtered out using Samtools software (Li et al., 2009). Variant datasets (SNPs/INDELs) were obtained by applying stringent filters to exclude variants with low quality or multi-allelic features (QUAL < 60, QD < 20, MQ < 30, FS > 13) using GATK software (McKenna et al., 2010; Chang et al., 2023). Additionally, heterozygous loci in the variant datasets of R-o-18 were excluded. Sequence variants of the MLPK gene between the reference genome Chiifu v3.0 and R-o-18 aligned reads were visualized using the Integrated Genome Viewer (IGV) software (Robinson et al., 2011). The analysis was conducted with an updated version of genome annotation (Chiifu v3.5) to ensure precision in interpreting the genomic variations (Zhang et al., 2022b).

KASP genotyping was employed to discern specific SNPs associated with the MLPK gene. Prior to designing the MLPK-specific marker, an extensive investigation into SNPs affecting MLPK functionality was conducted. Initially, SNPs within MLPK were extracted from the SNP dataset of 524 B. rapa genomes, as reported in a previous study (Cai et al., 2021). Subsequently, the precise SNP responsible for the G194R mutation within MLPK, which has been reported by Murase et al. (2004), was identified. Utilizing this SNP information, an MLPK-specific marker was designed using Primer3 software (Untergasser et al., 2012). Furthermore, accessions carrying the mm genotype were identified from the SNP dataset of 524 B. rapa genomes, with one such accession serving as a positive control for subsequent KASP genotyping.

Amplification and allelic discrimination for the KASP marker were performed using the Scientific QuantStudio™ 12K Flex Real-Time PCR System (Thermo Fisher Scientific, USA). The KASP genotyping protocol followed the methodology outlined in a previous study (Zhang et al., 2022a). Details regarding the KASP reaction volume and thermal cycling conditions for SC-MLPK are provided in Supplementary Tables S1 and S2 , respectively.

The CI value was determined following a previously published protocol (Xiao et al., 2019). Approximately 20 flowers per individual were manually self-pollinated one or two days after reaching full bloom. At the podding stage, CI was calculated using the formula: CI = number of seeds/number of pollinated flowers. Levels of SI, moderate self-compatibility (MSC), or SC were categorized as follows: SI (CI < 1); MSC (1 ≤ CI ≤ 4); SC (CI > 4). B. rapa accessions with the mm genotype, Chiifu, R-o-18, BC₃S₁ individuals with the mm genotype, were assessed of CI values. Additionally, from the BC₂ generation, plants with Mm genotypes were self-pollinated to generate BC₂S₁ populations. BC₂S₁ individuals with the mm genotype were assessed of CI values. All pollinations were conducted in the experimental greenhouse of IVF-CAAS, Beijing, China, at around 10:00 a.m. on sunny days.

Genome-wide INDELs were selected from variant datasets derived between Chiifu and R-o-18. The Primer3 software was employed to design INDEL markers, utilizing 150 bp sequences flanking each INDEL locus (Untergasser et al., 2012; Chen et al., 2021). The designed INDEL markers produced amplicons ranging from 130 bp to 260 bp, with primer lengths between 18 bp and 23 bp. Melting temperature (Tm) predictions for each amplicon were carried out using the nearest-neighbor model (Breslauer et al., 1986; SantaLucia, 1998), setting parameters such as the concentration of sodium ions and reaction temperature to 1 M and 37 °C, respectively. ΔTm for each INDEL amplicon was predicted, and INDEL markers with a ΔTm greater than 0.5°C were retained as HRM markers.

To address the limitations of physical distance in reflecting genetic recombination rates, HRM markers were anchored to the constructed genetic linkage map in B. rapa based on their physical positions (Liu et al., 2016). Genome-wide HRM markers, evenly distributed across 10 chromosomes of B. rapa, were chosen for HRM genotyping, guided by both physical and genetic maps. Amplification and allelic discrimination for HRM genotyping were performed following a previously reported method (Chen et al., 2021).

WGR data of representative individuals (BC₃-1S₁-6, BC₃-6S₁-15, BC₃-7S₁-14, and R-o-18) were mapped to assembly of Chiifu. Firstly, whole-genome SNP datasets were generated using the Variant Calls method. Loci with missing genotypes were filtered out from the SNP datasets. Genotypes obtained from the SNP datasets were anchored to the physical map based on their respective physical positions. Graphical genotyping analysis was conducted using our Python script available at https://github.com/clc-CAAS/Genetic_map_construction.

PR was calculated after genotyping using the formula: PR = (Recurrent parent allele count * 2 + Heterozygous allele count * 1 + Donor parent allele count * 0)/Total allele count * 2. HGB was calculated after genotyping using the formula: HGB = Homozygous allele count/Total allele count. To distinguish genotyping by HRM markers and WGR data, PR and HGB were designated as PR_m and HGB_m, respectively, when calculated after genotyping by HRM markers, and as PR_s and HGB_s, respectively, when calculated after genotyping by WGR data.

Chinese cabbage flowers from Chiifu, R-o-18, BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14 were manually pollinated with a consistent amount of self-pollen grains one or two days after reaching full bloom. Ten flowers were selected for each individual. After 6 hours, stigmas were cut 2 mm away from the stigmatic surface for aniline blue staining to visualize pollen tubes following a published approach (Zhang et al., 2021). Stigmas were fixed in Canoy’s fixative (methanol: acetic acid = 3:1), softened in 10 M NaOH, and stained with 0.1% aniline blue. Pollen tubes were visualized under UV fluorescence using a Leica MSV269 Stereo Microscope, and images were captured with a Leica DFC450 digital camera. SI, MSC or SC was assessed by the number of pollen tubes penetrating stigmas (NPT). The levels of SI, MSC or SC were classified as follows: SI (NPT < 10), MSC (10 ≤ NPT ≤ 25), SC (NPT > 25).

To identify the SNP responsible for inducing G194R within MLPK, we examined SNP datasets from 524 B. rapa genomes (Cai et al., 2021). We identified a total of 26 SNPs within MLPK ( Supplementary Table S3 ). A single nucleotide substitution (G to C) at position 1277 bp (G1277C) in the fourth exon of MLPK was found to induce the G194R substitution ( Figure 2A ; Supplementary Table S4 ). This G1277C mutation was present in only five out of 524 B. rapa genomes, while the MLPK/MLPK (MM) genotype observed in 509 genomes and missing genotype data in 10 genomes. Additionally, sequence variants of MLPK between Chiifu and R-o-18 were visualized, as illustrated in Supplementary Figure S1 . We developed a specific KASP marker named SC-MLPK based on the G1277C mutation ( Figure 2B ; Supplementary Table S5 ). We optimized the thermal cycling conditions of SC-MLPK for KASP genotyping, determining that 32 cycles after touchdown in the PCR program were superior to the native 26 cycles ( Supplementary Figure S2 ; Supplementary Table S2 ).

To identity accessions with mm genotypes, we genotyped 841 B. rapa accessions by SC-MLPK marker assay. Eight accessions with the mm genotype were identified, including five yellow sarson, one pak choi, two caixin. Notably, the mm genotype was not detected in any accessions of Chinese cabbage ( Figures 2C, D ). To further confirm the self-compatibility of the mm plants, we accessed CI values for these eight accessions ( Table 1 ), revealing CI values ranged from 4.25 to 16.3 among them. These results demonstrated that all of the eight accessions with mm genotypes were self-compatible. Thus, the mm genotype emerged as a potential factor for overcoming self-incompatibility in B. rapa.

To facilitate the introgression of mm allele into Chiifu through MABC, we developed HRM markers which discriminate genome backgrounds between Chiifu and R-o-18. By aligning 10 Gb WGR data of R-o-18 to the reference genome Chiifu v3.0, a total of 665,095 INDELs were identified. Among these, 188,289 INDELs meeting stringent criteria of homozygosity and a length of ≥ 3 bp were selected for marker design. Furthermore, 131 co-dominant HRM markers ( Supplementary Table S6 ) evenly distributed across 10 chromosomes of the B. rapa ( Table 2 ), were successfully developed for MABC. The ΔTm (melting temperature difference) values of INDEL amplicons between Chiifu and R-o-18 ranged from 0.8 °C to 2 °C for these HRM markers. The number of HRM markers ranged from 10 on chromosomes A07 and A08 to 17 on chromosome A03. The average genome-wide density of HRM markers was 2.28 Mb, ranging from 1.48 Mb on A10 to 3.01 Mb on A09. Furthermore, the average genetic interval across the genome was 7.29 cM, with a range from 5.15 cM on A10 to 9.04 cM on A07 ( Table 2 ). The HRM genotyping for a representative HRM marker on each chromosome is depicted in Supplementary Figure S3 . The 131 evenly distributed HRM markers across 10 chromosomes were used for further MABC of mm allele introgression from R-o-18 into Chiifu.

To facilitate the selection of individuals harboring the Mm genotype in successive backcrossing populations, the KASP assay for the SC-MLPK marker was applied to the BC₁, BC₂, and BC₃ populations. The distribution of individuals harboring the Mm genotype in each backcrossing population was as follows: BC₁ (6/12), BC₂ (n=534/1120), and BC₃ (n=390/768). An extensive population comprising 1120 individuals was utilized in the BC₂ population to expedite the selection process for plants with the Mm genotype and a high PR_m. The segregation ratio of MM: Mm in each backcrossing population aligned with the anticipated Mendelian ratio of 1: 1 (P > 0.05).

To streamline the genetic background analysis in the extensive BC₂ population, a two-step MAS strategy was implemented for individuals with the Mm genotype. In the first step, three HRM markers on chromosomes A01, A02, and A09 were used for preliminary selection ( Supplementary Figures S4A-D ). Subsequently, a total of 96 individuals selected from the first round screening were upon the second round selection using the left 128 HRM markers. Based on the result from assays by SC-MLPK and 131 HRM markers, the maximum PR_m in the BC₂ population was determined to be 85.1%. For background analysis of individuals with the Mm genotype in the BC₃ population, 39 markers with heterozygous genotypes from the BC₂ population were employed. The highest PR_m was 97.3% among the BC₃ population ( Supplementary Figure S4D ). These results elucidated that the genetic background of each backcrossing population progressively approached to that of Chiifu.

Three BC₃ individuals with Mm genotypes as well as high PR_m, including BC₃-1 (PR_m = 97.3%), BC₃-6 (PR_m = 95.1%), and BC₃-7 (PR_m = 94.7%), were retained for self-pollination to generate BC₃S₁ populations. From each BC₃S₁ population, ten mm genotype plants were selected to be genotyped using background HRM markers ( Figure 3A ). Genotyping revealed that PR_m for the BC₃S₁ population ranged from 93.9% to 98.5%. Three mm individuals with high PR_m and high HGB_m, namely BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14, were retained from BC₃S₁ population. Genotyping indicated that the PR_m of BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14 was 98.5%, 96.6%, and 95.8%, respectively. The HGB_m of BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14 was 100%, 96.2%, and 99.2%, respectively.

To further confirm the high similarity of genetic backgrounds between the three BC₃S₁ individuals to that of Chiifu, we conducted WGR for them. Genome-wide SNPs were extracted from aligned reads, utilizing the B. rapa genome Chiifu v3.0 as a reference. Loci with missing genotypes for each sample were filtered, and SNPs with homozygous genotypes in R-o-18 were retained, resulting in the identification of 1.6 million SNPs spanning chromosomes A01 to A10. These SNPs were anchored to both physical and genetic maps ( Figures 3B-D ). Based on the SNPs, boundaries of introgressed segments were delineated in BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14 ( Supplementary Table S7 ). Genotyping revealed that the PR_s of BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14 was 98.1%, 92.6%, and 88.5%, respectively. The HGB_s of BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14 was 99.6%, 97%, and 96%, respectively ( Table 3 ). The individual plant (BC₃-1S₁-6) with the highest PR_s of 98.1% in the BC₃S₁ population displayed a near-identical match with the recurrent parent Chiifu genotype, except for a specific region (approximately 5.2 Mb, from 12059865 to 17244940) of the MLPK gene linkage on chromosome A03 retaining the donor parent R-o-18 genotype ( Figure 3B ). This result indicated successful introgression of mm allele from R-o-18 into Chiifu.

To evaluate the phenotypic traits of individual plants within the BC₃-1S₁, BC₃-6S₁, and BC₃-7S₁ populations relative to Chiifu, we conducted a detailed analysis encompassing agronomic traits, pollen tube visualization, and assessment of CI value ( Figure 4 ). For agronomic trait analysis, 10 plants with mm genotype were selected from each BC₃S₁ population. Before leafy heads were ready for harvest, head height and width were evaluated and compared with Chiifu ( Table 4 ). The head height and width ranged from 23.5 to 25.1 cm and 10.7 to 13.0 cm, respectively. Notably, both the average head height and width in BC₃-6S₁ population were not significantly different from those of Chiifu. In each BC₃S₁ population, all individuals exhibited heading, hairy, and crinkle traits similar to Chiifu, while the R-o-18 line displayed non-heading, hairless, and smooth traits. Three individuals including BC₃-6S₁-15, BC₃-1S₁-6, and BC₃-7S₁-14 exhibiting both genotypic and phenotypic similarity to Chiifu were selected from each BC₃S₁ population ( Figures 4A-C ; Supplementary Figure S5 ).

Upon flowering, BC₃-6S₁-15, BC₃-1S₁-6, and BC₃-7S₁-14, together with Chiifu, R-o-18, were investigated for pollen tube elongation and CI. Fluorescence microscopy revealed that R-o-18, BC₃-6S₁-15, and BC₃-7S₁-14 plants had more than 25 pollen tubes, whereas no pollen tubes were observed in BC₃-1S₁-6 and Chiifu ( Figures 4D-F ; Supplementary Figure S6 ). During the seed-filling stage, seed setting was evaluated ( Figures 4G-I ). Notably, BC_3-6S₁-15 and BC₃-7S₁-14 plants exhibited significant SC, achieving CI values of 7.5 and 7.1, respectively. In contrast, BC₃-1S₁-6 exhibited SI with a CI of 0.1, similar to that of Chiifu (CI = 0.1). Consequently, the two newly developed plants (BC₃-6S₁-15, BC₃-7S₁-14), which exhibited phenotypic and genetic similarities to Chiifu while demonstrating SC, were self-pollinated to generate two novel lines designated as SC-Chiifu Lines (SC-Chiifu-a, SC-Chiifu-b, respectively). The CI values of the SC-Chiifu Lines compared to the CI value of Chiifu (CI = 0.1) increased tenfold. However, the CI values of the SC-Chiifu Lines remained lower than that observed for R-o-18 (CI = 12). These results highlight the successful transfer of self-compatibility from R-o-18 into Chiifu, leading to the development of two novel SC-Chiifu Lines.

In the BC₂S₁ population, 63 individuals with the mm genotype were selected for the assessment of CI values. The CI values varied significantly among these individuals, spanning from 0.2 to 13 ( Supplementary Table S8 ). Subsequent CI analysis delineated 29 individuals exhibiting SC (CI > 4), 22 individuals displaying MSC (1 ≤ CI ≤ 4), and 12 individuals manifesting SI (CI < 1). This observation indicated a pronounced segregation of SC within the mm plants of the BC₂S₁ population.

Furthermore, variations in SC were also observed among the three selected individuals (BC₃-6S₁-15, BC₃-7S₁-14, and BC₃-1S₁-6) with the mm genotype. Given that all three plants in the BC₃S₁ population possessed the same mm genotype and similar genetic backgrounds, the observed segregation suggested the involvement of gene(s) beyond the MLPK locus in influencing self-compatibility in B. rapa. We postulated that the gene(s) responsible being located at the segments shared by BC₃-6S₁-15 and BC₃-7S₁-14 individuals are distinct from those in BC₃-1S₁-6. We examined S haplotypes in BC₃-1S₁-6, BC₃-6S₁-15, BC₃-7S₁-14, Chiifu and R-o-18 using S-specific primers (Oikawa et al., 2011). We found that BC₃-1S₁-6, BC₃-6S₁-15, and BC₃-7S₁-14, along with Chiifu, had S⁶⁰ homozygotes, while R-o-18 exhibited an S⁵⁴ homozygote. There was no correlation between CI value and S haplotype in BC₃-1S₁-6, BC₃-6S₁-15 and BC₃-7S₁-14. Based on our analysis on the introgressed segment boundaries in the three selected individuals, we identified four such segments ( Figure 5 ). Specifically, four shared segments were located on chromosomes A03 (approximately 0.8 Mb, from 16451403 to 17244940), A06 (approximately 11.3 Mb, from 8101576 to 19449508), A07 (approximately 0.8 Mb, from 14002870 to 14798378), and A09 (approximately 0.7 Mb, from 12255709 to 12997755). These findings suggested that these shared segments emerged as candidate segments potentially influencing SC.