Molecular basis of high temperature-induced bolting in lettuce revealed by multi-omics analysis

High temperature induces early bolting in lettuce (Lactuca sativa L.), which affects both quality and production. However, the molecular mechanism underlying high temperature-induced bolting is still limited. We performed systematical analysis of morphology, transcriptome, miRNAs and methylome in lettuce under high temperature treatment. Through a comparison of RNA-Seq data between the control and the high temperature treated lettuces at different time points totally identified 2944 up-regulated genes and 2203 down-regulated genes, which cover three floral pathways including photoperiod, age and gibberellin (GA) pathways. Genome wide analysis of miRNAs and methylome during high temperature treatment indicated miRNAs and DNA methylation might play a role controlling gene expression during high temperature-induced bolting. miRNA targets included some protein kinase family proteins, which potentially play crucial roles in this process. Together, our results propose a possible regulation network involved in high temperature-induced bolting.


Introduction
Plants undergo different developmental stages in their life cycle [1]. Among them, the transition from vegetative to reproductive growth is a crucial step. In this process, the shoot apical meristem (SAM) plays a vital role. During vegetative stage, the SAM produces leaves from the peripheral zone. Upon translation, the SAM elongates and converts into the inflorescence meristem (IM), generating flowers [2]. This floral transition is tightly controlled by developmental, hormonal and environmental signals. Among them, the photoperiod and temperature are the most important environmental cues [3]. In Arabidopsis, previous studies have shown that six major genetic pathways, including photoperiod, ambient temperature, age, gibberellin (GA), autonomous, and vernalization pathways, control this transition [4]. The flowering locus T (FT) and suppressor of overexpression of constans 1 (SOC1/AGL20) are among the most important floral pathway integrators that control the transition to flowering [5]. FT acts as a long-distance signal molecular, moving from the leaves to the apical meristem Open Access † Jinghong Hao and Junwei Yang contributed equally to this work. Full list of author information is available at the end of the article through the phloem and promoting flowering at the SAM [6,7]. Previous studies have shown that many transcription factors are involved in the transcription regulation of FT. The photoperiod-associated factors constans (CO) and phytochrome-interacting factor 4 (PIF4) activate FT expression in response to long-day (LD) conditions and increasing temperature [8,9], while vernalization-dependent gene flowering locus C (FLC) and APETALA 2 (AP2) -like genes repress FT expression [10,11]. The transcription level of another floral integrator gene SOC1 is mainly regulated by CO and FLC [12]. CO promotes SOC1 expression mainly through FT, while FLC represses SOC1 transcript by binding to its promoter. Subsequently, FT activates the downstream floral meristem identity genes AP1 and its paralog CAULIFLOWER (CAL) [13], and SOC1 activates LEAFY (LFY), leading to floral initiation [14].
Although the genetic network of flower transition is well-studied in Arabidopsis, this knowledge in most crops is still lacking. Lettuce (Lactuca sativa), belonging to Compositeae family, is an important leafy vegetable worldwide [15]. Cool temperature is suitable for its cultivation, and high temperature induces its bolting (rapid stem elongation), which is irreversible and seriously affects the quality of harvestable crops and flavor [16]. So, studying the mechanism of heat-induced bolting is important for lettuce breeding. Previous studies showed that knockdown the homolog of the Arabidopsis SOC1 or FT in lettuce led to a bolting delay with reduced LFY and AP1 expression and insensitivity to high temperature [17,18], indicating that FT, SOC1, and their downstream genes AP1 and LFY likely play important roles in lettuce floral transition. In addition, GA and related genes were reported to be associated with bolting in response to high temperatures [19][20][21]. These studies suggested that the function of these central genes in floral translation seems to be conserved between lettuce and Arabidopsis, but FT and SOC1 integrate which floral pathways to regulate high temperature-induced bolting in lettuce remain largely unknown.
To further investigate the relevant signaling pathway involved flowering transition in lettuce, RNA-Seq, miR-NAs and DNA methylome were performed in lettuce tips during high temperature treatment. RNA-Seq data showed that three floral pathways, including photoperiod, age and gibberellin (GA), are potential involved in high temperature-induced bolting. Moreover, epigenome analysis indicated that miRNAs and DNA methylation may be also involved in bolting by regulating floral pathways indirectly. Together, our study provided new molecular insights into the mechanism of high temperature-induced bolting in lettuce.

