Combination of quantitative trait locus (QTL) mapping and transcriptome analysis reveals submerged germination QTLs and candidate genes controlling coleoptile length in rice

Abstract

Submerged germination ability is a key trait of rice (Oryza sativa L.) in the success of water direct seeding. Therefore, mining submerged germination-related loci and seeking high tolerant rice germplasms are important for developing rice direct seeding breeding strategies. Here, we treated 272 recombinant inbred lines (RILs, F₁₅) from a cross from Luohui 9 (indica) and RPY geng (japonica) with simulated submergence for 7 days. We successfully bred several superior homozygous RILs with high submergence resistance, namely, R53, R137, R148, R176, and R180, via indica and japonica hybridization strategies. Then, the QTL mapping of coleoptile length trait was performed by R/qtl and IciMapping softwares based on the high-density bin marker genetic map and two QTLs responsible for coleoptile length were detected in chromosomes 1 (qCL-1.1) and 3 (qCL-3.1). The haplotype results of RILs showed that the aggregation of elite alleles of these two QTLs in an individual was beneficial to improving the rice submergence tolerance. According to contrasting genotype of qCL-1.1 and qCL-3.1 with different coleoptile lengths, RIL07 (R07) and R180 were selected to further RNA-seq analysis for analyzing the possible molecular mechanisms of rice under submergence stress. Finally, LOC_Os01g04430, LOC_Os01g04530, and LOC_Os03g22720 for qCL-1.1 and qCL-3.1 were characterized to affect submergence stress by the overlapped analysis of two independent sets of RNA-seq data, qRT-PCR verification, and gene sequence alignments. In this study, the submerged germination-related QTLs/genes and high submergence-tolerant RILs provided new genetic resources for breeding rice varieties suitable for direct seeding via molecular breeding strategies in the immediate future.

INTRODUCTION

Among the various planting methods of rice, water direct seeding is gradually becoming an ideal method for high-efficiency cultivation of rice due to its advantages of simplicity, water-saving, and labor-saving (Jin et al., 2018; Yang, Sun, et al., 2019). However, the seeds of water direct seeding are exposed to submergence during germination and are subjected to hypoxia stress, which will cause no germination, seedling death, and poor crop seedling establishment (Ismail et al., 2009; Jeong et al., 2020). The rapid elongation of coleoptile in a short period is considered to be a witty “escape” strategy used by rice young seedlings to tolerate submergence during the germination phase, commonly known as submerged (anaerobic) germination (Hsu & Tung, 2017; Yang, Yang, et al., 2019). The faster the coleoptile grows, the sooner the rice seedlings can get rid of hypoxia stress, which will greatly increase the survival rate of rice under submergence conditions (Hsu & Tung, 2015; Kuya et al., 2019). However, not all rice varieties can germinate under submerged conditions and different rice accessions show different submerged germination abilities (Doley et al., 2018; Ghosal et al., 2020; Rauf et al., 2019; Zhang et al., 2017).

Under submerged conditions, coleoptile elongation is a typical quantitative trait, which varies greatly among different varieties, and has now become a typical trait for studying rice submerged germination and screening of rice germplasms with high tolerance to submergence (Ma et al., 2020; Magneschi et al., 2009). So far, several quantitative trait loci (QTLs) related to coleoptile elongation have been identified in different genetic or natural populations (Kuroha & Ashikari, 2020; Ma et al., 2020; Pucciariello, 2020). For example, Angaji et al. (2010) reported five QTLs associated with tolerance of submergence during germination based on SSR markers, namely, qAG-1-2, qAG-3-1, qAG-7-2, qAG-9-1, and qAG-9-2, respectively. However, only qAG-9-2 was cloned based on the QTL mapping of coleoptile elongation trait to date, encoding T6P phosphatase, and was involved in sugar signaling regulation and the linking of trehalose metabolism to starch mobilization (Kretzschmar et al., 2015). This QTL has been applied to crop improvement and has been successfully validated in the field (Chamara et al., 2018). Septiningsih et al. (2013) identified eight QTLs related to tolerance of submergence during germination, namely, qAG2, qAG5, qAG6, qAG7.1, qAG7.2, qAG7.3, qAG9, and qAG12. Benefitting from the development of second-generation sequencing technology coupled with the construction of high-density genetic bin maps, the AG2 (qAG7.1) interval was successfully reduced from 7 Mb to less than 0.7 Mb (Tnani et al., 2021). Recently, many submerged germination QTLs controlling coleoptile length were reported by high-density genetic map or genome-wide association study (GWAs) analysis (Hsu & Tung, 2015, 2017; Kuya et al., 2019; Yang, Sun, et al., 2019; Zhang et al., 2017). For instance, Yang, Sun, et al. (2019)) conducted QTL mapping analysis using a high-density bin map of a recombinant inbred line (RIL) population derived from a cross between YZX and 02428 and obtained 25 loci related to anaerobic germination tolerance.

