Mapping of a Pale Green Mutant Gene and Its Functional Verification by Allelic Mutations in Chinese Cabbage (Brassica rapa L. ssp. pekinensis)

Leaves are the main organ for photosynthesis, and variations in leaf color affect photosynthesis and plant biomass formation. We created two similar whole-plant pale green mutants (pem1 and pem2) from the double haploid (DH) Chinese cabbage line “FT” through ethyl methanesulfonate (EMS) mutagenesis of seeds. Photosynthetic pigment contents and net photosynthetic rates were significantly lower in the mutants than in the wild-type “FT,” and the chloroplast thylakoid endomembrane system was poor. Genetic analysis showed that the mutated phenotypes of pem1 and pem2 were caused by a single nuclear gene. Allelism tests showed that pem1 and pem2 were alleles. We mapped Brpem1 to a 64.25 kb region on chromosome A10, using BSR-Seq and map-based cloning of 979 F2 recessive individuals. Whole-genome re-sequencing revealed a single nucleotide polymorphism (SNP) transition from guanine to adenosine on BraA10g021490.3C in pem1, causing an amino acid shift from glycine to glutamic acid (G to E); in addition, BraA10g021490.3C in pem2 was found to have a single nucleotide substitution from guanine to adenosine, causing an amino acid change from E to lysine (K). BraA10g021490.3C is a homolog of the Arabidopsisdivinyl chlorophyllide a 8-vinyl-reductase (DVR) gene that encodes 3,8-divinyl protochlorophyllide a 8-vinyl reductase, which is a key enzyme in chlorophyll biosynthesis. Enzyme activity assay and chlorophyll composition analysis demonstrated that impaired DVR had partial loss of function. These results provide a basis to understand chlorophyll metabolism and explore the mechanism of a pale green phenotype in Chinese cabbage.


INTRODUCTION
Leaf color mutations that are widespread in the plant kingdom have been identified in many plants, such as maize (Lonosky et al., 2004), rice (Oster et al., 2000;Jung et al., 2003;Lee et al., 2005), barley (Preiss and Thornber, 1995;Rudoi and Shcherbakov, 1998), and Arabidopsis (Carol et al., 1999;Kumar et al., 2000). The typical variation in plant leaf color is yellowing, which has become an important material for studying the mechanism of photosynthesis, chlorophyll biosynthetic pathways, chloroplast development, and its genetic control mechanisms (Hansson et al., 1999;Zhao et al., 2008). The color of plant leaves is affected by genetic and environmental factors, and is determined by various pigments, including chlorophyll, carotenoids, and anthocyanins . Chlorophyll, the main pigment in photosynthesis, can absorb most of the red and purple light, and thus plays a key role in light absorption during photosynthesis and affects plant yield to some extent (Fromme et al., 2003). However, chlorophyll biosynthesis is a complex process involving many compounds. The chlorophyll biosynthetic pathway in higher plants consists of 16 enzymatic reactions involving 15 enzymes (Beale, 1999;Wang et al., 2010b). Chlorophyll is produced from L-glutamate in a sequential manner through the formation of protoporphyrin IX, Mg-protoporphyrin IX, protochlorophyllide a, and, finally, chlorophyll a and chlorophyll b via the tetrapyrrole biosynthetic pathway (Tanaka et al., 1998;Beale, 1999;Lange and Ghassemian, 2003). In Arabidopsis thaliana, the chlorophyll biosynthetic pathway is encoded by 27 genes; all of which have been cloned (Beale, 2005;Nagata et al., 2005). Many yellow mutants are generated by mutations in the genes involved in the chlorophyll biosynthesis pathway, which encode key enzymes, such as glutamyl tRNA reductase (GluTR), Mg-chelatase, and chlorophyll a oxidase (CAO) (Oster et al., 2000;Kim et al., 2005;Lee et al., 2005;Gao et al., 2016;Fu et al., 2019).
Classified on the number of vinyl side chains, there are two types of chlorophyll: monovinyl chlorophyll (MV-Chl) and divinyl chlorophyll (DV-Chl) (Nagata et al., 2005;Wang et al., 2010b), which are both accumulated in various organisms (Kim and Rebeiz, 1996;Rebeiz et al., 2003;Islam et al., 2008). To date, five DV-Chl derivatives have been identified, namely, DV-Mg-Proto (Kim and Rebeiz, 1996), DV-Mg-Proto IX ester (MPE) (Kolossov et al., 2006), DV-Pchlide a (Tripathy and Rebeiz, 1988), DV-Chlide a Rebeiz, 1992, 1995), and DV-Chl a (Adra and Rebeiz, 1998). The DV-Chl derivative is formed first, and then the MV-Chl derivative is converted by 3,8divinyl protochlorophyllide a 8-vinyl reductase (DVR) to the corresponding DV-Chl derivative (Nagata et al., 2005;Wang et al., 2010b). Almost all photosynthetic aerobic organisms in the biosphere use MV-chlorophyll to absorb light energy for photosynthesis (Chisholm et al., 1992;Porra, 1997). Only a few species (Prochlorococcus marinus) with a distinct photosynthetic antenna system can use DV-chlorophyll for photosynthesis (Partensky et al., 1999;Ito et al., 2008). Phylogenetic and genomic findings suggest that the progenitors of Prochlorococcus lost the DVR gene and acquired the ability to utilize DV-chlorophyll during evolution (Nagata et al., 2005). The DVR gene has been identified in only a few species to date, including A. thaliana, rice, maize, cucumber, and green sulfur bacteria. Nakanishi et al. (2005) and Nagata et al. (2005) independently identified DVR in chlorophyll content-reduced mutants as AT5g18660, which is the last identified gene required for chlorophyll biosynthesis in higher plants. However, a DVR gene in Brassica plants has not been identified.
Chinese cabbage (Brassica campestris L. ssp. pekenensis) is the main winter vegetable in northern Asia. Leaf color is related to photosynthesis and product quality of Chinese cabbage, so breeders have always focused on its leaf color variations. However, there are few reports about the cloning and identification of the leaf color gene in Chinese cabbage.
We identified two similar pale green mutants (pem1 and pem2) from a Chinese cabbage mutant library through ethyl methanesulfonate (EMS) mutagenesis of germinating seeds. Genetic and photosynthetic characteristics of the mutants pem1 and pem2 had been analyzed. BraA10g021490.3C, a homolog of Arabidopsis DVR, encoding 3,8 divinyl protochlorophyllide a 8vinyl reductase, was considered as a responsible gene of Brpem1 and Brpem2. The mutants could be used as material to investigate chlorophyll metabolism and explore the mechanism that leads to the pale green phenotype in Chinese cabbage.

