Solanum lycopersicum “Micro-Tom” (TOMJPF00001) and “Moneymaker” cultivars were obtained from the National Bioresource Project (MEXT, Japan) through the TOMATOMA database (Saito et al., 2011). Fruits were harvested at the mature green stage (before onset of autocatalytic ethylene production), washed in commercial bleach (1:10 dilution of sodium hypochlorite) and sorted to ensure uniform size, color and absence of defects or damage. This timing of harvest was well thought out to avoid the effect of large amounts of ethylene, which are produced at later maturity stages, during cold stress tests. For each cultivar, two groups of 50 fruits were used to characterize the ripening behavior at 20°C; the first group were treated (2–3 times a week) with 2 µLL^-1 1-MCP for 12 h, while the second group were a non-treated control. For cold stress tests, two groups of 50 fruits were also used; the first group was treated with 1-MCP as described above while the other group was a non-treated control. “Micro-Tom” fruits were stored at 5°C for 14 d before being transferred to 20°C for up to 21 d. Three separate storage trials were carried out on “Micro-Tom” fruits. “Moneymaker” fruits were also stored at 5°C but for 21 d; after every 7 d, 10 fruits were transferred to 20°C to observe CI symptoms. 1-MCP treatments were carried out to keep the fruits insensitive to ethylene. To release 1-MCP gas, SmartFresh™ powder (AgroFresh, PA, USA) was dissolved in water and soda lime was added in the sealed treatment containers to reduce CO₂ accumulation. For DNA methylation inhibitor treatment, “Moneymaker” fruits at the mature green stage were used. The fruits were injected with 100 µL of 50 mM 5-azacytidine aqueous solution into the columella at harvest and every 7 d during storage at 5°C. The negative controls were injected with the same amount of distilled water. In all treatments, pericarp samples of three replicate fruits were collected, frozen in liquid nitrogen and stored at -80°C for future analysis.

Peel color measurements were carried out on four evenly distributed equatorial sites using a Konica Minolta Color Reader CR-10 (Konica Minolta, Tokyo, Japan). The Hunter lab parameter a*, which is a measure of greenness or redness, was recorded and then expressed as the mean of six replicate fruits.

Individual fruits were incubated at the respective storage temperatures for 1 h. Headspace gas (1 mL) was then withdrawn and injected into a Shimadzu 5890 series gas chromatograph (Shimadzu, Kyoto, Japan) equipped with a flame ionization detector (200°C) and an activated alumina column (80°C). Ethylene production rates were expressed as the mean of 10 replicate fruits.

Fruits were assessed visually for severity of CI symptoms based on a five-point scale (0 = no injury; 1 = < 10%; 2 = 11 to 25%; 3 = 26–40%, and 4 = > 40%) consisting of three parameters: surface pitting, uneven ripening, and decay (Vega-García et al., 2010; Albornoz et al., 2019). The CI index was then calculated by determining the average of the injury levels of surface pitting, uneven ripening, and decay. For each time point, 10 fruits were evaluated individually, and the CI indexes were averaged.

Pericarp samples of “Micro-Tom” fruits collected at harvest (0 d), after 14 d storage at either 5°C or 20°C, and after 21 d at 20°C or 14 d at 5°C followed by 7 d at 20°C (three fruit per treatment) were used for RNA-Seq analysis. Total RNA was extracted from 100 mg samples using the RNeasy^® Plant Mini Kit (Qiagen, Hilden, Germany), treated with DNase I (Nippon Gene, Tokyo, Japan) to remove genomic DNA contamination and further purified with FavorPrep after Tri-Reagent RNA Clean-up Kit (Favorgen Biotech. Co., Ping-Tung, Taiwan). Paired-end libraries were then constructed using NEBNext^® Ultra^™ II Directional RNA Library Prep Kit (New England Biolabs), and sequencing was performed on an Illumina Novaseq 6000 platform (Illumina, Inc.).

The sequenced reads were analyzed primarily on the Galaxy online platform ¹ . Trimming was first carried out to exclude both low quality sequences and adapter sequences. The trimmed reads were then mapped to the reference S. lycopersicum genome (SL4.0), and mapped reads were counted using the featureCounts tool. Gene expression levels were then normalized as transcripts per kilobase million (TPM) reads. Differentially expressed genes (DEGs) were obtained by comparing the expression levels in samples after 14 d of storage at either 5°C or 20°C with those of at-harvest (0 d) samples on the iDEP (v. 0.96) web-based toolkit (Ge et al., 2018). Three criteria were used to detect DEGs: (i) TPM ≥ 1.0 in either of the three replicate samples, (ii) false discovery rate ≤ 0.01, and (iii) two-fold increase or increase in expression levels. Weighted gene co-expression network analysis (WGCNA) method (Langfelder and Horvath, 2008) was then employed to generate clusters of highly correlated genes with 8 and 0.15 as thresholding power and tree-cut parameters, respectively. Significantly enriched gene ontology (GO) terms and KEGG pathways were established using the ShinyGO (v. 0.77) web-based toolkit (Ge et al., 2020). The cut-off for significantly enriched terms was P < 0.05.

cDNA was synthesized from 1 µg of clean DNase I-treated RNA (same samples as those used for library construction) using the SuperScript^® III First-Strand Synthesis SuperMix for qRT-PCR kit (Invitrogen). Gene-specific primers ( Supplementary Table 1 ) were designed using the Primer3 online software (version 0.4.0 ² ). Gene expression of three biological replicates was examined on a Stratagene Mx3005P Real-Time QPCR System (Agilent Technologies, Santa Clara, CA, USA) using KOD SYBR^® qPCR Mix (Toyobo, Osaka, Japan). SlActin (Solyc03g078400) was used as the housekeeping gene after examining its constitutive expression pattern from the RNA-Seq data. Relative expression values were calculated using the 2^-ΔΔCt method with at-harvest (0 d) samples calibrated as 1.

Genomic DNA was extracted from the pericarp of “Moneymaker” fruits (three replicates for each time point) using the Nucleospin^® Plant II kit (Takara, Shiga, Japan).The genomic DNA (1 µg) was then digested at 37°C overnight with McrBC (Takara), before performing qPCR analysis as described in section 2.7 with 20 ng digested DNA as a template. Undigested gDNA samples were used as controls. Relative methylation was then calculated as 2^{Ct(digested) - Ct(undigested)}, such that higher relative McrBC-qPCR signals correspond to higher methylation levels.

Three-week old “Moneymaker” seedlings were transferred to 5°C or 20°C for up to 7 d. Leaves were collected at 0, 3, 5, and 7 d for RNA extraction and qPCR analysis as described in section 2.7.

Data obtained in this study were subjected to statistical analysis using R v.3.4.0 ³ . Differences in ethylene production rates, peel color, CI index, methylation levels and gene expression levels were determined using ANOVA followed by Tukey’s post hoc tests.

Non-treated tomato fruits started to lose their green peel color after 7 d at 20°C, attaining a uniform red color after 21 d ( Figure 1A ; Supplementary Figure 1 ). This was evidenced by an increase in parameter a* from negative values to positive values ( Figure 1C ; Supplementary Figure 1 ). The changes in peel color were consistent with ethylene production rates ( Figure 1D ), which increased from 7 d in fruits during storage at 20°C and peaked at 14 d in a typical climacteric pattern. For the 1-MCP-treated fruits, both peel color changes and ethylene production rates did not change significantly in the first 14 d at 20°C; although ethylene production rates increased later on, the peak was 58% lower than that of non-treated fruits, and the fruits never attained a red color throughout the 35 d storage period.

During storage at 5°C, both non-treated and 1-MCP-treated fruits did not show any noticeable change in peel color ( Figures 1B, C ), and the ethylene production rates also did not change significantly ( Figure 1D ). On rewarming at 20°C, ethylene production rates increased, with a lag of about 6 d in 1-MCP-treated fruits ( Figure 1D ). The rewarmed fruits, however, failed to ripen normally, developed pitted surfaces and showed decay symptoms as evidenced by the high CI index values ( Figures 1B, E, F ; Supplementary Figure 1 ).

RNA-Seq analysis was then performed to identify gene expression changes during normal fruit ripening (at 20°C) and in response to cold stress (triggered by storage at 5°C for 14 d) in “Micro-Tom” fruits. By comparing samples collected after 14 d storage at 20°C and 5°C with at-harvest (0 d) samples, we identified 10,042 DEGs. Interestingly, we found more upregulated and downregulated genes at 5°C (8,406) than at 20°C (4,814) ( Figure 2A ). Subsequent hierarchical clustering of the DEGs against temporal expression patterns in at-harvest, stored and rewarmed samples outlined 9 major modules ( Figure 2B ; Supplementary Tables 2 - 10 ). Among them, modules I, II, III, IV, VII and VIII comprised genes which were differentially expressed both at 20°C and 5°C, and 1-MCP treatment affected their expression levels at both temperatures; this indicated that they are regulated by either ethylene or low temperature. Genes in modules V and VI were up- or down-regulated during storage at 20°C, a change that was reversed by 1-MCP treatment; these expression changes were suppressed during storage at 5°C, but they recovered (at least to some extent) upon rewarming. Interestingly, module IX genes displayed insignificant expression changes during storage at 20°C, but they were up- or down-regulated at 5°C in both 1-MCP-treated and non-treated samples with a reversion to original levels upon rewarming at 20°C.