Since different rice varieties have different elite submerged germination tolerance genes, the aggregation of elite alleles from different varieties can significantly improve the submerged germination tolerance of rice (Pucciariello, 2020; Yang, Yang, et al., 2019). Thus, it is necessary to apply multi-omics methods, the well-designed populations, and the high-density genetic map to locate new submerged germination QTLs facilitating precise rice breeding in the future. In this study, a high-generation RIL population containing 272 lines from Luohui 9 (indica) × RPY geng (joponica) was developed, sequenced, and genotyped of 4758 bin blocks. We simulated the submergence treatment of 272 RILs and their parents at the germination stage and analyzed the variation of coleoptile length trait in this population. Then, we used R/qtl and IciMapping softwares to QTL mapping the coleoptile length trait, respectively. Next, we analyzed the effect of QTLs in the regulation of coleoptile length using RILs. Finally, candidate genes in QTLs were fine mapped based on a comprehensive analysis of QTL mapping and RNA-seq.

MATERIALS AND METHODS Plant materials and submergence treatment

A total of 272 RILs (F₁₄), Luohui 9, and RPY geng were planted in Ezhou, Hubei, China, during summer, and the seeds were harvested after maturity for submergence simulation experiment in an artificial climate adjustment room. The harvested seeds (F₁₅) were treated at 40°C for 5 days to break the dormancy. Sixty filled and healthy seeds were selected for each line to be sterilized with 5% sodium hypochlorite and then soaked in distilled water for three times. Then, the seeds were placed in a conical flask with a water depth of 5 cm for 7 days submergence treatment with a daily photoperiod of 14-h light 30°C/10-h dark 28°C. All RILs were placed randomly with three biological replications. For each replication, the coleoptile of 20 randomly selected seeds was measured by Image J software (https://imagej.nih.gov/ij/download.html) (Figure S1).

In addition, the coleoptile samples of R07 and R180 after 7 days submergence treatment were collected and immediately frozen in liquid nitrogen and stored at −80°C for further RNA-seq.

Linkage map construction

Young leaves were collected from the seedlings of 272 RILs, and genomic DNA was extracted using the minor modified CTAB (Healey et al., 2014). Genomic DNA libraries of 272 RILs were sequenced by paired-end (2 × 150 bp) method on Illumina HiSeqTM 2500 (Illumina). The Nipponbare MSU v7.0 (http://rice.plantbiology.msu.edu/index.shtml) was the reference genome during genetic map construction. We carried out data filtering, SNP calling, and bin map construction of these 272 RIls according to the method described by Huang et al. (2009).

QTL analysis

The QTL mapping of the average coleoptile length of three biological replications was analyzed by R/qtl (Arends et al., 2010) and IciMapping (Meng et al., 2015). In R/qtl, the CIM interval mapping method was adopted, the threshold was set by 1000 times of PT test, and the LOD values of 0.1 and 0.01 significance levels were taken as the critical value of QTL significance. The confidence interval was calculated with the function “lodint” (Dupuis & Siegmund, 1999), and the drop value was set to 1.5. In IciMapping, QTLs were identified by inclusive composite interval mapping combined with additive mapping (ICIM-ADD) method, LOD value was set to 3, and p-value was set to 0.05.

RNA-seq analysis

Total RNA was isolated from coleoptile samples from R07 and R180 using CTAB method, and cDNA library was prepared NEBNext^® Ultra™ RNA Library Prep Kit for Illumina^® (New England Biolabs) and sequenced on an Illumina sequencing platform (HiSeq 2500). All raw data were filtered by fastp (Chen et al., 2018) and mapped to the Nipponbare (NIP) genome (MSU v7.0, http://rice.plantbiology.msu.edu/index.shtml) using hisat2 (Kim et al., 2015). Mapped reads were counted by featureCounts (Liao et al., 2014), and differentially expressed genes (DEGs) were identified by DEseq2 with |log₂^{fold change}| ≥ 2 and a false discovery rate (FDR) < 0.01 (Kong et al., 2020).

In addition, we downloaded 18 RNA-seq data from NCBI (accession number: SRP193577) containing 2, 3, and 4 days sequencing data of two rice (YZX, indica, and 02428, japonica) with opposite submergence tolerance using the same data processing flow (Yang, Sun, et al., 2019). DEGs in QTLs were also identified by DEseq2 with |log₂^{fold change}| ≥ 2 and a false discovery rate (FDR) < 0.01 (Kong et al., 2020).

qRT-PCR verification of RNA-seq data

To confirm the expression of the submergence-responsive genes in coleoptile samples, qRT-PCR was performed (Primers in Table S1). The qRT-PCR reaction (10 μl) was formulated using the 2 X SYBR Green qPCR Master Mix (US Everbright^® Inc., Suzhou, China). All qRT-PCRs were carried out on a CFX96 Touch™ Real-Time PCR Detection System (Bio-Rad). The UBI gene was used as internal controls for qRT-PCR analysis. The gene expression fold change was calculated by the 2^−ΔΔCT method from three biological replicates.

Sequence analysis of candidate genes for qCL-1.1 and qCL-3.1

To characterize the sequence variation of candidate genes between parents, we obtained the coding sequences (CDSs) of candidate genes from the high-quality chromosomal level genomes of parents (RPY geng and Luohui 9, unpublished results) and NIP. These CDSs were verified by the sequencing results of PCR products (Primers in Table S1). Sequence alignments of candidate genes were performed with DNAman (https://dnaman.software.informer.com/).

RESULTS Significant differences of coleoptile length in RILs

The parents of RILs showed significant differences in submergence tolerance, and Luohui 9 had a longer coleoptile than RPY geng (33.5 mm vs. 22.5 mm) (Figure 1a). 272 RILs showed differentiation of the average coleoptile length from 21.37 mm to 42.57 mm (Figure 1b). We noticed that the average coleoptile length of 271/272 RILs exceeds RPY geng and 113/272 RILs exceeds Luohui9, and the coleoptile length of 5 RILs exceeds 40 mm, namely, R53, R137, R148, R176, and R180. These results indicated that the hybridization of indica X japonica rice can produce rice varieties with stronger submergence tolerance.

FIGURE 1. Coleoptile length of RPY geng, Luohui 9, and 272 recombinant inbred lines (RILs)

QTL detection of coleoptile length

In this study, we obtained 1,339,300 high-quality SNPs based on 3× resequencing of RILs and divided these SNPs into 4758 bin markers among 12 linkage groups with a total map distance of 2356.41 cM according to linkage relationship (Figure S2). The QTL mapping of R/qtl identified two coleoptile length-related QTLs distributed on chromosome 1 and chromosome 3 (Figure 2). The QTL mapping of IciMapping completely overlapped R/qtl results, which verified the accuracy of our QTL results (Table 1). We took the overlapped regions of the two software as the confidence intervals of qCL1.1 (2002462 bp−2026897 bp) and qCL-3.1 (12604120 bp −12940183 bp). qCL1.1 and qCL-3.1 had 7.14% and 8.81% phenotypic variation explained rates, respectively, and opposite additive effects (Table 1 and Figure 2). These results proved that coleoptile length is a complex quantitative trait controlled by multiple minor genes.

FIGURE 2. QTLs mapping results in the high-density bin map of the Luohui 9 X RPY geng. Note: The black line and the dashed line represented the LOD value of 0.1 and 0.01 significance levels, respectively

TABLE 1 QTL mapping results of coleoptile length

[Table omitted. See PDF]

Confirmation of the qCL1.1 and qCL-3.1 for coleoptile length

To further confirm the function of the newly identified genetic locus of qCL1.1 and qCL-3.1, all RILs were divided into RPY geng genotyping (AA) and Luohui 9 (BB) RILs based on our bin genetic map genotyping results. qCL-1.1 (BB) RILs had longer coleoptiles than qCL-1.1 (AA) RILs (Figure 3a), while qCL-3.1 (AA) RILs had longer coleoptiles than qCL-1.1 (BB) RILs (Figure 3b). It was worth to mention that the synergistic genotypes of these two QTLs were completely opposite. As expected, the RILs of qCL1.1 and qCL-3.1 synergistic genotype aggregation showed the longest coleoptile among four different combinations (Figure 3C). For example, R180 (qCL-1.1^{Luohui 9} qCL-3.1^{RPY geng}, black point in Figure 3c) had the longest coleoptile (42.57 mm), while R07 (qCL-1.1^{RPY geng} qCL-3.1^{Luohui 9}, black triangle in Figure 3c) had only 25.52 mm coleoptile. These results indicated that the identified QTLs had stable effects on coleoptile growth under submergence stress.

FIGURE 3. Confirmation of the qCL1.1 and qCL-3.1 for coleoptile length. (a) qCL-1.1 (BB) RILs had longer coleoptiles than qCL-1.1 (AA) RILs. (b) qCL-3.1 (AA) RILs had longer coleoptiles than qCL-1.1 (BB) RILs. (c) The longest coleoptile was found in the four genotypes. Note: AA stands for RPY geng genotype, and BB stands for Luohui 9 genotype

RNA-seq analysis and DEGs annotation

To better understand the transcriptomic response to submergence stress, we performed RNA-seq analysis of coleoptile samples from R180 (qCL-1.1^{Luohui 9} qCL-3.1^{RPY geng}, tolerant) and R07 (qCL-1.1^{RPY geng} qCL-3.1^{Luohui 9}, sensitive) after submergence stress for 7 days. A total of 19.99 and 18.38 Gb clean data of R07 and R180 were obtained, and the Q30 of each sample exceeded 94% (Table S2). The mapping results showed that 93.92%–95.37% of clean reads in all 6 samples can be mapped to the MSU 7.0 genome (Table S2). Pearson correlation results showed that there was an excellent correlation between three biological duplications, while the correlation between R07 and R180 was less than 0.5, indicating that R180 and R07 have different reaction mechanisms under submergence stress (Figure 4a). Then, 1226 up-regulated DEGs and 1027 down-regulated DEGs were identified in R180 versus R07 (Figure 4b and Table S3). GO annotation (corrected p-value < 0.05) revealed that DEGs involved multiple stress-related biological processes, including defense response, response to stress, response to stimulus, lipid localization, and lipid transport (Table S4). Heat map of DEGs clustered all DEGs into two groups, namely high expression group in R180 (up-regulated DEGs) and high expression group in R07 (down-regulated DEGs) (Figure 4c). We next performed KEGG enrichment (corrected p-value < 0.05) of up-regulated and down-regulated DEGs, respectively, and the results showed that metabolism of other amino acids, starch and sucrose metabolism, glutathione metabolism, metabolism, and metabolism of terpenoids and polyketides were enriched by both up-regulated and down-regulated DEGs. Down-regulated DEGs were unique to plant–pathogen interaction, cytochrome P450 pathways, and diterpenoid biosynthesis, while up-regulated DEGs involved in 15 pathways, namely amino acid metabolism, biosynthesis of other secondary metabolites, carbohydrate metabolism, carbon fixation in photosynthetic organisms, cyanoamino acid metabolism, energy metabolism, environmental adaptation, etc. (Figure 4d).

FIGURE 4. RNA-seq results of R180 versus R07. (a) Pearson correlation heatmap of all six samples. (b) Volcano Plot of all differentially expressed genes (DEGs), black and red dots mean non-DEGs and DEGs, respectively. (c) Heat map of DEGs. (d) KEGG enrichments of Up/down-regulated DEGs

Candidate genes prediction

In this study, we found that qCL-1.1 and qCL-3.1, respectively, contained 26 genes and 184 genes involved in multiple biological processes, which is difficult to find reliable candidate genes (Table S5). To narrow down the candidate genes, a total of three and three DEGs were identified in qCL-1.1 and qCL-3.1 (Figure 5) after filtering low-expressed genes (FPKM < 2). Of them, LOC_Os01g04430 and LOC_Os01g04530 showed higher expression level (>five-fold changes) in R07 than R180, while four remaining DEGs showed higher expression level (>five-fold changes) in R180 than R07 (Figure 5a). A BlastP homology search of these six DEGs sequences to Arabidopsis thaliana revealed that LOC_Os01g04430 and LOC_Os01g04530 were the homologous genes of AtLRK10L1.2 and LOC_Os03g22210 and LOC_Os03g22720 were homologous genes of AT3G16660 coding Pollen Ole e 1 allergen and extension family protein precursor. Interestingly, LOC_Os01g04500 and LOC_Os03g22720 that only expressed in R180 did not find homologous genes in A. thaliana. In this study, LOC_Os01g04430, LOC_Os01g04530, LOC_Os03g22210, and LOC_Os03g22720 were verified by RT-PCR and the results were completely consistent with RNA-Seq (Figure 5b).

FIGURE 5. Expression profile of DEGs in in qCL-1.1 and qCL-3.1. (a) The heatmap of DEGs' expression level between R07 (sensitive) and R180 (tolerant). (b) RT-PCR results of four selected DEGs. (c) The heatmap of DEGs’ fold change in qCL-1.1 and qCL-3.1 between YZX (sensitive) and 02428 (tolerant). Note: The red genes represented overlapped DEGs of two independent RNA-seq data sets

Interestingly, we also found LOC_Os01g04430, LOC_Os01g04530, and LOC_Os03g22720 as DEGs in the published RNA-seq of YZX (sensitive) and 02428 (tolerant) (Figure 5c, Figure S3). The expression profiles of these three DEGs from public data were completely consistent with our RNA-seq result that LOC_Os01g04430 and LOC_Os01g04530 had higher expression level in sensitive varieties, and LOC_Os03g22720 has higher expression level in tolerant varieties, which suggested that LOC_Os01g04430, LOC_Os01g04530, and LOC_Os03g22720 were functionally conservative under submergence stress and they were accepted as the most likely candidate genes for qCL-1.1 and qCL-3.1.

Next, we aligned the nucleotide sequences of candidate genes for qCL-1.1 and qCL-3.1 from NIP, RPY geng, and Luohui 9. There were nine synonymous mutations, seven non-synonymous mutations, one insertion mutation (at amino acid sequence position 13), and one complex mutation of LOC_Os01g04430 in NIP and RPY geng vs Luohui 9, and Luohui 9 had a longer sequence of LOC_Os01g04430 gene compared with NIP and RPY geng (Figure 6a). Similarly, we identified 11 synonymous mutations, 26 non-synonymous mutations, one 3-bp deletion mutation (at amino acid sequence position 29), and one complex mutation of LOC_Os01g04530 in sequence alignment of NIP and RPY geng vs Luohui 9 (Figure 6b). Notably, LOC_Os03g22720 was not found in Luohui 9, which explained why LOC_Os03g22720 was not expressed (FPKM = 0) in R07 (qCL-1.1^{RPY geng} qCL-3.1^{Luohui 9}) (Figure 6c). The major differences in the structure of these genes further supported LOC_Os01g04430, LOC_Os01g04530, and LOC_Os03g22720 as candidate genes for qCL-1.1 and qCL-3.1.

FIGURE 6. Gene sequence alignment of candidate genes for qCL-1.1 and qCL-3.1 between NIP, RPY geng, and Luohui 9. Note: The red arrows mean non-synonymous, insertion, deletion, and complex mutations

DISCUSSION

Coleoptile elongation was considered as the strategy of rice escaping from submergence stress (Hsu & Tung, 2015) and was a very effective trait for screening submergence-tolerant rice varieties and identifying QTLs related to submergence tolerance (Gao et al., 2020; Kuroha & Ashikari, 2020; Kuya et al., 2019). In this study, we evaluated the submergence tolerance of 272 RILs and their parents by using the coleoptile length as an important tolerant trait. Luohui 9 with longer coleoptiles was more resistant to submergence than RPY geng. We speculated that this difference may be related to the stored amylose content between indica and japonica rice. Xu et al. (2013) reported that the amylose content of indica rice was higher than that of japonica rice varieties. On the contrary, our experiments demonstrated the feasibility of indica-japonica hybridization to breed stronger submergence-tolerant varieties and successfully selected several high submergence-tolerant RILs. The excellent performance of these high submergence-tolerant RILs can be derived from the synergistic genotype aggregation of the parental QTLs.

Based on the high-density genetic map of the O. nivara introgression lines (IL) population, Liu et al. (2021) identified one coleoptile length-related QTLs (qCLN1) in chromosome 1 under anaerobic conditions. The peak marker (bin150) of qCLN1 was located between 38,174,447 and 38,387,308 bp harboring the “green evolution” gene semi-dwarf1 (sd1). Kuroh et al. (2018) also reported sd1 responsible for submergence-induced internode elongation via GWAs and the SD1 protein direct increased synthesis of gibberellins, which promoted internode elongation. Su et al. (2021) found several QTNs in chromosomes 1 and 3 using GWAs of 209 natural rice populations. Hsu and Tung (2015), Nghi et al. (2019), and Nishimura et al. (2020) also reported submerged germination QTLs controlling coleoptile length in chromosomes 1 and 3. The distribution of these previous QTLs in different segments of chromosomes 1 and 3 suggested that these two chromosomes carry multiple submerged germination-related QTLs. The physical location comparison of our identified QTLs with the previous QTLs confirmed that we identified two new submerged germination QTLs controlling coleoptile length. The aggregation of elite alleles can significantly improve rice submergence tolerance, suggesting that these two QTLs are hopeful to further improve rice submergence tolerance through aggregation of submergence-tolerance genes that have been widely used, such as Sub1A and qAG-9-2.

The regulation network of submergence tolerance is very complicated and involves many factors, such as carbohydrate metabolism, fermentation, reactive oxygen species, lipid peroxidation, hormone induction, cell division and expansion, Ph reduction, and ferrous ions (Liu et al., 2021; Yang, Sun, et al., 2019). Of them, insufficient energy supply under low-O₂ conditions caused by submergence is a major bottleneck for seed germination and seedling survival (Miro & Ismail, 2013). The better submergence tolerance of R180 compared with T07 may be attributed to the more efficient energy supply under hypoxic conditions, because Up-DEGs (highly expressed in R180) were significantly enriched in energy metabolism (Figure 4d).

Earlier studies demonstrated that DEGs analysis within QTLs can effectively locate candidate genes (Kong et al., 2020; Yang, Sun, et al., 2019; Zhang et al., 2017). RNA-seq analysis showed that LOC_Os01g04430, LOC_Os01g04530, and LOC_Os03g22720 were DEGs between sensitive and tolerant varieties and had the same expression bias in public RNA-seq data. These results suggested that these three genes are important candidate genes for qCL-1.1 and qCL-3.1 and may be functionally conserved in rice submergence tolerance among different rice varieties. The potential role of LOC_Os01g04430 and LOC_Os01g04530 in abiotic stresses was supported by homologous genes in other species. The AtLRK10L1.2, an Arabidopsis ortholog of LOC_Os01g04430 and LOC_Os01g04530, was GO annotated as cellular response to abscisic acid stimulus, cellular response to water deprivation, photoperiodism, and flowering. Lim et al. (2015) reported that the AtLRK10L1.2 is involved in ABA-mediated signaling and drought resistance. Notably, LOC_ Os01g04430 and LOC_ Os01g04530, and LOC_Os03g22720 showed completely opposite expression bias between sensitive and tolerant varieties. LOC_ Os01g04430 and LOC_ Os01g04530 showed higher expression in sensitive varieties than in tolerant varieties (Figure 6a), suggesting that they may play a negative role in rice submergence tolerance, and further experiment was required for validating this speculation.

CONCLUSION

In this report, we developed a high-generation RILs population of indica and japonica rice and constructed a high-density genetic map of this population. We found two novel submerged germination QTLs, qCL-1.1 and qCL-3.1, co-localized from two linkage mapping software. The integration of QTL mapping and RNA-seq analysis assisted us to narrow down the region within the identified QTLs, and three candidate genes were accepted as the most likely candidate genes, which provided highly reliable targets for further gene cloning and breeding submergence-tolerant rice cultivars suitable for direct seeding. On the contrary, this multi-omics comprehensive analysis provided new insights into the genetic basis controlling submergence tolerance in rice, which laid a foundation in rice improvement.

ACKNOWLEDGMENT

We thank Professor Qian Qian (State Key Laboratory of Rice Biology, China National Rice Research Institute, Hangzhou, Zhejiang, 310006, People’s Republic of China) for proving us RPY geng materials.

CONFLICT OF INTEREST

The authors declare that they have no competing interests for this research.

AUTHOR CONTRIBUTIONS

W.K. and Y.L. conceived and designed the experiments. W.K. performed all of the experiments, analyzed the data, prepared the figures and tables, and wrote the paper. S.L. performed all of the experiments and completed the preliminary analysis of the data. C.Z. and Y.Q. performed parts of the experiments, figures, and tables. All authors read and approved the final version of the manuscript.

DATA AVAILABILITY STATEMENT

The RNA-seq data are publicly available from the NCBI (BioProject ID PRJNA760818, https://dataview.ncbi.nlm.nih.gov/object/PRJNA760818?reviewer=mpd74fr5iqcfua99okpbo1l2ja).