Plant Materials
The mutant pem1 and pem2 lines were obtained from the "FT" line, which is a double haploid (DH) Chinese cabbage generated via microspore cultures after treating the seeds of the "FT" Chinese cabbage line with EMS. The mutants were crossed with wild-type "FT" to construct F 1 , F 2 , and BC 1 populations for genetic analysis, and with K23, a Pak choi inbred line, which has great difference with the mutant in characteristics, to construct F 2 segregation populations for gene mapping. All plants were grown in a greenhouse or open ground at Shenyang Agriculture University.

Determination of Photosynthetic Pigments
Chlorophyll was extracted from fresh fifth true leaves (0.1 g), using 80% (V/V) ethanol in acetone in darkness for 24 h. Thereafter, the absorbance of extract was measured at 470, 645, and 663 nm, using a DU800 Ultraviolet Spectrophotometer (Beckman Coulter, La Brea, CA, USA). The control was 80% (V/V) ethanol in acetone. The total chlorophyll, chlorophyll a, chlorophyll b, and carotenoid contents in the plants were calculated as described by Holm (1954). Three biological repeats were determined for each material, and three times the technical repeats were performed in each biological repeat.

Determination of Photosynthetic Characteristics
The photosynthetic characteristics of the fifth true leaves of mutant and wild-type (WT) plants at the same growth stage were assessed during a clear morning. Net photosynthetic rates (PN; µmol CO 2 m −2 s −1 ), stomatal conductance (GS; mol H 2 O m −2 s −1 ), intercellular carbon dioxide concentrations (CI; µmol CO 2 mol −1 ), and transpiration rates (TS; mmol H 2 O m −2 s −1 ) were measured, using an Lincoln,NE,USA). Three biological repeats were determined for each material, and three times the technical repeats were performed in each biological repeat.

Observation of Chloroplast Ultrastructure
Fifth true leaves of 6-week-old plants (∼5 × 1 mm) were fixed in 3% glutaraldehyde at 4 • C overnight, washed four times with a 1% phosphoric acid buffer, and then placed in osmium tetroxide (1%) for 2 h. After three washes with the 1% phosphoric acid buffer, the tissues were dehydrated, using graded series of 50, 70, and 80% ethanol, followed by 90 and 100% acetone for 24 h, and then epoxy resin Epon812 was used for encapsulation and polymerization. The tissues were sectioned, using a Leica EM UC7 Ultra-Thin Microtome (Leica Microsystems GmbH, Wetzlar, Germany) and stained with uranyl acetate and lead citrate. Chloroplast ultrastructure was examined, using an HT7700 transmission electron microscope (Hitachi Ltd., Tokyo, Japan) at Analysis and Test Center of Shenyang Agricultural University.