To further understand the molecular changes that are induced during cold response versus fruit ripening, we performed GO term and KEGG pathway-based enrichment analyses of the DEGs in each of the modules identified above. GO analysis confirmed that genes in module VI were associated with ‘fruit ripening’ ( Figure 2C ), especially the ‘carotenoid biosynthesis pathway’ ( Figure 2D ). On the other hand, module IX genes were dominated by functions such as ‘gene expression’, ‘translation’, ‘ribosome biogenesis’, ‘proteasome’ and ‘NcRNA metabolism’ ( Figures 2C, D ). Other functions that were uniquely enriched in module IX genes included ‘spliceosome’, ‘organonitrogen biosynthesis’, and ‘linolenic acid metabolism’.

All in all, our results indicated that there was a divergence in the molecular adjustments triggered by cold response and those involved in fruit ripening in tomatoes. Module VI genes, in particular, were mostly involved in changes that lead to fruit ripening and they were under the influence of ethylene signaling (i.e., affected by 1-MCP treatment). Conversely, module IX genes were regulated by low temperature while ethylene signaling had little or no significant influence on their expression levels (i.e., unchanged by ethylene at 20°C and unaffected by 1-MCP at 5°C.

Genes that have been commonly associated with fruit ripening were selected from among the DEGs for further examination and to assert the inhibitory effect of 1-MCP. Most of these genes belonged to module VI; they were upregulated at 20°C in accordance with ethylene production rates and 1-MCP treatment inhibited this change ( Figure 3A ). However, one gene (SlSGR, Solyc12g056480), which encodes magnesium dechelatase that is involved in chlorophyll degradation (Shimoda et al., 2016), was upregulated only at 5°C and 1-MCP treatment did not inhibit this expression change, indicating that it responded only to low temperature and was independent of ethylene signaling.

Six of the genes in Figure 3A were then selected for validation by qPCR analysis. These included ripening inhibitor (SlRIN) and non-ripening (SlNOR) which encode established molecular regulators of fruit ripening, 1-aminocyclopropane-1-carboxylic acid synthase 2 (SlACS2) which encodes an ethylene biosynthetic enzyme, polygalacturonase 2a (SlPG2a) which is involved in fruit softening, phytoene synthase 1 (SlPSY1) which is involved in carotenoid metabolism and the aroma volatile production-associated alcohol acyl transferase (SlAAT). As expected, the expression levels of these six genes increased markedly in non-treated fruits during storage at 20°C ( Figure 3B ), which corresponded with the ethylene production patterns observed earlier ( Figure 1D ). In 1-MCP-treated fruits at 20°C, little or insignificant changes in expression levels were registered during the first 14 d; increases however occurred at 21 d as ethylene production levels rose. The expression levels of these ripening-associated genes (except for SlSGR) remained mostly unchanged during storage at 5°C and even after 1 d of rewarming at 20°C; significant transcript accumulation occurred after post-cold storage for 7 d at 20°C. Together, these results reaffirmed that ethylene signaling regulates the ripening process in tomatoes and 1-MCP treatment effectively blocks ethylene-mediated regulation of the genes involved.

Cold stress results in increased levels of reactive oxygen species which cause oxidative damage in fruits and vegetables (Sevillano et al., 2009; Valenzuela et al., 2017). To avoid or tolerate this oxidative damage, several enzymatic reactions catalyzed by lipoxygenases (LOX), peroxidases (POX), alternative oxidases (AOX), catalases (CAT), and superoxide dismutases (SOD) are generated in the affected plants. We therefore examined further the expression patterns of genes encoding these antioxidant defense enzymes from among the DEGs. Interestingly, out of the 18 genes that we found ( Figure 4A ), a majority (13 genes, 72%) belonged to module IX as they showed insignificant changes in expression levels during storage at 20°C, but they were differentially expressed at 5°C with or without 1-MCP treatment. Furthermore, only 2 out of these 13 genes were downregulated while the rest (11 genes) were upregulated in fruits stored at 5°C.

We then studied the expression of SlAOX1a (Solyc08g075540) and SlPOX65 (Solyc08g007150) in both “Micro-Tom” and “Moneymaker” fruits during storage at 5°C and 20°C via qPCR analysis. As indicated in Supplementary Figure 1 , “Moneymaker” fruits also exhibited CI symptoms during rewarming at 20°C following storage at 5°C for 7, 14, and 21 d. SlAOX1a and SlPOX65 transcript levels notably increased in both “Micro-Tom” and “Moneymaker” fruits during storage at 5°C with or without 1-MCP treatment, and decreased immediately after the fruits were rewarmed to 20°C ( Figures 4B, C ). During storage at 20°C, the expression levels remained mostly unchanged except for a relatively slight but relevant increase in SlAOX1a transcript levels in “Moneymaker” fruits. Therefore, it appears that the adaptive responses to cold stress-triggered ROS accumulation in tomato fruits are typically not influenced by ethylene signaling.

Term enrichment analysis of our RNA-Seq data ( Figures 2C, D ) indicated that one of the molecular responses to low temperature alone (with respect to fruit ripening) was adjustments among ribosomal proteins and the ribosome biogenesis pathway in addition to ‘NcRNA metabolism’. Visualization of enriched pathways suggested that low temperature triggers comprehensive changes in both large and small ribosomal subunits as well as ribosomal biogenesis factors ( Figure 5A ). After a closer look at the DEGs, we found 30 genes belonging to module IX ( Figure 5B ). Their expression levels did not change significantly during storage at 20°C, but 18 of them were upregulated while the remaining 12 were downregulated during storage at 5°C both in 1-MCP-treated and non-treated fruits. qPCR analysis further confirmed the upregulation of SlZFP622/REIL1 (Solyc08g006470), encoding a ribosome biogenesis factor, during storage at 5°C with transcript levels peaking at 7 d in both “Micro-Tom” and “Moneymaker” tomatoes ( Figure 5C ).

Another category that stood out among the DEGs responding to low temperature alone (module IX) was ‘spliceosome’ ( Figure 2D ), a protein complex in which introns are removed from immature mRNAs to generate uninterrupted open reading frames for translation or to produce different splicing variants of the same gene and hence potentially increase the total number of proteins in the cell (Reddy et al., 2013; Wang et al., 2022). Further visualization of the spliceosome pathway confirmed that several components were highly enriched among the module IX DEGs ( Figure 6A ). Twelve genes associated with the spliceosome were then isolated from the total DEGs ( Figure 6B ), of which 11 were upregulated only in fruits stored at 5°C for 14 d with or without 1-MCP treatment while 1 gene was downregulated. This included SlSmEb (Solyc06g072280) whose ortholog in Arabidopsis has been associated with proper alternative splicing events during chilling stress responses (Wang et al., 2022). By qPCR analysis, we confirmed that SlSmEb transcripts accumulate in both “Micro-Tom” and “Moneymaker” tomatoes during storage at 5°C irrespective of 1-MCP treatments, but little or insignificant changes in expression occur at 20°C ( Figure 6C ). The cold-specific upregulation of SlSmEb is further supported by the swift drop in expression levels (within 1 d in “Micro-Tom” tomatoes) when fruits were transferred from 5°C to 20°C.

As is the case with other biotic and abiotic stresses, the physiological and biochemical alterations which occur in plants under cold stress are controlled by several regulatory hubs that include various transcription factors, hormones, chaperones and signaling peptides (Kosová et al., 2018; Mehrotra et al., 2020). We therefore scrutinized our DEGs particularly in module IX to see how genes associated with these regulatory hubs are affected by cold stress. This led to the identification of 99 DEGs of which 60 encoded various transcription factor families, 14 were associated with jasmonate and auxin signaling, 10 encoded heat shock proteins, 13 were related with calcium signaling and 2 encoded clavata3/embryo surrounding region-related (CLE) signaling peptides ( Supplementary Table 11 ). qPCR analyses were then performed on a selected nine genes to validate the RNA-Seq data and to evaluate their expression patterns in response to cold stress in “Micro-Tom” and “Moneymaker” tomatoes. In both cultivars, transcripts of SlERF13, SlNAC22, SlHSP, SlMADS43, SlJAZ2, SlBEL3 and SlCLE2 increased dramatically, while those of SlPIF3 decreased during storage at 5°C regardless of 1-MCP treatment; changes at 20°C were insignificant ( Figure 7 ). Transcript levels in fruits at 5°C also dropped following transfer to 20°C. However, while SlCBF2 transcripts increased during storage at both 5 and 20°C in “Micro-Tom” fruits ( Figure 7A ), increased expression levels were registered only at 5°C in “Moneymaker” fruits especially after 14 d ( Figure 7B ).

The finding that ‘NcRNA metabolism’ was one of the highly enriched GO terms among module IX DEGs ( Figure 2C ) prompted us to also look into the expression patterns of genes associated with epigenetic modifications. Out of the 20 genes that we identified, 12 were associated with DNA methylation while 8 were associated with chromatin formation/remodeling ( Figure 8A ). Among the DNA methylation-related genes, 6 encoded short interfering RNA-producing dicer-like enzymes of which 4 were upregulated while 2 were downregulated only during storage at 5°C for 14 d regardless of 1-MCP treatments. Interestingly, we also found 2 genes that encode RNA-directed DNA methylation 1 (RDM1) protein; SlRDM1-like 1 (Solyc10g045190) was slightly downregulated while SlRDM1-like 2 (Solyc09g082480) was highly upregulated in both 1-MCP-treated and non-treated fruits after 14 d at 5°C. For the genes associated with chromatin formation/remodeling, 5 were upregulated while only 3 were downregulated in both 1-MCP-treated and non-treated fruits after 14 d at 5°C. All of these genes (both DNA methylation- and chromatin formation/remodeling-associated) showed insignificant expression changes during storage at 20°C and the expression changes induced at 5°C were reversed when fruit were rewarmed to 20°C.

We then selected three genes which showed high and clear expression changes based on their TPM values, that is SlRDM1-like 2, dicer-like 2d (SlDCL2d, Solyc11g008530) and chromodomain helicase 9 (CHD9)-like (SlCHD9-like, Solyc04g074680), for further analysis by qPCR. In addition, we also examined the expression of Demeter-like 2 (SlDML2, Solyc10g083630) which is associated with active DNA demethylation in tomatoes (Lang et al., 2017). In both “Micro-Tom” and “Moneymaker” tomatoes, SlDML2 transcripts increased significantly during storage at 20°C and peaked after 7 d, but this increase was inhibited in 1-MCP-treated fruits ( Figures 8B, C ). However, the expression of SlDML2 did not change significantly during storage at 5°C, but increases were observed only after rewarming at 20°C. Transcripts of SlRDM1-like 2, SlDCL2d, and SlCHD9-like showed little or no significant changes during storage at 20°C, but they accumulated highly during storage at 5°C irrespective of 1-MCP treatments ( Figures 8B, C ). The increased expression registered at 5°C was reversed within 1 d (in “Micro-Tom” fruits) and after 7 d (in “Moneymaker” fruits) of rewarming at 20°C.

After examining gene expression patterns, we determined the involvement of DNA methylation in cold response and CI development in “Moneymaker” fruits using 5-azacytidine, an inhibitor of DNA methylation (Huang et al., 2019; Zhu et al., 2020) and McrBC, a restriction enzyme that cuts methylated DNA. Fruits injected with 5-azacytidine developed relatively severe symptoms (dark spots on the skin and darkening of columella) compared to those injected with distilled water (mock) ( Figure 9A ). Furthermore, fruits injected with 5-azacytidine were generally firmer at the end of cold storage than mock fruits ( Supplementary Figure 2 ). We then assayed the DNA methylation levels of SlDML2 and SlRDM1-like 2 using McrBC-qPCR analysis. Genomic DNA was first extracted from the pericarp of fruits stored at either 5°C or 20°C for 7 d and then digested with McrBC before performing qPCR using both digested and undigested DNA as templates. Results revealed that methylation levels of both SlDML2 and SlRDM1-like 2 were higher in samples at 5°C than at 20°C, and treatment with 5-azacytidine significantly lowered the methylation levels at 5°C with respect to mock samples ( Figure 9B ). To further investigate the influence of 5-azacytidine treatment on gene expression, we examined the transcript levels of 4 genes induced by low temperature alone which we previously identified by RNA-Seq analysis. As shown in Figure 9C , 5-azacytidine treatment decreased the expression of SlRDM1-like 2, SlCBF2, and SlAOX1a only after 7 d. However, it is intriguing that the inhibitor completely decreased the expression of SlERF13 throughout the 21 d of storage at 5°C ( Figure 9C ). These results show that DNA methylation has an important role in cold responses and CI development in tomato fruits.