The mutation of ent-kaurene synthase, a key enzyme involved in gibberellin biosynthesis, confers a non-heading phenotype to Chinese cabbage (Brassica rapa L. ssp. pekinensis)

The presence of a leafy head is a vital agronomic trait that facilitates the evaluation of the yield and quality of Chinese cabbage. A non-heading mutant (nhm1) was identified in an ethyl methanesulfonate mutagenesis population of the heading Chinese cabbage double haploid line FT. Segregation analysis revealed that a single recessive gene, Brnhm1, controlled the mutant phenotype. Using MutMap, Kompetitive allele-specific PCR, and cloning analyses, we demonstrated that BraA07g042410.3C, which encodes an ent-kaurene synthase involved in the gibberellin biosynthesis pathway, is the nhm1 mutant candidate gene. A single-nucleotide mutation (C to T) in the fourth exon of BraA07g042410.3C caused an amino acid substitution from histidine to tyrosine. Compared to that of the wild-type FT, BraA07g042410.3C in the leaves of the nhm1 mutant had lower levels of expression. In addition, gibberellin contents were lower in the mutant than in the wild type, and the mutant plant phenotype could be restored to that of the wild type after exogenous GA3 treatment. These results indicate that BraA07g042410.3C caused the non-heading mutation. This is the first study to demonstrate a relationship between gibberellin content in the leaves and leafy head formation in Chinese cabbage. These findings facilitate the understanding of the mechanisms underlying leafy head development in Chinese cabbage.


Introduction
Chinese cabbage (Brassica rapa L. ssp. pekinensis) is a widely cultivated and economically important vegetable species in China and Southeast Asia. The presence of a leafy head is an essential morphological trait for evaluating the yield and quality of Chinese cabbage. Leafy heads of Chinese cabbage are rich in vitamins, dietary fiber, and nutrients 1,2 .
The formation of leafy heads in Chinese cabbage is biologically complex and usually spans four stages of development: the seedling stage, rosette stage, folding stage, and heading stage 3 (Fig. 1). Numerous factors, such as auxin concentrations, light intensity, temperature, and the carbohydrate-to-nitrogen ratio, may influence this process 4 . Several genes associated with the formation of Chinese cabbage leafy heads have been reported in previous studies. For instance, indoleacetic acid (IAA) induces transcriptional expression of the BcpLH (LEAFY HEADS) gene, indirectly promoting leafy head formation 1 . Mao et al. 5 suggested that the expression levels of BrpTCP (TEOSINTE BRANCHED 1, cycloidea, and PCF transcription factor) may affect the shape of leafy heads. Wang et al. 6 revealed that BrpSPL9 (SQUAMOSA PRO-MOTER BINDING-LIKE 9) genes could shorten the seedling and rosette stages to make the formation of leafy heads occur earlier. Moreover, Yu et al. 7 revealed the BrAN3 (ANGUSTIFOLIA 3) gene controls the formation of leafy heads. Other previous studies have also shown that multiple genes control leafy head formation 8,9 . Therefore, further investigation of the genes associated with leafy head formation could enhance our understanding of leafy head development.
Auxins are the most important hormones that affect the formation of leafy head in Chinese cabbage. Li 10 suggested that changes in the concentrations of auxin in the adaxial and abaxial surfaces of leaves cause the leaves to curl inwards to form a leafy head. Similarly, He et al. 3 suggested that auxin concentrations and distribution could influence the formation of leafy heads. Using a genome convergence approach to analyze the genomes of B. rapa and B. oleracea, Cheng et al. 2 found that various processes, including auxin-mediated signaling, may result in similar leafy head formation in Chinese cabbage and cabbage. Although some studies have shown that the synthesis, transport, and signal transduction of phytohormones, especially auxins, regulate the formation of Chinese cabbage leafy heads, the underlying molecular mechanisms and the genetic basis of leafy head development remain unclear.
In the present study, we identified a non-heading mutant, nhm1, whose mutation is controlled by a single recessive gene, from a Chinese cabbage ethyl methanesulfonate (EMS) mutagenesis population. Based on Mut-Map and Kompetitive Allele-Specific PCR (KASP) analyses, we demonstrate that BraA07g042410.3C, which encodes an ent-kaurene synthase (KS) participating in gibberellin (GA) biosynthesis, is the nhm1 candidate gene. As a key enzyme, KS plays an important role in GA synthesis (Fig. 2), and its mutation should affect plant development. By determining endogenous GA contents and applying exogenous GA 3 treatment, we demonstrated that nhm1 was a GA-deficient mutant. In contrast to previous studies, our study is the first to report that GA contents in leaves are related to the formation of leafy heads in Chinese cabbage. These results provide novel insights that could facilitate our understanding of the mechanisms of leafy head formation in Chinese cabbage.