BSR-Seq (Bulked Segregant RNA-Sequencing)
The F 2 population derived from Chinese cabbage pale green mutant pem1 crossed with Pak choi K23 was used to construct mutant-phenotype pool (MPpool) and wild-type phenotype pool (WPpool) for BSR-Seq. Total RNA was, respectively, extracted from the third true leaves of pale green and green plants, using the plant total RNA extraction kit (Tiangen Biotech, Beijing, China).
The two mixed RNA pools were sequenced, using an Illumina HiSeq (Illumina Inc., San Diego, CA, USA). The quality of the original transcriptome sequencing data was evaluated, using fastQC (version 0.10.1) and preprocessed, using Cutadapt (version 1.9.1) to screen out low-quality data and remove contaminants and joint sequences. Subsequently, screened data were mapped on the reference genome (http://brassicadb.org/ brad/datasets/pub/Genomes/Brassica_rapa/V3.0/) using HISAT (v2.0.14). Transcripts predicted by StringTie (v1.0.4) were classified by ASprofile (v1.0.4), and their expression was analyzed. Finally, transcriptional assembly and single nucleotide variants (SNV) were analyzed. Mutant loci with sequencing coverage greater than 3X in the MP-pool and WP-pool were screened, and the Euclidean distance 5 (ED 5 ) values of SNV were calculated to eliminate background noise. According to the ED 5 value to rank, the top 1% was taken as the threshold to identify chromosome regions linked to the target traits (Tan et al., 2019). The raw sequence data of the BSR-Seq are available from the Sequence Read Archive at NCBI under the accession No. SRR14339734.

Whole-Genome Re-Sequencing
The total genomic DNA of wild-type "FT" and mutant pem1 was extracted, using a DNA Secure Plant Kits (Tiangen, Beijing, China). A DNA library was constructed with an insert fragment of 400 bp, and the paired-end (PE) of the library was sequenced, using next-generation sequencing (NGS) based on a NovaSeq 6000 System sequencer (Illumina, San Diego, USA). The coverage of whole-genome re-sequencing was 25x, and the sequencing mode was paired-end, 2 × 150 bp. Highquality data were generated from the raw data by removing 3' end joint contamination, using Adapter Removal (version 2); quality filtering, using the sliding window method; and length filtering. By mapping to the reference genome, we detected SNP (single nucleotide polymorphism) and INEDL (insertion and deletion), CNV (copy number variation), SV (structural variation of chromosome), using GATK (Zhu et al., 2015), CNVnator v0.2.7 (Abyzov et al., 2011), Breakdancer 1.3.7 (Fan et al., 2014), respectively. The SNP/INDEL loci were annotated, using ANNOVAR (Wang et al., 2010a). The raw sequence data of the whole-genome re-sequencing are available from the Sequence Read Archive at NCBI under the accession numbers SRR14340669 and SRR14340670.

DNA Extraction, PCR, and Marker Development
Total genomic DNA extracted from plant leaves, using the modified CTAB method, was amplified, using Mastercycler R nexus GSX1 (Eppendorf AG., Hamburg, Germany). The primers used for PCR amplification were shown in Supplementary Table S2. The reaction volume was 10 µL, and the conditions were 95 • C, 5 min, followed by 35 cycles of 95 • C for 30 s; 56 • C for 30 s, then final extension at 72 • C for 5 min.
The sequences of BSR-Seq target regions were downloaded from the Brassica database. The SSR and INDEL primers were designed, using Primer Premier 5.0 software (Premier Inc., Charlotte, USA) and synthesized by GENEWIZ (Tianjin, China). Polymorphisms in the primers were screened between the parent mutants pem1 and K23, using polyacrylamide gel electrophoresis.

Linkage Analysis and Genetic Map Construction
The recombinants from recessive individuals in the F 2 segregation population were screened, using polymorphic primers between parents to analyze linkage. According to the recombination frequency between polymorphic markers and mutant genes, a genetic linkage map between molecular markers and mutant genes was constructed, using Mapmaker V. 3.0 (Whitehead Institute, Cambridge, MA, USA).

Clone Sequencing
The full-length and promoter sequences of candidate genes were obtained by cloning and sequencing. Primers were designed according to the gene sequence information and are shown in Supplementary Table S5. PCR products were purified, using a Gel Extraction Kit (Omega Bio-Tek Inc., Norcross, GA, USA), then ligated to pGEM R -T Easy Vector (Promega, USA) and transformed into Top 10 competent cells (CWBIO, Beijing, China). The positive clones were screened as blue and white colonies in the Xgal-IPTG medium and sequenced by GENEWIZ (Tianjin, China), using 3730 × lDNA Analyzer (Applied Biosystems, Foster City, CA, USA). Sequences were aligned, using DNASTAR (DNASTAR, Inc., Madison, WI, USA).

Bioinformatic Analysis of BrDVR
The sequence and the information of BrDVR were downloaded from the Brassica database. The conserved domain of BrDVR was analyzed online at NCBI (https://www.ncbi.nlm.nih.gov/ Structure/cdd/wrpsb.cgi). Protein localization was predicted, using the TargetP-2.0 (http://www.cbs.dtu.dk/services/TargetP/) server and chloroplast transit peptides (cTP) in protein sequences, and the locations of potential cTP cleavage sites were predicted, using ChloroP (http://www.cbs.dtu.dk/services/ ChloroP/). Three-dimensional structures of the proteins were built, using SWISS-MODEL (https://swissmodel.expasy.org/ interactive/). The amino acid sequence and accession numbers of BrDVR homologs were obtained from NCBI BLAST. Sequence alignment was performed, using DNAMAN V6 (Lynnon BioSoft, Canada). A phylogenetic tree was constructed with MEGA-X (Kumar et al., 2018), using Clustal W and neighbor joining based on 1,000 bootstrap replications.

Verification of Gene Expression Using qRT-PCR
The expression of candidate genes was detected, using quantitative real-time (qRT)-PCR. Total RNA was extracted from the leaves of mutant and wild-type plants, using the plant total RNA extraction kit (Tiangen, Beijing, China). Single-stranded RNA was reverse transcribed into cDNA, using FastQuant RT Super Mix (Tiangen, Beijing, China); gene-specific primers designed using Primer Premier 5.0 (Supplementary Table S6) were used for qRT-PCR with SYBR Green PCR Master Mix (Takara Bio Inc., Kusatsu, Japan). The reaction mix (20 µL) included 10 µL of 2 × UltraSYBR mixture, 0.8 µL of diluted cDNA, 0.4 µL of forward and reverse primers, and 8.4 µL of double-distilled H 2 O. The ACTIN gene was used as internal standard, and the reaction proceeded under the following conditions: 95 • C for 10 min, followed by 40 cycles at 95 • C for 15 s, 60 • C for 1 min, 95 • C for 15 s, 60 • C for 1 min, 95 • C for 15 s, and 60 • C for 15 s. All reactions included three biological and two technical replicates.

Enzyme Activity Assays
DVR activity was evaluated, using the Plant DVR R ELISA Kit (Meimian Biotech Co., Ltd., Jiangsu, China) via a double antibody sandwich method. Fresh leaf samples of 6-weekold (0.1 g) were homogenized in 1 ml of PBS (pH 7.4) and centrifuged at 12,000 × g for 20 min; the supernatant was then collected. The kit assay DVR activity in the sample used Purified Plant DVR antibody to coat microtiter plate wells, made solid-phase antibody, then added a DVR sample to the wells, combined antibody, which, with HRP labeled, became 1 | Content of photosynthetic pigments in pem1 and pem2 compared with wild-type "FT" (mgg-1).

Material
Chla (ratio of mutant to "FT") Chlb (ratio of mutant to "FT") Total Chl (ratio of mutant to "FT") Car (ratio of mutant to "FT") antibody-antigen-enzyme-antibody complex, and, after washing completely, added TMB substrate solution. TMB substrate became blue color at HRP enzyme-catalyzed, reaction was terminated by the addition of a sulfuric acid solution, and the color change was measured, using spectrophotometer at a wavelength of 450 nm. The enzyme activity of DVR in the samples was then determined by comparing the O.D. of the samples to the standard curve. Three biological repeats were determined for each material, and three times the technical repeats were performed in each biological repeat.

HPLC Analysis for Chlorophyll Composition
Chlorophyll used for HPLC analysis was extracted from fresh leaf samples of 6-week-old pem1, pem2, and "FT" plants, using 100% acetone. After centrifuged at 12,000 × g for 15 min, the supernatant was dried under N 2 gas, redissolved in acetone, and filtered with 0.22 µm organic filter membrane. HPLC analysis was performed, using an Alliance e2695 instrument with a C8 column (4.6 mm i.d. × 250 mm long; 5 µm; Agilent) according to the method of Zapata et al. (2000), using an isocratic system, consisting of methanol: acetonitrile: acetone at 1:3:1 at 40 • C and a flow rate of 1 m min −1 . Eluted samples were monitored, using a 2,489 UV/Vis detector at 660 nm. Chl a and b standards (Sigma, St. Louis, USA) were used as controls. Three biological repeats were determined for each material, and three times the technical repeats were performed in each biological repeat.

Statistical Analysis
The significance difference between the mutants and wild type "FT" was analyzed by Analysis of Variance (ANOVA) of software SPSS (IBM, Armonk, USA) at the significance level of 0.05.

Mutant Creation and Phenotypic Characterization
In order to further study the functional genome of Chinese cabbage, DH line "FT, " which was generated by microspore culture as the WT to create a Chinese cabbage mutant library. Treated 7,800 "FT" germinated seeds in 0.8% EMS for 12 h, to obtain 701 mutant M2 lines, including 286 yellow mutants. Among which, seven mutants with similar phenotype were investigated, and their allelism test (see as the section of "The Allele Test of Brpem1 and Brpem2") revealed that the mutant genes of pem1 and pem2 were allele, and pem1 and pem2 were chosen for further investigation. The mutant lines pem1 and pem2 exhibited the plant pale green phenotype since the cotyledon stage persisted until the heading stage (the leafy head formation stage) (Figure 1). The pale green phenotype was evident in all the aboveground organs. The mutant plants were weak and grew so slowly that they could not head during the heading stage; however, pollination and seed setting were normal.

Analysis of Photosynthetic Pigments
The contents of chlorophyll a, chlorophyll b, total chlorophyll, and carotenoids were significantly decreased in pem1 and pem2 compared with "FT, " and the total Chl contents were, respectively, decreased by 57.49 and 54.54% ( Table 1). The ratios of chlorophyll a content to chlorophyll b content (Chla/b) and of carotenoids content to total chlorophyll content (Car/Chl) were significantly increased (Table 1). Therefore, the chlorophyll content showed a higher decrease than the carotenoid content, and the decrease was greater for chlorophyll b content than chlorophyll a content. Previous studies have shown that chlorophyll is responsible for the green pigmentation in plants, and the presence of carotenoids affects leaf color (Wu et al., 2018;Sun et al., 2020). We speculated that this significant change in the photosynthetic pigment was the direct cause of the pale green phenotype in pem1 and pem2.

Defective Chloroplast Structure and Inefficient Photosynthesis
Chloroplast development in pem1 and pem2 was investigated, using transmission electron microscopy (TEM). The chloroplasts in "FT" had a clear internal structure, with a well-developed somatic membrane system, and many closely arranged stacks of grana with rich lamellae (Figure 2A), whereas those in pem1 and pem2 were slightly smaller and had a poor internal structure, singular grana lamellae, and no distinct grana stacks (Figures 2B,C). Therefore, we speculated that the mutation affected the chloroplast development.
Photosynthetic parameters significantly differed between "FT" and the two mutant lines. Compared with "FT, " the net photosynthetic rate of pem1 and pem2 decreased significantly by 55.65 and 52.67%, respectively ( Figure 2D). The stomatal conductance and transpiration rates of the mutants did not significantly change, while the intercellular CO 2 concentrations of pem1 and pem2 significantly increased by 69.34 and 65.04%, FIGURE 2 | Ultrastructure of chloroplasts and photosynthetic characteristics in pem1 and pem2 compared with "FT" in Chinese cabbage. Mesophyll cells of "FT" (A), pem1 (B), and pem2 (C) were visualized with TEM (magnification, × 50,000). cp, chloroplast; gl, grana lamella; is, intercellular space; m, mitochondria; pg, plastoglobule; s, stroma; v, vacuole. Scale is shown at the bottom right. Photosynthetic characteristics: net photosynthetic rates (D), stomatal conductance (E), intercellular carbon dioxide concentrations (F), and transpiration rates (G) of pem1 and pem2 compared with the wild-type "FT." *indicate significant differences at 5% (ANOVA).
respectively (Figures 2E-G). This indicated that the mutant had a weak ability to fix CO 2 , which might be related to the lower chlorophyll content.

Analysis of Genetic Characteristics
We constructed seven generation populations to clarify the genetic characteristics of mutant pem1. Table 2 shows stable pale green and the green phenotype in pem1 and "FT, " respectively. The pem1 mutant and "FT" reciprocal cross F 1 plants were green, indicating that the pale green mutant features were due to an inherited nuclear recessive gene. The F 2 population showed character separation, and the segregation ratio was 3:1 for green:pale green plants (χ 2 < χ 0.052 = 3.84), which corresponded to the Mendelian laws. Moreover, the BC 1 (F 1 × "FT") population was green, whereas the BC 1 (F 1 × pem1) population had the character separation of 1:1 for green:pale green plants (χ 2 < χ 0.052 = 3.84). These results confirmed that the FIGURE 3 | Morphological appearance of allele mutant lines (pem1 and pem2). Left to right: 1 and 2, "FT"; 3 and 4, pem1; 5 and 6, pem1 × pem 2; 7 and 8, pem2 × pem1; 9 and 10, pem2. mutation phenotype of pem1 was controlled by a single nuclear recessive gene. Table 2 also shows the population phenotype of pem2. The results indicate that the pale green phenotype of pem2 was caused by the recessive mutation of a single nuclear gene.

The Allele Test of Brpem1 and Brpem2
Among 286 yellow mutants in the mutant library, we found seven mutants with a similar phenotype and inheritance. In order to test its allelism, we carried out semi-diallel crossing experiments. If the phenotype of the hybrid between two parents remains the mutant phenotype, the two parents are allele mutation. It was found that the hybrid between pem1 and pem2 showed pale green (Figure 3), which confirmed that the mutant genes of pem1 and pem2 were allele.

Brpem1 Location on Chromosome A10 Determined by BSR-Seq
We constructed two mixed pools with extreme phenotypes (MP-pool and WP-pool) to identify possible regions of pem1, using BSR-Seq. We mapped 110,869 SNV to the Chiifu reference genome that differed between MP-pool and WP-pool. A distribution map of the ED 5 value of different SNV on chromosomes revealed a peak in chromosome A10 (Figure 4).
Two intervals were obtained on chromosome A10 (A10: 7.55-9.56, 10.71-15.35 Mb; Supplementary Table S1) according to the top 1% ED 5 value threshold. Based on the reference genome and the BSR-Seq candidate region, 40 pairs of SSR markers were developed. Among these, 13 markers were polymorphic between mutant line pem1 and Pak choi K23 and produced stable, clear amplification bands (Supplementary Figure S1). Thirty recessive plants were identified in the F 2 population to verify whether chromosome A10 is linked to the pale green trait. We found that SSRA1-1 and SSRA2-16 were linked to and located on both sides of Brpem1 (Supplementary Table S2 and Supplementary Figure S2). Therefore, Brpem1 was mapped to an approximately 1.7-Mb region on chromosome A10 between SSRA1-1 and SSA2-16 at estimated genetic distances of 8.27 and 1.48 cm, respectively.

Fine Mapping of Brpem1 by Linkage Analysis
The mapping populations were expanded to 979 individuals, and more SSR and INDEL markers between SSRA1-1 and SSA2-16 were developed to map Brpem1 more precisely. The distribution of recombinant individuals showed that SSRA5-1, SSRA7-1, INDEL-D5, INDEL-N2, and INDEL-N14 were located on the same side of Brpem1 with SSRA1-1, but SSRA6-21, SSRA6-18, and INDEL-I8 were located on the same side as SSRA2-16 ( Figure 5). In addition, both INDEL-N14 and INDEL-I8 were more tightly linked with Brpem1 at a genetic distance of.05 cm. We mapped Brpem1 on a 64.25-kb region, containing 13 genes between the markers INDEL-N14 and INDEL-I8.

Identification of the Candidate Gene for Brpem1 Based on Whole-Genome Re-Sequencing
Limited by mapping population size, traditional map-based cloning could not precisely locate candidate genes. The whole-genome re-sequence was conducted to detect variations between "FT" and pem1. After filtering out heterozygous sites, four non-synonymous SNPs were identified in the candidate region between INDEL-N14 and INDEL-I8. All the four were located on BraA10g021490.3C (Supplementary Table S4). Clone sequencing showed that three of the four SNPs were not stable at pem1 and only one SNP invariably existed (Supplementary Figure S3 and Figure 6A). This SNP located on the 1,070th nucleotide of BraA10g021490.3C in pem1 was GA (GGG/GAG), which caused an amino acid transformation from Gly to Glu (G to E) on amino acid residue 357 of BraA10g021490.3C.
According to the genetic map of Brpem1, two recombined individuals between the two most closely linked markers, INDEL-N14 and INDEL-I8, were used to clone sequencing. Sequencing alignment of BraA10g021490.3C in two recombined individuals proved that this SNP co-segregated with the pale green trait ( Figure 6B). The predicted three-dimensional structure of the protein showed that the change in the amino acid in pem1 led to a difference in the three-dimensional protein structures in the mutant site between pem1 and "FT" (Figures 6C,D). Hence, we speculated that BraA10g021490.3C was the candidate gene for Brpem1. Figure 6A and Supplementary Figure S3 show the sequencing results of BraA10g021490.3C in pem2. A single nucleotide variation (SNV) from GA (GAG/AAG) was found in nucleotide 1054 of BraA10g021490.3C in pem2, which led to a transformation from Glu to Lys (E to K) in amino acid residue 352 (Figure 6A). The three-dimensional structure of the protein is shown in Figures 6E,F. The amino acid transformation caused by the SNV also led to differences in the three-dimensional structure of the protein at the mutation site. These results suggest that BraA10g021490.3C was probably responsible for the pale green trait of pem2.
biosynthesis. Therefore, we designated BraA10g021490.3C as BrDVR. The single-copy BrDVR in the Chinese cabbage genome has only one 1,233-bp exon, and it encodes a mature protein of 410 amino acids. TargetP and ChloroP predicted that BrDVR has a section chloroplast transfer peptide of 58 amino acids (Supplementary Tables S7, S8). Online NCBI analysis revealed a conserved domain of BrDVR between amino acids 34 and 410 (Supplementary Figure S4A). The amino acid transformations in pem1 (GGG/GAG: G-to-E) and pem2 (GAG/AAG: E-to-K) were both located in the conserved domain of BrDVR.
The amino acid sequences of 11 species were aligned to understand the relationship between BrDVR and other species homologs. The phylogenetic tree showed that BrDVR has significant similarity to B. napus, B. oleracea, R. sativus, and A. thaliana, with identities of 93.04, 89.10, 86.08, and 80.51%, respectively. It was also similar to S. lycopersicum (66.82%), Capsicum annuum (66.36%), Citrus sinensis (68.31%), Populus trichocarpa (67.67%), Oryza sativa (55.68%), and Zea mays (56.84%) (Supplementary Figures S4B, S5). These results indicated that BrDVR is relatively conservative in 11 species rang from dicotyledons to monocotyledons, including fruit trees, vegetables, and crops and more closely related to the DVR of Brassica than to other species.

Expression of BrDVR
We applied qRT-PCR to investigate the spatiotemporal expression of BrDVR in total RNA extracted from the various organs and leaves at different stages of "FT." The results showed that BrDVR was expressed mostly in the leaves, followed by the stems, pods, flowers, buds, and roots ( Figure 7A). The amount of BrDVR expressed in leaves differed among developmental stages, being maximal during the sixth true leaf stage and minimal during cotyledon stages ( Figure 7B). The BrDVR expression level was increased in pem1 and pem2 plants but was not significantly different from that in "FT" plants ( Figure 7C).

Enzyme Activity Assay
DVR activity was measured via enzyme-linked immunosorbent assay (ELISA) to establish whether the SNV in BrDVR affects protein function. The activity of DVR did not significantly differ between pem1 and pem2, but the activity in both mutants was significantly reduced to 46.94 and 48.60%, respectively, that in the wild-type "FT" (Figure 8A).

Chlorophyll Composition Analysis
HPLC analysis found that the chlorophyll (Chl) of pem1 and pem2 differed from those of wild-type "FT" (Figures 8B-D). The peaks corresponding to Chl a and Chl b of the mutants were eluted a little earlier than those of the wild type. Three peaks were observed in both pem1 and pem2, and peak 2' was eluted earlier than peak 2 (chlorophyll a) by approximately 15 s. The characteristics of these peaks are consistent with those of divinyl-chlorophyll (Nagata et al., 2005;Nakanishi et al., 2005;Wang et al., 2010b). The finding indicated that mutants accumulated DV-Chls and a small amount of MV-Chls.

DISCUSSION
Chinese cabbage takes leaves as photosynthetic and product organs, whose leaf color is an important agronomic character of Chinese cabbage. Here, we identified a pair of allele pale green mutants (Figure 3), pem1 and pem2, in a 286 yellow Chinese cabbage mutant library. Based on the BSR-Seq and map-based cloning result, whole-genome re-sequencing indicated that BrDVR (BraA10g021490.3C), a homolog of Arabidopsis DVR, is a candidate gene for Brpem1. Co-segregation analysis and allele mutation analysis results verified that a single (G to A) nucleotide substitution at 1070 of pem1 and 1054 of pem2 caused an amino acid transformation in the conserved domain of the BrDVR protein. This is the first DVR gene identified in Brassica plants. Genetic, phenotypic, physiological characteristics and sequence of pem1 and pem2 will help to understand the biosynthesis of chlorophyll and the FIGURE 8 | DVR activity assay and chlorophyll composition analysis among the wild-type "FT," pem1 and pem2. (A) DVR activity of "FT" and pem1 and pem2. The asterisk represents significance. Chlorophyll composition analysis by HPLC: wild-type "FT" (B), pem1 (C), pem2 (D). Peak 1, chlorophyll b; peak 2, chlorophyll a, peak 1 ′ , chlorophyll b-like; peak 2 ′ , chlorophyll a-like. The abscissa axis is the elution time. molecular mechanism underlying the plant pale green trait in Brassica plants.
Chlorophyll biosynthesis in eukaryotic organisms consists of 16 enzyme-catalyzed reactions, with DVR being the most recently discovered enzyme in chlorophyll synthesis. DVR catalyzes the reduction of the 8-vinyl group to an ethyl group on the tetrapyrrole ring and is encoded by a single copy gene in model plants, A. thaliana and Oryza sativa. Nagata et al. (2005) first identified DVR (AT5G18660) as indispensable for the synthesis of MV-Chls as mutation in DVR resulted in an accumulation of DV-Chls. At the same time, mutant DVR caused a pale green phenotype with a lower Chls content and higher Chl a/b ratios. Nakanishi et al. (2005) discovered that single nucleotide substitutions in DVR led to chloroplasts lacking obvious grana stacks and accumulated DV-Chls with only a small account of MV-Chls. Nine nucleotide deletions in OsDVR resulted in the mutant presenting a yellow-green phenotype, with reduced total Chl levels, sole accumulation of DV-Chls, arrested chloroplast development, and decreased output (Wang et al., 2010a). The expressive level of BrDVR in mutants was not significantly different from that of the wild-type "FT, " and the promoter sequences showed no variation; however, DVR activity was decreased by approximately 50%. The marked decline in DVR activity inhibits chlorophyll synthesis and lowers chlorophyll content. Similarly, mutants of BrDVR lead to a partial loss of protein function, with an accumulation of mainly DV-Chls and a small amount of MV-Chls, with high Chl a/b ratios (Nakanishi et al., 2005). The higher Chl a/b ratios caused by the absence of functional BrDVR are attributed to the unique chemical structure (Nakanishi et al., 2005). The DVR catalyzes the conversion of the vinyl group at C8 of divinyl protochlorophyllide a 8 to the ethyl group, whereas CAO catalyzes the conversion of Chla to Chlb by converting the methyl group at C7 to the formyl group. Hence, the vinyl group side chain at C8 might interfere with the catalytic reaction of CAO at C7.
Single nucleotide substitutions in genes involved in plant biology may or may not result in large phenotypic changes. For example, Li et al. (2020) found that a single nucleotide substitution in sedoheptulose 1,7-bisphosphatase (SBPase) severely affects plant growth and grain yield in rice. A single nucleotide mutation (G664A) in the TraesCS7A01G480700.1 (chlI of Mg-chelatase) resulted in the mutant protein losing the ability to interact with the normal protein TaCHLI-7A . In A. thaliana, a single nucleotide substitution in DVR leading to one amino acid substitution in the protein did a partial loss of DVR function. However, another single nucleotide substitution in DVR resulted in early termination of translation in pbc2, which eventually led to a complete loss of protein function. In the present study, mutant BrDVR in pem1 and pem2 was caused by single nucleotide substitutions (G to A), and our results show that there is a reduction in DVR activity by approximately 50% with a partial loss of DVR function. This was observed by the reduction in chlorophyll content, accumulation of DV-Chls, and higher Chl a/b ratios in the mutants compared with those in the wild type, and these differences in the pigment contributed to the pale green phenotype. Phylogenetic tree and sequence alignment illustrated that DVR is conserved in a broad range of species. The amino acid variations at 357th of pem1 and 352th of pem2 were located in a highly conserved region in those species (Supplementary Figure S5). Therefore, we speculated that this region is associated with protein effective function.
The allelic mutant can be used to verify the function of the mutant gene as shown in lettuce (Lactuca sativa) by Huo et al. (2016); sorghum (Sorghum bicolor) by Jiao et al. (2017) and Chinese cabbage (Brassica rapa ssp. pekinensis) by Fu et al. (2020) and so on. Fu et al. (2020) mapped two allelic early bolting mutants to BrSDG8 in Chinese cabbage whereas Jiang et al. (2015) used two allelic mutants, nal1-2 and nal1-3, to identify that NAL1 plays a role in leaf growth and development. Zhang et al. (2018) effectively identified a round leaf shape gene in cucumber in two round-leaf (rl) allelic mutants. We selected seven similar pale green mutants from 286 yellow mutants in our library to analyze their allelism. It was found that the mutant genes of pem1 and pem2 were allelic. Mapping and sequence analysis found that a single non-synonymous nucleotide substitution (G to A) at position 1070 of pem1 and 1054 of pem2 in BrDVR. BrDVR activity assay and chlorophyll composition analysis demonstrated that impaired BrDVR led to partial loss of function, causing a pale green phenotype. Thus, we considered that BrDVR was responsible for the pale green phenotype of pem1 and pem2.
Our experiment revealed that BrDVR is the caused gene for the pale green phenotype of pem1 and pem2, and partial loss of its protein function is due to the single nucleotide substitution in BrDVR. The two allelic pale green phenotype mutants and the present study provide a good foundation for further investigation into DVR functions.

DATA AVAILABILITY STATEMENT
The raw sequence data of the BSR-Seq has been uploaded to the Sequence Read Archive at NCBI under the accession number SRR14339734. The raw sequence data of the whole genome re-sequencing has been uploaded to the Sequence Read Archive at NCBI under the accession number SRR14340669 and SRR14340670.

AUTHOR CONTRIBUTIONS
HF and SH conceived and designed the study. YZhao and MZ performed the experiments. YZhao wrote, reviewed, and edited the manuscript. YZhao and YZhan analyzed the data. All authors read and approved the manuscript. Supplementary Table S1 | Localization of chromosomal loci related to the pale green trait.