The effect of high temperature on lettuce bolting
To study high temperature-induced floral transition in lettuce, the seedlings with seven leaves were transferred to the plant incubators with 14 h light at 33 °C and 10 h darkness at 25 °C, while the control was cultured at 20 °C on the day and 13 °C at night. We observed that the high temperature treatment could induce a longer stem compared to the control group at the 8th day (C8) ( Fig. 1a, b). After 24 day treatment, the lettuce stem has elongated to approximately 7 cm (Fig. 1a). To further analyze the morphological changes of lettuce stem tips, histological sections were performed to observe the morphological features at different time points. As shown in Fig. 1, compared to control at 24th day (C24), the SAM of the lettuce treated by high temperature at 24th day (HT24) is becoming from flat to domed, marking the transition from vegetative to reproductive growth with elongated SAM (Fig. 1c). So, the lettuce at C0, C2, C8, C24, HT2, HT8 and HT24 were selected for further study.

Transcriptome reprogramming during high temperature-induced bolting
To explore which pathways were involved in high temperature induced bolting in lettuce, we performed RNA-Seq of the same control and high-temperature treated lettuces as in morphological analysis. Hiseq-PE150 sequencing produced about 6 Gb data for sample. The principal component analysis (PCA) revealed homogeneity between replicates at each time point (Fig. S1). We identified 2944 up-regulated differentially expressed genes (DEGs) and 2203 down-regulated DEGs totally by comparison between the control and the high temperature treated lettuces ( Fig. 2a- Gene Ontology (GO) enrichment analysis was performed at different time points to probe the pathway involved in high temperature-induced bolting. For those up-regulated DEGs under high temperature treatment, the significant enrichment pathways included 'transcription factors and regulators' and 'trans-membrane signaling receptor' (Fig. S2-S4). For those down-DEGs under high temperature treatment, the 'heat shock protein binding' was shown to be enriched significantly ( Fig. S2-S4). Moreover, the Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis showed that the plant hormone signal transduction is the significant enrichment pathway after high temperature treatment for 8 days (Fig. S2-S4). It is noteworthy that those pathways may be highly related to bolting in lettuce.

Dynamic expression of genes related to flower transition during high temperature-induced bolting
To study the DGE expression pattern, these DEGs that may be involved in the flowering pathway were selected for analysis. The result showed that the expression levels of most these DGEs were induced by high temperature treatment (Fig. S5). Previous studies reported that the floral pathway integrators SOC1 and FT play important roles in lettuce bolting [17,18], we first analyzed the expression levels of SOC1 and FT based on the RNA-Seq data. Unfortunately, the expression of FT (Lsat_1_v5_ gn_2_17881) was too low to be detected. The expression of SOC1 (Lsat_1_v5_gn_9_53801) was began to upregulated after high temperature treatment at 24th day (Fig. S6a), indicating that SOC1 may play an important role in lettuce bolting. These results prompted us to further investigate which floral pathway involved in flower transition regulation from the RNA-Seq data. Among these genes, most of them are COs and PIFs, which are involved in photoperiod pathways. Meanwhile, the gene GA 20-OXIDASE 1 (GA20OX1) responsible for GA biosynthesis was significantly up-regulated in high temperature-treated lettuce at 8th day (Fig. 2c). Moreover, some AP2-like genes in age pathway were down-regulated after high temperature treatment at 24th day (Fig. 2). Apart from them, vernalization is also a major determinant of flowering time, which leads to down regulation of FLC to promote flower [3]. To study whether vernalization pathway is involved in high temperature-induced flowering in lettuce, we analyzed the expression levels of vernalization pathway genes VIN3 and FLC after high temperature treatment. The expression level of VIN3 was not changed, which FLC was down-regulated after high temperature treatment (Fig. S6), indicating that vernalization pathway may be also required for heat induced flowering in lettuce. Taken together, these data suggested that multiple pathways may be involved in the regulation of high temperature-induced bolting in lettuce.

Dynamic changes in microRNAs during high temperature-induced bolting
MicroRNAs (miRNA) play important roles in gene expression regulation by targeting mRNAs [22]. To identify whether miRNAs potentially regulate high temperature-triggered bolting, small non-coding RNA sequencing was performed in different lettuce samples (i.e., C0, C2, C8, C24, HT2, HT8 and HT24). In total, 5295 unique miR-NAs were identified in lettuce. As expected, the expression levels of some miRNAs were increased gradually in lettuce with increasing high temperature treatment (Fig.  S7). To study the potential role of miRNAs in regulating gene expression during high temperature induced bolting, the differential analysis was performed in HT2 vs C2, HT8 vs C8 and HT24 vs C24 and identified 39, 139 and  Table S1). These DE miRNAs potentially targeted 40, 188 and 167 coding genes in HT2 vs C2, HT8 vs C8 and HT24 vs C24, respectively (Fig. 3b, Table S1). Among them, 3, 22 and 21 genes were overlapped between the DEGs and these coding genes potentially targeted by miRNAs during high temperature treatment at 2, 8 and 24 day, respectively (Fig. 3c, Table S2-4), suggesting that the expression of some DEGs were potentially affected by miRNAs (Fig. 3c), although no gene was reported to be involved in the known floral pathway. Among these potentially targeted coding genes, there are some kinase protein family proteins, which may be involved in bolting regulation (Fig. 3c, Table S2-4). These results suggested that the involvements of miRNAs during high temperature treatment possibly affect coding gene expression and bolting in lettuce.

