The cucumber is one of the most important vegetables worldwide and is used as a research model for study of phloem transport, sex determination and temperature-photoperiod physiology. The shoot apex is the most important plant tissue in which the cell fate and organ meristems have been determined. In this study, a series of whole-genome small RNA, degradome and transcriptome analyses were performed on cucumber shoot apical tissues treated with high vs. low temperature and long vs. short photoperiod.
A total of 164 known miRNAs derived from 68 families and 203 novel miRNAs from 182 families were identified. Their 4611 targets were predicted using psRobot and TargetFinder, amongst which 349 were validated by degradome sequencing. Fourteen targets of six miRNAs were differentially expressed between the treatments. A total of eight known and 16 novel miRNAs were affected by temperature and photoperiod. Functional annotations revealed that “Plant hormone signal transduction” pathway was significantly over-represented in the miRNA targets. The miR156/157/SBP-Boxes and novel-mir153/ethylene-responsive transcription factor/senescence-related protein/aminotransferase/acyl-CoA thioesterase are the two most credible miRNA/targets combinations modulating the plant’s responsive processes to the temperature-photoperiod changes. Moreover, the newly evolved, cucumber-specific novel miRNA (novel-mir153) was found to target 2087 mRNAs by prediction and has 232 targets proven by degradome analysis, accounting for 45.26–58.88% of the total miRNA targets in this plant. This is the largest sum of genes targeted by a single miRNA to the best of our knowledge.
These results contribute to a better understanding of the miRNAs mediating plant adaptation to combinations of temperature and photoperiod and sheds light on the recent evolution of new miRNAs in cucumber.
Temperature and photoperiod are two of the most important factors affecting plant development [1,2,3]. The cucumber is an important vegetable worldwide and is extensively produced via protected cultivation. Temperature and photoperiod can affect plant architecture, floral bud differentiation, sex expression, stress tolerance, fruit productivity, phytochemicals accumulation, flavours and nutrient formation in cucumbers [4,5,6,7]. Therefore, revealing the molecular mechanisms underlying the developmental modification in response to temperature and photoperiod changes will contribute substantially to cucumber production.
The plant genome responds to environments via epigenetic, transcriptional and post-transcriptional regulations . The epigenetic modifications include chromosome organization, DNA methylation and histone modifications [9, 10]. DNA methylation changes have been revealed in cucumber shoot apexes grown in differential temperatures in previous studies . The DNA methylation patterns also changed in male gametophytes of cucumber plants infected by viroids . Transcriptional regulation of genes encoding photosynthetic enzymes, H+-ATPase, leaf volatile compounds, Argonaute, Dicer-like, RNA-dependent RNA polymerase, NAC transcription factors and superoxide dismutase (SOD) have been proven to play roles in suboptimal temperature and light adaptation and abiotic stress tolerance in cucumbers [13,14,15,16,17,18,19]. Transcriptome sequencing has been performed on cucumbers to analyse the transcripts responsible for sex expression, fruit length, waterlogging reaction and Phytophthora capsici resistance [20,21,22,23]. From the post-transcriptional control level, non-coding small RNA is the best-studied mechanism till now. Small RNAs not only play fundamental roles in post-transcriptional regulation via mRNA cleavage and translation inhibition but also link up these three levels of regulation in the way of epigenetic reprogramming [24, 25].
MicroRNAs (miRNAs) are a class of well-studied small non-coding RNAs. In plants, miRNAs are generated from primary miRNA transcripts (pri-miRNAs) containing a distinctive hairpin structure, which are first trimmed by DICER-LIKE1 (DCL1) to generate miRNA precursors (pre-miRNAs) in the nucleus and then transported to the cytoplasm and further processed to generate ~ 21 nt mature miRNAs [26, 27]. Subsequently, the mature miRNAs are loaded onto the RNA-induced silencing complex (RISC) and guide transcript cleavage (Bartel, 2004; Baulcombe, 2004) or translational repression of their target mRNAs [25, 28]. Plant miRNAs are highly complementary to their targets, which make it possible to predict and verify their targets using bioinformatics and degradome sequencing methods.
miRNAs play important roles in a wide range of plant developmental processes, including plant architecture , leaf morphology , root development , tuberization , the vegetative-reproductive phase change , floral organ identity [28, 34], flowering time , self-incompatibility , cytoplasmic male sterility , plant nutrient homeostasis , and response to biotic and abiotic stresses [39,40,41]. In cucumbers, the miRNAs derived from leaf, root, stem and phloem exudate were analysed in three independent high-throughput sequencing studies [42,43,44]. However, neither miRNAs in the shoot apex nor the miRNAs at different temperature and photoperiod changes were profiled in cucumber plants.
In cucumbers, approximately 10–15 axillary buds or primordia are formed at the shoot apex in seedlings. The environmental conditions of the seedling stage determine the lateral organs such as male or female flowers in the formation of the next 10–15 nodes. Thus, the performances of genes in the cucumber shoot apexes are particularly significant from a horticulture perspective.
In this study, the miRNAs in cucumber shoot apexes were identified. The expression profiles from short and long photoperiods, low and high temperatures were compared, by a series of small RNA sequencing experiments. The targets of these miRNAs were predicted and proven by degradome sequencing. The expression levels of the predicted targets were also profiled by transcriptome sequencing. Based on these data, the roles of miRNAs on the cucumber adaptation to environmental temperature and photoperiod changes were investigated. A newly evolved cucumber specific novel-miRNA and its evolutionary path were highlighted.
sRNA sequencing and miRNA identification
To reveal the roles of miRNAs in plant adaptation to temperature-photoperiod environments, small RNAs were pyrosequenced from the shoot apical tissues of cucumber plants grown in four differential temperature-photoperiod treatments (16 h light at 28 °C /8 h dark at 25 °C, HL; 8 h light at 28 °C /16 h dark at 25 °C, HS; 16 h light at 20 °C /8 h dark at 15 °C, LL; 8 h light at 20 °C /16 h dark at 15 °C, LS). The 16 h/8 h and 8-h/16 h of light/dark photoperiods simulated long- and short-day seasons. The 28 °C/25 °C and 20 °C/15 °C conditions were used to simulate the high and low temperature climates, respectively. A total of 209.25 M reads were produced from the 12 samples (four treatments, three biological replicates). Approximately 16.63 ~ 17.65 M high quality clean reads were obtained for each sample (Additional file 1: Table S1). The majority of the tags ranged in size between 18 and 28 nt, with the 24 nt and 21 nt lengths dominating (Additional file 2: Figure S1). The total reads were converted to non-redundant unique tags and then were aligned to GenBank, Rfam and the cucumber genome. The reads annotated as potential miRNAs and those derived from intronic regions and exon antisense regions as well as the unannotated reads were used for miRNA mining based on the secondary structure (hairpin folding) of their precursors and other criteria defined in the Methods section. The miRNA candidates were aligned to all known plant miRNAs in miRBase 21.0, allowing two-base mismatch and a two-nucleotide overhang. The matched miRNAs were termed as known miRNAs and named after their best-match families. The members within the same miRNA family were alphabetical labelled according to their genome location. A total of 68 known miRNA families containing 164 members were identified (Table 1). The largest family was miR7129, which contains 13 members, followed by miR169, which harbours 12 members. Both 5p and 3p miRNAs were detected for 76 precursors, while the remaining 88 precursors had only one mature miRNA for each. Detailed information, including nomenclature, genome location, mfe, sequences and counts of 5p and 3p mature miRNAs, precursor sequences and their secondary structures are listed in (Additional file 1: Table S2). The remaining 203 miRNAs that could not be incorporated into the known families were termed novel miRNAs. Amongst these, 31 novel-mirs were merged into 10 novel families and the other 172 were single copy miRNA genes (Additional file 1: Table S3). Out of these, 19 novel miRNAs contained both 5p and 3p sequences (Table 2).
Comparative expression patterns of miRNAs in the cucumber shoot tips
The expression levels of known and novel miRNAs were profiled based on tag counts. For multiple precursors sharing the same mature miRNA sequences, only one member of the mature miRNA was used for expression profiling. miR166, miR167 and miR157 were the three highest expressed known miRNAs, which had 200,167-485,949, 91,482–378,695 and 68,360–369,907 copies, respectively, in the cucumber shoot apical tissues (Additional file 1: Table S4). The two highest abundant novel miRNAs were novel-mir166 and novel-mir175, which had 7139–10,992 and 4605–8714 copies, respectively (Additional file 1: Table S5). Another 16 novel miRNAs also had more than 1 thousand copies for each samples. The detailed expression levels of the known and novel miRNAs are listed in Additional file 1: Tables S4 and S5.
The expression level between the treatments was compared using normalized tag counts from the three biological replicates. A total of 8 known and 16 novel miRNAs were differentially expressed (more than two-fold change, p < 10− 3) between at least one pair in an environmental comparison (Fig. 1). Among them, three known miRNAs (miR398, miR530, and miR860) and two novel miRNAs (novel-mir47 and novel-mir87) were suppressed by higher temperature in both short and long photoperiod conditions. One known miRNA (miR7741) and six novel miRNAs (novel-mir3, novel-mir14, novel-mir15, novel-mir125, novel-mir140 and novel-mir321) were suppressed by higher temperature in the long photoperiod but did not fulfil the p-value cut off level in short-photoperiod conditions. Three known and four novel miRNAs were induced by higher temperature. Among these, miR827, novel-mir97 and novel-mir200 were induced specifically in the short photoperiod; miR8125 and novel-mir250 were stimulated only in the long photoperiod, while miR399 and novel-mir55 were significantly up-regulated in both conditions. In low temperatures, no miRNA expression was affected by the photoperiod. However, in high temperatures, two known miRNAs (miR399 and miR827) and three novel miRNAs (novel-mir170, novel-mir173 and novel-mir206) were suppressed by long-day conditions. Interestingly, one novel miRNA (novel-mir51) was down-regulated by both high temperature and long photoperiod considering the fold-changes. However, neither of them fulfilled the p-value < 10− 3 cut off level. When comparing the high-temperature, long-photoperiod to low-temperature, short-photoperiod conditions, both the fold-change and the p-value were improved to being past the cut off level. When comparing the high-temperature, short-day condition (as often happens in natural light greenhouses in the winter season) to the low-temperature, long-day condition (as often happens in low latitude and high elevation regions), the miRNA expression resembled that of the high vs. low temperature in the short-day condition. Comparison of the high-temperature, long-day condition with the low-temperature, short-day condition resembled the summer vs. winter season, showing similar miRNA expression patterns to high vs. low temperature in the long-day condition.
Target prediction by degradome analysis
Using the psRobot  miRNA target prediction tool, 2027 and 6423 targets were assigned to the 66 known and 184 novel miRNAs (Additional file 1: Tables S6 and S7). Using another tool, TargetFinder , we predicted 2775 targets for known miRNAs and 8291 targets for novel miRNAs (Additional file 1: Tables S8 and S9). Assignments of 1007 of the known miRNA targets and 3604 of the novel miRNA targets were the same between the two methods (Additional file 1: Tables S10 and S11).
In plants, miRNAs function mainly by degrading their target mRNAs. Thus, we performed a degradome sequencing analysis to validate the miRNA targets. A total of 17.05 M clean tags were generated from the mixed samples of cucumber shoot tips cultivated in the four temperature-photoperiod conditions. Among them, 12.34 M (72.38%) tags were mapped to the reference cucumber (Chinese Long) genome (Cucumis sativus L. var. sativus cv. 9930 _v20) . A total of 7.70 M (45.13%) tags were mapped to cDNA sense chains and were used for target prediction. A total of 349 mRNAs (394 target sites) were identified as miRNA targets, including 102 genes targeted by 37 known miRNA families and 256 genes targeted by 45 novel miRNAs. Amongst these, novel-mir153 had 232 targets on 201 genes. The cut sites are shown in (Figs. 2, 3 and 4) and the annotation of the targets are list in (Additional file 1: Table S12). A large proportion of miRNAs, including 31 known and 137 novel, had no targets detected, possibly due to the target mRNAs being expressed at low levels that the cut ends were not detected by our pipeline.
The targets were classified by gene ontology (GO) annotation, with the top two terms being “metabolic process” and “cellular process” in the “biological process” cluster. “Cell” and “organelle” were the two most common “cellular component” targets. “Binding” and “catalytic activity” were the two most common terms in the “molecular function” cluster. The “cellular nitrogen compound biosynthetic process”, “macromolecule biosynthetic process”, “aromatic compound biosynthetic process”, “translation”, “gene expression”, “nucleic acid binding” and “binding” were the over-represented GO terms for the targets (Fig. 5, Additional file 1: Table S13). These targets were further categorized by KEGG, and the “Plant hormone signal transduction” pathway was the only significantly over-represented pathway (Additional file 1: Table S14). Among the 629 genes of this pathway, 27 were predicted as targets of miRNAs in cucumber shoot tips.
Target expression profiles based on transcriptome sequencing
To profile the expression of the targets, a transcriptome sequencing project was performed using the same samples as those used for miRNA sequencing. A total of 25.35–27.60 M clean reads were generated from each of the 12 samples (three biological replicates for four treatments). For each sample, 19,045–19,242 genes were expressed (Additional file 1: Table S15). A total of 544 genes were differentially expressed between at least one pair of the six comparisons. The expression levels of the target genes were extracted from the transcriptome data. Amongst these, 14 targets of six miRNAs were differentially expressed among at least one pair of comparison. There were seven targets (one target of novel-mir97, one target of novel-mir132 and five targets of novel-mir153) were up-regulated and the six targets (three targets of miR156/157, one target of miR164 and two targets of novel-mir153) were down-regulated by high temperature, respectively. One gene targeted by miR7741 were down-regulated by long-day condition (Fig. 6, Additional file 1: Table S16). Amongst, novel-mir153 targeted to seven genes, which composed 50% of the total differentially expressed targets. However, the expression of novel-mir153 among the treatments displayed more than two-fold changes but did not reach the p-value< 10− 3 cut-off level. The accumulation of miR164 and novel-mir132 showed no significant different. Only three miRNAs, including miR156/157, miR7741, and cas-novel-mir97, were differentially expressed simultaneously with their targets (Additional file 1: Table S16).
miRNA-target modulating the adaptation to temperature-photoperiod changes
The roles of the miRNA-target network in modulating the plant adaptation to temperature-photoperiod changes were investigated based on the functional annotations of the expressional affected miRNAs and targets. Of the 24 temperature and photoperiod affected miRNAs, 12 have 32 targets detected by degradome analysis (Additional file 1: Table S17). Five out of the 32 targets were differentially expressed (Additional file 1: Table S17). However, only miR156/157 and its three SBP-box targets have shown negative correlated expression pattern (Fig. 7a). MiR156/157 is well studied and negatively controls plant juvenile to adult phase changes, flowering, and other developmental processes via SBP-box targets [29, 48, 49]. By down-regulation of miR156/157, the long-day condition promoted the juvenile to adult phase changing, flower initiation and other developmental and metabolic processes in cucumber.
Though did not achieve the significantly different cut-off level, the expression of novel-mir153 among the treatments displayed negative correlation pattern with its five high-temperature induced targets (Fig. 7b). These targets include an aminotransferase, an ethylene-responsive transcription factor, a senescence-related protein and an acyl-CoA thioesterase coding gene. As ethylene is the central regulator signal for flower sex development in cucumber, the novel-mir153 has the potential to play roles in sex ratio regulation.
A newly evolved novel miRNA targets a large number of genes
Noticeably, a novel miRNA (novel-mir153) targeted 2087 mRNA targets, predicted by the intersection of the psRobot and TargetFinder tools, which accounted for 45.26% of the total targets in this plant. Amongst, 232 cut sites were proved by degradome analysis, which accounted for 58.88% of the total degradome results. This is the largest sum of genes targeted by a single miRNA in plants. By BLAST on NCBI, we did not find a similar hairpin structure sequence from other organisms. When this novel-mir153 precursor was traced back to the cucumber (Chinese Long 9930) genome, we noticed that it was generated from the third intron of the alcohol dehydrogenase-like 6 gene (Csa1G044860). We then retrieved its homologue genes from seven other cucurbit genomes, including two cucumbers (Gy14 and PI183967), melon (DHL92), two watermelons (97,103 and Charleston Gray), Cucurbita maxima (Rimu) and Cucurbita moschata (Rifu) . The coding sequences were highly conserved between these plants (Additional file 2: Figure S2). However, the intron regions have diverged between these species. The 3′-end part of the third intron was conserved between these five species, indicating that these introns were generated from a common ancestor. The direct repeat sequence at the borders indicates the intron was probably derived from transposable element insertion. The 5′-end of the third introns was lost in pumpkins or gained by the progenitor of cucumber, melon and watermelon after divergence from pumpkins. After this, two regions highly diverged between Cucumis and Citrullus. One of them evolved the novel-mir153 in cucumber (9930 and wild) by at least three steps (Fig. 8a and b). The secondary structure of the precursor of novel-mir153 showed perfect hairpin folding (Fig. 8b). This novel-miRNA is very young, evolved from the ancestors of melon and cucumber to cucumber cultivar Gy14 and then to cucumber cultivar Chinese long by two rounds of “AGA” duplications. The novel-mir153 and its targets matches are shown in Figs. 4 and 8c, perfectly meeting the criterion for regulation.
We checked the functions of these targeted genes, 182 out of the 201 genes were functionally annotated. These targets covered quite wide biological functions, including RNA transcription, protein translation, protein, lipid and carbohydrate metabolic, growth regulation, environmental response, transporters, etc. (Additional file 1: Table S12). It is worth mentioning that 23 transcription factors, nine kinase and seven cyclin-like genes were putative targets of novel-mir153, indicating this novel miRNA plays an important role in cucumber development.
Three previous studies by Ling et al., (2017); Mao et al., (2012) and Martinez et al., (2011) have identified 25, 29 and 31 known and 7, 2 and 49 novel families of miRNAs in cucumber plants (leaf, root and stem tissues) [42,43,44]. Among these previously identified miRNAs, 27 families of known miRNAs were also detected in our present study. The other eight known miRNA families (miR170, miR319, miR858, miR894, miR1030, miR1515, miR2911 and miR2916) identified in these previous reports were not detected in our study. However, our present study detected 41 new known miRNA families that have never been identified in previous reports (Table 1). These differences present the differential tissue patterns of miRNAs between the leaf, root, and stem in the previous reports and the shoot apex in our present study.
The expressional profiling showed that the miRNAs were primarily affected by temperature rather than by photoperiod. The temperature had an epistatic effect on photoperiod (Fig. 1). When cultivated in low temperature, the change of photoperiod did not affect the expression levels of miRNAs in these cucumber shoot tips. However, when temperature rose above a suitable level, the photoperiod significantly influenced the expression of many miRNAs. These results give guidance for cucumber production, especially for greenhouse-based cultivation, in that temperature has a greater impact than photoperiod and merits greater attention. When cucumbers are growing in a greenhouse during the winter season, if the temperature is not high enough, the extension of the photoperiod will be insufficient. Only when the temperature are high enough will a prolonged photoperiod stimulate flowering and other developmental processes. However, cucumber has many ecotypes. The “Chinese long (9930)” used in this study is a cultivar successfully used in greenhouse cultivation in winter and early spring seasons. The above findings should be cautiously used if extended to other ecotype cucumbers such as the “Xishuangbanna” variety.
The functional annotations of the differentially expressed miRNAs provide some answers for the performance of cucumber within differential temperature-photoperiod environments (Additional file 1: Table S18, Additional file 2: Figure S3). For example, down-regulation of miR7741 and seven other miRNAs releases the activity of targets involved in transcription initiation, elongation, mRNA processing, degradation, and protein translation in high temperature. This can explain the fact that these plants are more vigorous in higher temperatures. Under miRNA regulation, the expression or translation of proton-dependent transporter genes and vacuolar transport genes were biased towards high or low temperatures, respectively. This can contribute to the phenomena of mineral element absorbance and transport being more efficient in high temperatures, while the carbohydrate and alkaloid storage in vacuoles is more efficient in low temperatures. It is well known that low temperatures can cause fungal diseases in cucumber. The suppression of the resistance and defence related genes via miRNAs can partially explain this phenomenon. However, due to the limited sequencing depth, a lot of predicted targets were not proven in the degradome analysis. Therefore, deeper degradome sequencing is needed to prove these functional deductions.
One of the important outputs of the temperature-photoperiod influence is the flower sex ratio changes. DNA methylation has been proven to play a role in this process . DNA methylation achieves its functions by transcriptional regulation of genes. miRNA is a more rapid and direct regulator of gene expression. In this study, miR156/157 differentially expressed and regulated the SBP-box transcription factors. SBP-box genes regulate flowering time, male fertility and gynoecium patterning in Arabidopsis [48, 49] but did not play a strong role in sex determination, according to the studies that were carried out using hermaphrodite flowering plants. Thus, we cannot exclude the possibility of these SBP-box transcription factors playing roles in flower sex ratio regulation in cucumber. Beside the significant differentially expressed miRNAs, those miRNAs such as novel-mir153 target Ethylene-responsive transcription factors whose expression was significantly affected by temperature and photoperiod. As ethylene is the central regulator signal for flower sex development in cucumber, these genes and the corresponding miRNAs have the potential to play roles in sex ratio regulation.
There are three models that describe miRNA origin and evolution: origination from a random formation of hairpins in non-coding sections such as introns or intergenic regions; origination from inverted duplications of protein-coding sequences; and origination from the miniature inverted-repeat transposable elements [51, 52]. In this study, we found a solid example of a new miRNA (novel-mir153) that originated from an intron. The comparison of the genes from pumpkins, watermelons, melons and cucumbers displayed the step-by-step evolution of this miRNA. Because novel-mir153 is a quite newly evolved miRNA and absent in other plants, it is not an indispensable element for the plant life cycle. However, it does target a surprisingly large number of target genes. The detailed function of novel-mir153 and its contribution to cucumber evolution is worthy further investigation.
Materials and methods
Plant materials and RNA extraction
The cucumber inbred line 9930 was used in this study. Seedlings were grown in artificial climate chambers with conditions set as follows: 16 h light in 28 °C /8 h dark in 25 °C (high temperature, long day, HL), 8 h light in 28 °C /16 h dark in 25 °C (high temperature, short day, HS), 16 h light in 20 °C /8 h dark in 15 °C (low temperature, long day, LL), and 8 h light in 20 °C /16 h dark in 15 °C (low temperature, short day, LS). A 70% humidity was used for all of these plants. When four true-leaves were unfolded, shoot apical tissues (1 mm) were dissected under a microscope and snap-frozen in liquid nitrogen and kept at − 80 °C for further use. In each experiment, more than 500 shoot apices and three biological replicates were used. Total RNA was extracted using the TRIzol reagent (Invitrogen, USA). DNase (Promega, USA) was used to remove potential DNA contamination.
Small RNA library construction and sequencing
The RNA samples were quantified and equalized. A total of 30 μg of RNA was resolved on denatured polyacrylamide gels. Gel fragments with the size range of 18–30 nt were excised and recovered. These small RNAs were ligated with 5′ and 3’RNA adapters using T4 RNA ligase. The adapter-ligated small RNAs were subsequently transcribed into cDNA by Super-Script II Reverse Transcriptase (Invitrogen) and amplified using primers specific for the ends of the adapters. The amplified cDNA products were purified and finally sequenced using Illumina sequencing technology (BGI, Shenzhen, China).
Identification of known and novel miRNAs
The adapter sequences, impurities, and sequences with less than 18 nt or more than 30 nt were filtered out from the raw sequence reads. For miRNA identification, the total reads were converted to non-redundant unique tags to filter out duplicates. The unique tags were aligned to GenBank, Rfam and the cucumber genome  and therefore annotated as rRNA, miRNA, snRNA, snoRNA, tRNA, repeat, exon sense, exon antisense, intron sense, intron antisense, and unannotated sequences. The putative miRNA tags, unannotated reads and those derived from the intron regions and exon antisense regions were used for miRNA mining based on the criteria as follow: hairpin pre-miRNAs that can fold into secondary structures and mature miRNAs that were present in one arm of the hairpin precursors; the 5-p and 3-p mature miRNAs present 2-nucleotide 3′ overhangs; hairpin precursors lacking large internal loops or bulges; the secondary structures of the hairpins were steady, with a free energy of hybridization lower than or equal to − 18 kcal/mol; and the copy number of mature miRNAs with predicted hairpins should be greater than 5 in the alignment result. The candidates were first aligned to all known plant miRNAs in miRBase 21.0 allowing a two-base mismatch and a two-nucleotide overhang. The matched miRNAs were termed as known miRNAs and named after their best-match families. The members within the same miRNA family were alphabetical labelled according to the genome location of their precursors. The rest of the miRNAs that could not be incorporated into known families were termed as novel miRNAs.
The expression of miRNAs was produced by summing the total count of tags with no more than 3 mismatches on the 5′ and 3′ ends and no mismatches in the middle from the alignment result. The differentially expressed miRNAs were calculated with the following procedures. First, the expression of miRNAs was normalized to obtain the expression of transcript per million (TPM). The normalization formula was as follows: Normalized expression = Actual miRNA count/Total count of clean reads× 1,000,000. Second, the fold-change and P-value were calculated from the normalized expressions.
Degradome library construction and target identification
Total RNA was extracted from the same samples as those used for sRNA sequencing. Approximately 200 μg of total RNA was polyadenylated using the Oligotex mRNA mini kit (Qiagen). A 5′ RNA adapter was added to the cleavage products (which possessed a free 5′-monophosphate at their 3′ termini) using the T4 RNA ligase (Takara). Then, the ligated products were purified using the Oligotex mRNA mini kit (Qiagen) for reverse transcription to generate the first strand of cDNA using an oligo dT primer via SuperScript II RT (Invitrogen). After the cDNA library was amplified for 6 cycles (94 °C for 30 s, 60 °C for 20 s, and 72 °C for 3 min) using Phusion Taq (NEB), the PCR products were digested with restriction enzymes and further ligated to double-stranded DNA adapters. The ligated products were subjected to PCR amplification and gel purification and finally used for high-throughput sequencing with the Illumina HiSeq 2000.
Low-quality sequences and adapters were removed, and the unique sequence signatures were aligned to the cucumber (9930) genome  using SOAP software . CleaveLand was used to detect potentially cleaved targets . The tags mapped to cDNA sense strands were adopted to predict cleavage sites. The miRNA-mRNA pairs and p-values were analysed using PAREsnip, p-values less than 0.05 were used for the t-plot . All alignments with scores not exceeding 4 and possessing the 5′ end of the degradome sequence coincident with the tenth and eleventh nucleotides of the miRNA complementary were retained. The targets were GO annotated using Blast2GO .
Target gene profiling
Total RNAs (10 μg) were subjected to poly-A selection, fragmentation, random priming and cDNA synthesis with the Illumina Gene Expression Sample Prep kit (CA, USA). The cDNA fragments were end repaired, ligated to adapters and then enriched by PCR. The fragments were purified with 6% TBE PAGE gel electrophoresis. After denaturation, the single-chain fragments were fixed onto the Solexa Sequencing Chip (Flowcell) and consequently grown into single-molecule cluster sequencing templates through in situ amplification . Double-end pyrosequencing was performed on the Illumina HiSeq 2000 with the read lengths of 90 bp for each end. The clean reads were aligned to cucumber reference genes . Gene expression levels were calculated using the FPKM. Statistical comparison between treatments was performed using the data from three biological replicates. Unigenes were considered to be differentially expressed when the mean FPKM between treatments displayed a more than two-fold change with a probability of more than 0.8. The expression profiles of the target genes were extracted from the transcriptome data.
Schultz TF, Kay SA. Circadian clocks in daily and seasonal control of development. Science. 2003;301:326–8.
Savvides A, Dieleman JA, van Ieperen W, Marcelis LF. A unique approach to demonstrating that apical bud temperature specifically determines leaf initiation rate in the dicot Cucumis sativus. Planta. 2016;243:1071–9.
Castellano M, Martinez G, Marques MC, Moreno-Romero J, Kohler C, Pallas V, et al. Changes in the DNA methylation pattern of the host male gametophyte of viroid-infected cucumber plants. J Exp Bot. 2016;67:5857–68.
Bi HG, Wang ML, Jiang ZS, Dong XB, Ai XZ. Impacts of suboptimal temperature and low light intensity on the activities and gene expression of photosynthetic enzymes in cucumber seedling leaves. Ying Yong Sheng Tai Xue Bao. 2011;22:2894–900.
Zhang XM, Yu HJ, Sun C, Deng J, Zhang X, Liu P, et al. Genome-wide characterization and expression profiling of the NAC genes under abiotic stresses in Cucumis sativus. Plant Physiol Bioch. 2017;113:98–109.
Gan D, Zhan M, Yang F, Zhang Q, Hu K, Xu W, et al. Expression analysis of argonaute, Dicer-like, and RNA-dependent RNA polymerase genes in cucumber (Cucumis sativus L.) in response to abiotic stress. J Genet. 2017;96:235–49.
Zhou Y, Hu L, Wu H, Jiang L, Liu S. Genome-wide identification and transcriptional expression analysis of cucumber superoxide dismutase (SOD) family in response to various abiotic stresses. Int J Genomics. 2017;2017:7243973.
Deng J, Yu HJ, Li YY, Zhang XM, Liu P, Li Q, et al. Leaf volatile compounds and associated gene expression during short-term nitrogen deficient treatments in Cucumis seedlings. Int J Mol Sci. 2016;17:E1713.
Wang L, Cao C, Zheng S, Zhang H, Liu P, Ge Q, et al. Transcriptomic analysis of short-fruit 1 (sf1) reveals new insights into the variation of fruit-related traits in Cucumis sativus. Sci Rep. 2017;7:2950.
Zhang Y, Zhao G, Li Y, Mo N, Zhang J, Liang Y. Transcriptomic analysis implies that GA regulates sex expression via ethylene-dependent and ethylene-independent pathways in cucumber (Cucumis sativus L.). Front Plant Sci. 2017;8:10.
Xu X, Chen M, Ji J, Xu Q, Qi X, Chen X. Comparative RNA-seq based transcriptome profiling of waterlogging response in cucumber hypocotyls reveals novel insights into the de novo adventitious root primordia initiation. BMC Plant Biol. 2017;17:129.
Mansfeld BN, Colle M, Kang Y, Jones AD, Grumet R. Transcriptomic and metabolomic analyses of cucumber fruit peels reveal a developmental increase in terpenoid glycosides associated with age-related resistance to Phytophthora capsici. Hortic Res. 2017;4:17022.
Zhang X, Zou Z, Zhang J, Zhang Y, Han Q, Hu T, et al. Over-expression of sly-miR156a in tomato results in multiple vegetative and reproductive trait alterations and partial phenocopy of the sft mutant. FEBS Lett. 2011;585:435–9.
Wei X, Zhang X, Yao Q, Yuan Y, Li X, Wei F, et al. The miRNAs and their regulatory networks responsible for pollen abortion in Ogura-CMS Chinese cabbage revealed by high-throughput sequencing of miRNAs, degradomes, and transcriptomes. Front Plant Sci. 2015;6:894.
He H, Liang G, Li Y, Wang F, Yu D. Two young MicroRNAs originating from target duplication mediate nitrogen starvation adaptation via regulation of glucosinolate synthesis in Arabidopsis thaliana. Plant Physiol. 2014;164:853–65.
Mao W, Li Z, Xia X, Li Y, Yu J. A combined approach of high-throughput sequencing and degradome analysis reveals tissue specific expression of microRNAs and their targets in cucumber. PLoS One. 2012;7:e33040.
Fahlgren N, Howell MD, Kasschau KD, Chapman EJ, Sullivan CM, Cumbie JS, et al. High-throughput sequencing of Arabidopsis microRNAs: evidence for frequent birth and death of MIRNA genes. PLoS One. 2007;2:e219.
We thank the National Mid-term Genebank of Vegetable Germplasm Resources, the National Science & Technology Infrastructure, and the Beijing Research Station of Vegetables Crop Gene Resource & Germplasm Enhancement, Ministry of Agriculture, P.R. China.
This work was supported by grants from the National Key Research and Development Programme of China (2016YFD0100204), the National Natural Science Foundation of China (31171961), the National Key Technology R & D Programme of China (2013BAD01B04), and the Science and Technology Innovation Programme of the Chinese Academy of Agricultural Sciences (CAAS-ASTIP-IVFCAAS).
Availability of data and materials
The RNA sequencing datasets generated or analysed during the current study are available in NCBI Sequence Read Archive (SRA) with the accession numbers SRP104621 and SRP104378 and in NCBI Gene Expression Omnibus (GEO) with the accession number GSE117867.
Xiaohui Zhang, Yunsong Lai and Wei Zhang contributed equally to this work.
Authors and Affiliations
Key Laboratory of Biology and Genetic Improvement of Horticultural Crops, Ministry of Agriculture; Institute of Vegetables and Flowers, Chinese Academy of Agricultural Sciences, Beijing, 100081, China
Xiaohui Zhang, Yunsong Lai, Wei Zhang, Jalil Ahmad, Yang Qiu, Xiaoxue Zhang, Mengmeng Duan, Tongjin Liu, Jiangping Song, Haiping Wang & Xixiang Li
Institute of Pomology & Olericulture, Sichuan Agricultural University, Chengdu, 611130, China
XHZ, XL designed the experiments, analysed the data, and wrote the manuscript; WZ, YL, XHZ performed the experiments; XXZ, MD, TL, JS, HW, YQ participated in plant material growth and experiments; and JA participated in manuscript writing. All authors read and approved the final manuscript.
Table S1. summary of sRNA sequencing data. Table S2. summary of known miRNAs. Table S3. summary of novel miRNAs. Table S4. expression levels of known miRNAs in cucumber shoot tips grown in four temperature-photoperiod conditions. Table S5. expression levels of novel miRNAs in cucumber shoot tips grown in four temperature-photoperiod conditions. Table S6. targets of known miRNAs predicted by psRobot. Table S7. targets of novel miRNAs predicted by psRobot. Table S8. targets of known miRNAs predicted by TargetFinder. Table S9. targets of novel miRNAs predicted by TargetFinder. Table S10. targets of known miRNAs predicted by both psRobot and TargetFinder (intersection). Table S11. targets of novel miRNAs predicted by both psRobot and TargetFinder (intersection). Table S12. targets identified by degradome sequencing. Table S13. GO term enrichment of targets. Table S14. KEGG classification of miRNA targets. Table S15. statistics of transcriptome sequencing for each sample. Table S16. differentially expressed degradome confirmed targets. Table S17. targets (degradome proven) of the temperature and photoperiod affected miRNAs. Table S18. functional annotations (computer predicted targets) of the temperature and photoperiod affected miRNAs. (XLSX 1042 kb)
Figure S1. Length distribution of small RNAs in cucumber shoot tips. Figure S2. Alignment of exons of alcohol dehydrogenase-like 6 genes Figure S3. Schematic diagram of the miRNA-target network modulating the plant adaptation to temperature changes. (DOCX 590 kb)
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Zhang, X., Lai, Y., Zhang, W. et al. MicroRNAs and their targets in cucumber shoot apices in response to temperature and photoperiod.
BMC Genomics19, 819 (2018). https://doi.org/10.1186/s12864-018-5204-x