Phenotypic analysis and inheritance of nhm1
The leaves of the mutant exhibited geotropic growth throughout all stages of development, and leafy heads were not formed during the heading stage (Fig. 3). This is in contrast to the wild-type heading Chinese cabbage double haploid (DH) line FT, which exhibited leafy head formation. Based on analysis of scanning electron microscopy (SEM) images, there was no significant difference in the areas of the adaxial and abaxial epidermal cells between the wild-type FT and the nhm1 mutant at the rosette stage (Fig. 4).
The nhm1 phenotypic data are listed in Table 1. All the F 1 plants obtained from the reciprocal cross between the wild-type FT and nhm1 mutant showed the wild-type phenotype. The segregation ratio between the wild-type and mutant phenotypes in F 2 plants was close to 3:1 (χ 2 = 0.16). In addition, in the backcross (F 1 × nhm1 mutant) generation, the separation ratio was close to 1:1 (χ 2 = 0.10). Therefore, we concluded that the nhm1 mutant trait was controlled by a single recessive nuclear gene, which we named Brnhm1.

Identification of candidate genes using MutMap
Through genome resequencing of the wild-type FT, nhm1 mutant, and the DNA pool, we obtained 135,076,750, 187,153,604, and 170,387,334 high-quality reads, respectively. Subsequently, 96.6%, 97.97%, and 97.78% of the clean reads in the wild-type FT, nhm1, and When a single-nucleotide polymorphism (SNP) index of 0.86 was used as the threshold, we located a candidate region of 2.91 Mb (26,010,373-28,918,429) on chromosome A07 (Fig. 5a). We found two SNP mutations (SNP A07, 27,142,769; SNP A07, 28,395,204) that were localized to the exon and caused the nonsynonymous mutation (Table S1).

SNP genotyping through KASP
To confirm the causal SNP, we used F 2 individuals to conduct KASP genotyping analysis of two mutated SNPs. The genotypic assay showed that the SNP (A07, 28,395,204) of BraA07g042410.3C cosegregated with the mutant phenotype. The mutant-phenotype plants had a T:T genotype, and the wild-type-phenotype plants had C: T or C:C genotypes for this SNP (Table S2). However, a recombinant event was observed at another SNP (A07, 27,142,769); here, both the G:G and G:A genotypes were detected in the mutant-phenotype plants (Table S2). This further confirms that BraA07g042410.3C was the causal gene of the nhm1 mutant. BraA07g042410.3C is homologous to Arabidopsis AtKS1 (At1g79460) and encodes the key enzyme KS in the biosynthesis pathway of GA.