DNA methylation patterns during high temperature-induced bolting
To investigate the overall methylation patterns of lettuce response to high temperature, control and high temperature treated lettuces at different time points were selected for methylation analysis. Each sample produced about 25 Gb data. In the C0 sample, methylcytosine occurred at CHH sites (3.51%), followed by CHG (44.25%) and CG (52.24%) (Fig. 4a). To further explore DNA methylation distribution in different genomic regions, we analyzed the methylation profiles with genes and their flanking sequences. Our results showed that methylation occurs mostly at intergenic regions (Fig. 4b). Moreover, methylation at CGs, CHGs and CHHs within protein-encoding genes exhibited different distribution across sequence contexts. CG and CHG methylation primarily occurred in the genebodies, while methylation at CHHs preferentially occurred in both genebodies and promoters (Fig. 4b).
To examine the DNA methylation patterns in high temperature-induced bolting, we analyzed the changes of DNA methylation after high temperature treatment. It can be seen that with the increase of treatment time, the number of differentially methylated region was increased significantly (Fig. S8). To further explore whether any role of DNA  (Table S5).
We further identified differentially methylated regions (DMRs) using the R packages "DSS" and these DMRs were associated with 7, 16 and 25 protein-coding genes at 2d, 8d and 24d, respectively (Fig. 5b, Table S6). However, only 0, 2 and 5 genes were overlapped between the DMR-associated genes and DEGs at 2d, 8d and 24d during high temperature treatment (Fig. 5b, Table S7-8), indicating a potential role of DNA methylation in regulating gene expression during high temperature-induced bolting.