Candidate gene cloning and sequence analysis
The MutMap and KASP analyses supported BraA07g042410.3C being the candidate gene for nhm1. Gene annotations showed that BraA07g042410 was 3752 bp in length and consisted of 13 exons. We cloned the cDNA sequences of BraA07g042410.3C from wild type FT and the nhm1 mutant. There was only a single base substitution at position A07 (28,395,204; C to T) in nhm1, which caused an amino acid to change from histidine (H) to tyrosine (Y) (Fig. 5b, c).  We performed quantitative real-time PCR (qRT-PCR) to measure the expression level of BraA07g042410.3C in the leaves of the wild type FT and nhm1 mutant. BraA07g042410.3C expression levels were lower in the nhm1 mutant than in wild type FT in all four stages of leaf development (Fig. 6).
The data are the means of three replicates (±SDs). *, ** Significant differences in expression levels at P < 0.05 and P < 0.01, respectively (Student's t-test).
nhm1 mutant plants had lower levels of GA and could be rescued by GA 3 treatment Bioactive GA 1 was detected only in the wild-type plants, and the level of bioactive GA 4 was significantly lower in the nhm1 mutant than in the wild type. However, bioactive GA 3 and GA 7 were not detected. In addition, the GA content (GA 9 , GA 12 , GA 15 , GA 19 , GA 20 , GA 24 , GA 29 , GA 34 , and GA 51 ) was significantly lower in the nhm1 mutant than in the wild type (Fig. 7). Based on the differences in GA contents between the wild-type and mutant plants, we investigated the effects of exogenous GA 3 application on the leaves of nhm1 mutant plants. As illustrated in Fig. 8, the phenotype of the mutant plants treated with exogenous GA 3 could be restored to wildtype phenotype.

Discussion
In the present study, a Chinese cabbage EMS-induced non-heading mutant, nhm1, was identified. First, phenotypic   identification revealed that the leaves of the mutant appeared to exhibit geotropic growth during the entire growth period, resulting in non-heading. Second, following MutMap and KASP analyses, a nonsynonymous mutation was detected in the fourth exon of BraA07g042410.3C. This gene encodes a protein with KS activity and thereby catalyzes the second step in the cyclization of GGPP to ent-kaurene in the GA biosynthesis pathway; we found that this gene harbors an SNP on A07 (at 28,395,204). Third, cloning and sequence analyses revealed that this mutant SNP involved only a single base substitution (C to T), which caused an amino acid to change from H to Y. Fourth, the BraA07g042410.3C expression levels were lower in the leaves of the mutant than in the leaves of the wild-type plants. Consistent with the qRT-PCR results, the GA contents in the mutant were significantly lower than those in the wild type. Last, the mutant phenotype could be rescued after exogenous GA 3 treatment.
The results indicate that the mutation in BraA07g042410.3C, which encodes a key GA biosynthesis enzyme, was responsible for the mutant phenotype. GAs, which are class of diterpenoid plant hormones, play key roles in plant life cycles, such as in the promotion of cell elongation and division, hypocotyl and stem elongation, seed germination, root growth, and flowering 11,12 . To ensure normal plant growth and development, active GAs must be produced and must accumulate. However, the biosynthesis of active GAs requires various functional enzymes to catalyze different intermediates, which involves a complex and multistep process 13 . In the model plant species Arabidopsis thaliana and rice, the main genes that encode key enzymes in each step of the GA biosynthetic pathway have been identified [14][15][16] . In the GA biosynthetic pathway, copalyl pyrophosphate synthase catalyzes the cyclization of GGDP to CDP at an early step, whereas Ent-kaurenoic acid oxidase in the endoplasmic reticulum catalyzes the conversion of ent-kaurene acid to GA 12 12,17,18 . KS catalyzes the second step in the cyclization of GGPP into ent-kaurene in the GA biosynthetic pathway. Given that KS the second key enzyme in the GA synthesis pathway, mutations in the gene encoding this protein considerably affect plant development. At present, the gene encoding KS has been cloned in pumpkin, Arabidopsis thaliana, stevia, lettuce and other herbs [19][20][21][22] . The GA2 mutant of Arabidopsis thaliana exhibits a lack of KS catalytic function, and germination cannot occur  24 , and GA 9 precursors and is also converted into GA 1 through the precursors GA 53 , GA 53 , GA 44 , GA 19 , and GA 20 in the 13-hydroxylation pathway 26 . GA 1 , GA 3 , GA 4 and GA 7 are thought to function as bioactive forms in plants. In the present study, we measured the endogenous GA levels in the leaves of mutants and the wild type. As illustrated in Fig. 7, the level of bioactive GA 4 in the mutant was significantly lower than that in the wild type. We speculate that this could be because the endogenous GA 9 , GA 12 , GA 15 , and GA 24 levels were significantly lower in the mutant than in the wild type. The absence of bioactive GA 1 in the nhm1 mutant may have been due to the significantly lower endogenous GA 12 , GA 19 , GA 20 , GA 44 , and GA 53 levels in the mutant than in the wild type. In addition, bioactive GA 3 and GA 7 were not detected in either the mutant or the wild type. These results indicate that the KS mutation may affect GA biosynthesis and that the non-heading phenotype may be caused by a lack of bioactive GAs.
Previous studies have indicated that GA is one of the most important hormones in regulating plant height. Mutations in GA3ox resulted in dwarf-type phenotypes of maize 27 and rice 28 . A mutation in the sd1 allele, which encodes GA20ox, in the IR8 rice variety reduced endogenous GA levels and led to a short stature 29,30 . In oilseed rape, the semidwarf mutant gene ds-3, which encodes a DELLA protein that negatively regulates plant height via the GA signaling pathway, was identified 31 . Li et al. 32 suggested that GmDW1 encodes KS and was responsible for a dwarf phenotype in soybean; the phenotypes of mutant plants that were treated with exogenous GA 3 were restored to those of the wild type. In the present study, a non-heading nhm1 mutant was identified. All leaves of the nhm1 mutant exhibited geotropic growth during the entire growth period, and leafy heads did not form. The mutated gene-encoded KS is a key enzyme in GA synthesis. We have demonstrated that nhm1 is a GA-deficient mutant through the determination of endogenous GAs and the exogenous application of GA 3 . In addition, we have demonstrated that the mutation of KS is responsible for the non-heading phenotype of Chinese cabbage.
Many previous studies have suggested that leafy head formation is linked to auxins in Chinese cabbage. He et al. 3 reported that auxins participate in the regulation of leafy head formation. Gao et al. 33 conducted a whole-genome annotation and a bioinformatics analysis of the auxin transport genes BrAUX/LAX (AUXIN/LIKE AUXIN RESISTANT), BrPIN (PIN-FORMED), and BrPGP (P-GLYCOPROTEIN) in Chinese cabbage and found that they play important roles in leafy head development. The authors also observed that polar auxin transport and the uneven distribution of auxins in leaves influenced the formation of leafy heads. Li et al. 34 identified a nonheading mutant (fg-1) of Chinese cabbage and showed that auxin and ABA signaling pathways were involved in regulating the formation of the leafy head. Our study is the first to demonstrate that GA contents in the leaves are related to leafy head formation in Chinese cabbage. Although the formation of leafy heads in Chinese cabbage may involve interactions among multiple hormones, how such hormones interact and their influence on leafy head formation require further investigation.
In conclusion, a non-heading mutant of Chinese cabbage, nhm1, was obtained following EMS mutagenesis. We demonstrated that BraA07g042410.3C, which encodes KS, an enzyme involved in the GA biosynthesis pathway, conferred a non-heading phenotype to the plants. These findings could facilitate our understanding of the mechanisms underlying leafy head formation and provide a genetic resource for Chinese cabbage crop improvement studies and activities.

Plant materials and mutant screening
FT is a DH line that is derived from the Chinese cabbage variety Fukuda 50 and exhibits desirable characteristics such as heat resistance, early maturity, an overlapping heading type, a compact head and a short growth period (only 45 days). In total, 7800 FT seeds were germinated at 25°C for 12 h. The germinated seeds were subsequently immersed in 0.8% EMS and placed in a 50turn shaker for 12 h. The treated seeds were then washed thoroughly in flowing water for 12 h and sown into trays. The surviving plants (M 0 ) were self-pollinated, and the M 1 seeds were harvested. A total of 14 stably inherited non-heading mutant (nhm) plants were harvested by screening the M 1 generation and identification in the M 2 generation. We selected a non-heading mutant, named nhm1, as the study material.

Genetic analysis of nhm1
For genetic analysis, wild type FT (P 1 ) was crossed with the nhm1 mutant (P 2 ) to generate F 1 and F 2 populations. Phenotypic characterization was performed for each generation (P 1 , P 2 , F 1 , F 2 ), and the segregation ratios of the BC 1 and F 2 populations were analyzed using the chisquare (χ 2 ) test. All plants were germinated and sown in a greenhouse at Shenyang Agricultural University.

Scanning electron microscopy
To observe the sizes of the adaxial and abaxial epidermal cells, similar parts of the rosette leaves from the wild type FT and the nhm1 mutant were fixed in 2.5% glutaraldehyde solution at 4°C for 48 h, processed according to the methods of Lin et al. 35 , and then examined under SEM (Hitachi TM3030, Japan).

Identification of candidate genes using MutMap
The candidate gene was identified by the modified MutMap method 36 . Fifty mutant plants were selected from among the F 2 plants, and a DNA pool was constructed by pooling equal amounts of plant leaf tissue from each plant. DNA from the wild type FT, nhm1 mutant, and the DNA pool was extracted using a DNAsecure Plant Kit (Tiangen, Beijing, China) and was resequenced using a NovaSeq 6000 System sequencer (Illumina, San Diego, USA).
Raw reads, including reads with adapters, reads with more than 10% unknown bases, and low-quality reads, were filtered and removed. The clean reads were then aligned to the Brassica reference genome sequence 37 (brassicadb.org/brad/datasets/pub/Genomes/Brassica_rapa/V3.0/) using Burrows-Wheeler Aligner (BWA) 38 . We used SAMtools to transform the alignment files into SAM/BWA files and to perform SNP calling 39 . The detection of SNPs and insertion-deletions was conducted using the GATK software toolkit 40 , and the results were analyzed using ANNOVAR 41 . Circos 42 was used to plot variations in the genome. We determined the SNP indices for 10 Chinese cabbage chromosomes. The ΔSNP index across the chromosomes of the B. rapa genome was obtained using sliding-window analysis, with a 4 Mb window size and a 10 kb increment for each step. The process for SNP filtering was as follows: first, we filtered out SNPs with multiple genotypes; second, we filtered out parental heterozygous SNPs; third, we filtered out SNPs that were not typical EMS-mutagenized-type mutations (C to T or G to A); and fourth, we retained SNPs that were homozygous and those that differed between the parents.

SNP genotyping through KASP
KASP genotyping analysis was used to detect the cosegregation of each SNP and to identify nhm1 candidate genes. One hundred eighty-four F 2 individual plants, including 50 plants that exhibited a mutant phenotype and 134 plants that had the normal phenotype, were used for KASP genotyping. The experiment was conducted by the Vegetable Research Center of Agriculture and Forestry Academy in Beijing. The allele-specific primers used are listed in Table S3.
Cloning and sequencing of the candidate gene The coding sequences of candidate genes in wild type FT and the nhm1 mutant were amplified by full-length PCR primers (Table S4). We purified the PCR products with a Gel Extraction Kit (CWBIO, Beijing, China) and then introduced the purified products into a pGEM-T Easy Vector (Promega, USA). The products were sequenced at GENEWIZ (Suzhou, China). We analyzed the sequences using DNAMAN V6 (Lynnon BioSoft, Canada).

Candidate gene expression analysis
The expression patterns of the nhm1 candidate gene were analyzed via qRT-PCR. Total RNA from the cotyledon, first true leaf, third true leaf, sixth true leaf, rosette leaf, and heading leaf of both the wild-type FT plants and the nhm1 mutant plants was extracted using TRIzol reagent (Invitrogen, USA) following the manufacturer's instructions. First-strand cDNA was synthesized using FastQuant RT Super Mix (Tiangen). The cDNA was subsequently used as a template for qRT-PCR together with SYBR Green PCR Master Mix (TaKaRa, Dalian, China). Actin was used as an internal control, and we used the 2 −ΔΔCt method to calculate relative gene expression levels 43 . Each experiment was performed for three biological replicates, and Bio-Rad iCycler IQ5 Manager (Bio-Rad Laboratories, Inc., Hercules, CA) was used to analyze the data. The sequences of the primers used for qRT-PCR are listed in Table S5.

Endogenous GA determination and exogenous GA 3 treatment
Endogenous GA levels in the leaves from both the wild-type and the mutant plants were determined using liquid chromatography-tandem mass chromatography. Wild-type FT and nhm1 mutant plants were grown under the same conditions. The third true leaves were obtained, after which the leaves were frozen directly in liquid nitrogen, ground to fine powder, and extracted with 1.0 mL of 80% methanol (v/v) at 4°C for 12 h. The extraction and content determination of GAs (GA 1 , GA 3 , GA 4 , GA 7 , GA 8 , GA 9 , GA 12 , GA 13 , GA 15 , GA 19 , GA 20 , GA 23 , GA 24 , GA 29 , GA 34 , GA 44 , GA 51 and GA 53 ) were based on those of a method reported previously 44 , with minor modifications.
To explore the response of nhm1 to GA, seedlings were sprayed with exogenous GA 3 solution (500 mg/L) six times at one-day intervals once the cotyledons had fully expanded. The control plants were treated with similar volumes of distilled water. Each treatment had three replicates. GA 3 powder (Solarbio) was first dissolved in 5 mL of ethanol and then diluted to the final concentration with distilled water.