Discussion
In this study, we performed the transcriptome and epigenome analysis to study the molecular mechanism of high temperature-induced bolting. Our results unravelled that multiple floral pathways are involved in this process. It is well known that there are six major genetic pathways controlling flower transition in Arabidopsis and these flowering signals converge at two regulators, FT and SOC1 [4]. Previous studies showed that high temperature treatment could promote LsFT and LsSOC1 expression in lettuce. LsFT or LsSOC1 knockdown by RNA interference led to insensitive to high temperature and late bolting phenotype [17,18]. But which floral pathways are involved in LsFT and LsSOC1 expression regulation during high temperature-induced bolting is lacking. RNAseq analysis performed in our study revealed that thousands of genes were differentially expressed after high temperature treatment. Among them, most of the flowering time genes are PIF and CO genes. Consistently, recent study also showed that CO LIKE-9 (COL9) was identified as a candidate affecting flowering and bolting in lettuce by utilizing a lettuce mapping population [23], indicating that photoperiod/clock pathway may be a major pathway for high temperature-induced bolting. Consistent with previous study that high temperature treatment results in increased expression level of LsGA2ox1, LsGA3ox1 with corresponding enhanced GA8 and GA1 endogenous levels [19,21], our results showed the gene GA20OX1 responsible for GA biosynthesis was significantly up-regulated after high temperature treatment for 8 days. Moreover, some AP2-like genes in age pathway were down-regulated after high temperature treatment (Fig. 2), indicating that hormonal and age pathways are also involved in high temperature-induced bolting.
Given that epigenome reprogramming plays an important role in gene expression regulation, genome-wide analysis of miRNAs and DNA methylation was performed in this study. We found that a few of DEGs were potentially targeted by DE miRNAs and/or associated with DMRs. Our results showed that about 6% DEGs including some protein kinase family proteins were potential targets of miRNAs during high temperature-induced bolting, although none is known to play a role in at first glance. Previous studies have demonstrated that some protein kinase family proteins, such as CLAVATA1 (CLV1) and CLV2, regulate flower time [24,25]. So, it will be interesting to study whether these protein kinases are involved in high temperature-induced bolting in the future.

Conclusions
Together, our study provided new molecular insights into the mechanism of high temperature-induced bolting by multi-omics analysis. To characterize the functions of the candidate genes will promote the understanding this process and facilitate lettuce breeding in the future.

Plant materials and growth conditions
The heat-sensitive Lactuca sativa L. cultivar GB-30, an easy bolting variety, was used in this study. The lettuce was grown in greenhouse of Beijing University of Agriculture Experimental Station with 14 h light at 20 °C/10 h darkness at 13 °C. When the seedlings have seven true leaves, they were randomly divided into two groups and treated with different temperatures: 33 °C as the treatment and 25 °C as the control. These seedlings were treated for 24 days.

RNA-Seq analysis and enrichment analysis
RNA-Seq was performed by novogene using the Hiseq-PE150 platform. Tips of the lettuce at four different time points (0d, 2d, 8d, 24d) under normal and high temperature treatments were collected separately in three biological replicates and subjected to RNA-Seq analysis. Differentially expressed genes (DEGs) were identified based on lettuce RNA-Seq data from different treatments. All reads were mapped to the lettuce reference genome with default parameters using HISAT2 software (https:// daehw ankim lab. github. io/ hisat2, version 2.1.0), and assemble transcripts by StringTie (http:// ccb. jhu. edu/ softw are/ strin gtie/, version 2.0.3) from the read alignment data. DEGs were identified through the package 'DEseq2' in R software (version 3.6.0) for RNA-Seq data from two treatments with the cutoff: P < = 0.05 and fold change > = 2. The expression patterns of those candidate DEGs among the HT2 vs C0, HT8 vs HT2, and HT24 vs HT8 were displayed according their FPKM using R package 'pheatmap' . Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were conducted for both the upregulated and down-regulated genes using the software "TopGO" and "clusterProfiler" [26], respectively.

Identification of miRNAs and their targets
smallRNA-Seq (sRNA-Seq) was performed by novogene using the Nova-SE50 platform. The lettuce miRNAs were predicted by miR-PREFeR (miRNA PREdiction From small RNA-Seq data) (https:// github. com/ hange lwen/ miR-PREFeR. git, version 0.23.0) default pipline. Differentially expressed (DE) miRNAs between the control and high temperature treated lettuces were identified using the package 'DEseq2' , and their patterns were showed using R package 'pheatmap' . Then, targets of those differentially expressed (DE) miRNAs between the control and high temperature treated lettuces were predicted by TargetFinder (https:// github. com/ carri ngton lab/ Targe tFind er. git).

Learn more biomedcentral.com/submissions
Ready to submit your research Ready to submit your research ? Choose BMC and benefit from: ? Choose BMC and benefit from: