Identification of genes associated with the regulation of cold tolerance and the RNA movement in the grafted apple

Abstract

In grafted apple, rootstock-derived signals influence scion cold tolerance by initiating physiological changes to survive over the winter. To understand the underlying molecular interactions between scion and rootstock responsive to cold, we developed transcriptomics and metabolomics data in the stems of two scion/rootstock combinations, 'Gala'/'G202' (cold resistant rootstock) and 'Gala'/'M9' (cold susceptible rootstock). Outer layers of scion and rootstock stem, including vascular tissues, were collected from the field-grown grafted apple during the winter. The clustering of differentially expressed genes (DEGs) and gene ontology enrichment indicated distinct expression dynamics in the two graft combinations, which supports the dependency of scion cold tolerance on the rootstock genotypes. We identified 544 potentially mobile mRNAs of DEGs showing highly-correlated seasonal dynamics between scion and rootstock. The mobility of a subset of 544 mRNAs was validated by translocated genome-wide variants and the measurements of selected RNA mobility in tobacco and Arabidopsis. We detected orthologous genes of potentially mobile mRNAs in Arabidopsis thaliana, which belong to cold regulatory networks with RNA mobility. Together, our study provides a comprehensive insight into gene interactions and signal exchange between scion and rootstock responsive to cold. This will serve for future research to enhance cold tolerance of grafted tree crops.

Figure 1

RNA-seq experiment design and DEG analysis. (a) A scheme for cold tolerance signaling study. Stem vasculature tissues of ‘Gala’ apple scion grafted onto two contrasting rootstock cultivars showing a different level of cold tolerance (‘G202’, ‘M9’) were sampled at three winter stages: CA (early winter/cold acclimation, 6th June), DW (deep winter, 17th July), DA (late winter/cold de-acclimation, 26th August). (b,c) Principal component analysis of RNA-seq expression data. Top 5000 genes listed according to the variance rank were chosen as informative genes and used for PCA analysis: (b) rootstock (RS) and (c) scion (SC). (d) Venn diagram of DEGs between scion and rootstock with a criterion of |fold change|≥ 1.5 FDR ≤ 0.05. (e) RNA-seq validation by qRT-PCR. A high significant correlation (r ≥ 0.9) was made between RNA-seq result (RPKM) and qPCR.

Figure 2

Expression dynamics of DEG clusters, functional annotation and comparative analysis of DEGs to Arabidopsis cold regulatory network. (a) K-means clustering of rootstock 14,747 DEGs (k = 6). (b) Gene ontology enrichment (biological process) of DEGs from five rootstock clusters. All clusters except CR5 (All high in ‘M9’) showed a significant GO terms with a criteria of FDR 0.05. (c) K-means clustering of scion 1261 DEGs (k = 5). (d) Gene ontology enrichment (biological process) of scion DEGs from with a criteria of p < 0.05. (e–g) Comparative analysis of scion DEGs to the cold regulatory network data conserved in Arabidopsis. (e) Venn diagram between scion orthologs and genes involved in cold regulatory network. 830 non-duplicated putative orthologs of Arabidopsis were obtained from 1251 scion DEGs and used for cross-species comparison. 72 scion DEGs were commonly found consisting of 11 transcription factors (TFs) and 61 targets (TGs). (f,g) Two sub-networks were constructed from the 72 scion DEGs found in cold regulatory network in Arabidopsis. The expression data at CA stage was used for visualization, as the reference network data was constructed from the expression profile of Arabidopsis exposed to non-freeze cold temperature, which conditions were similar to those of CA stage in our RNA-seq data. Expression ratios in scions of genes in these sub-networks (‘Gala’/‘G202’ over ‘Gala’/‘M9’ at CA undergoing cold acclimation) are color-coded. CA cold acclimation, DW deep winter, DA cold de-acclimation, G202 ‘Gala’/‘G202’, M9 ‘Gala’/‘M9’.

Figure 2

Expression dynamics of DEG clusters, functional annotation and comparative analysis of DEGs to Arabidopsis cold regulatory network. (a) K-means clustering of rootstock 14,747 DEGs (k = 6). (b) Gene ontology enrichment (biological process) of DEGs from five rootstock clusters. All clusters except CR5 (All high in ‘M9’) showed a significant GO terms with a criteria of FDR 0.05. (c) K-means clustering of scion 1261 DEGs (k = 5). (d) Gene ontology enrichment (biological process) of scion DEGs from with a criteria of p < 0.05. (e–g) Comparative analysis of scion DEGs to the cold regulatory network data conserved in Arabidopsis. (e) Venn diagram between scion orthologs and genes involved in cold regulatory network. 830 non-duplicated putative orthologs of Arabidopsis were obtained from 1251 scion DEGs and used for cross-species comparison. 72 scion DEGs were commonly found consisting of 11 transcription factors (TFs) and 61 targets (TGs). (f,g) Two sub-networks were constructed from the 72 scion DEGs found in cold regulatory network in Arabidopsis. The expression data at CA stage was used for visualization, as the reference network data was constructed from the expression profile of Arabidopsis exposed to non-freeze cold temperature, which conditions were similar to those of CA stage in our RNA-seq data. Expression ratios in scions of genes in these sub-networks (‘Gala’/‘G202’ over ‘Gala’/‘M9’ at CA undergoing cold acclimation) are color-coded. CA cold acclimation, DW deep winter, DA cold de-acclimation, G202 ‘Gala’/‘G202’, M9 ‘Gala’/‘M9’.

Figure 3

Investigation of highly-correlated genes for seasonal flow. (a) A scheme for the detection of mobile mRNAs of DEGs with seasonal flow. (b,c) Pearson correlation analysis on 691 DEG set. (b) Highly-correlated DEGs (r ≥ 0.7) between rootstock and scion were selected for further investigation. (c) Seasonal flow DEGs. Highly-correlated DEGs were divided into six groups by two flow directions: upwardly mobile (rootstock to scion) or downwardly mobile (scion to rootstock) combined with the cold stages. (d) Venn diagram of up- and downwardly mobile mRNAs of DEGs. (e) Functional categorization of highly-correlated 544 DEG set based on their DNA sequence homology using Mercator v3.6 tool. CA cold acclimation, DW deep winter, DA cold de-acclimation.

Figure 4

SNP detection and functional assays to validate the mRNA mobility associated with cold. (a–d) SNP detection to screen the mobility of rootstock-derived mRNAs to ‘Gala’ scion associated with the seasonal cold tolerance. (a) Estimation of heterozygosity of ‘Gala’, ‘G202’, ‘M9’ genotypes based on k-mer distribution. Heterozygous and homozygous peaks are shown from resequencing reads from paired-end libraries. (b) Genome-wide SNP calling from ‘Gala’, ‘G202’, ‘M9’ genotypes using resequencing data. (c) SNP detection of 544 highly-correlated seasonal flow DEGs. SNPs unique to each genotype were detected with the parameters of coverage 10, count 10, and frequency 20% in each cold stage. (d) Expression heatmap of 31 mobile mRNAs that contain rootstock-derived SNPs pooled from all winter stages. (e–h) Tobacco agroinfiltration assay. (e,f) Expression profile of genes in four tissue types of transiently overexpressing tobacco including the infiltrated leaf, petiole and the adjacent stem above/below the infiltrated leaf. Tissues were sampled at three cold stages (early, deep, late) following three days after agroinfiltration. (e) MdTSJT1; (f) eGFP (control). There were six to eight biological replicates. NbFBOX (Niben.v0.3.Ctg24993647) was used as reference gene for qRT-PCR calculation. (g) Degree of cold-dependent mRNA mobility in tobacco stem. (h) Comparison of degree of cold-dependent mobility of MdTSJT1 between qRT-PCR of transiently overexpressing tobacco and RPKMs of RNA-seq data (upward > 0, downward < 0). (i,j) Arabidopsis grafting assay and the investigation of cold-dependent AtVNI2 mRNA movement. (i) A scheme for Arabidopsis grafting for cold treatment. Six-day-old seedlings of wild-type or mutant medium-grown Arabidopsis plants were grafted. After post-grafting recovery, twenty-day-old grafted seedlings (14 days after grafting) were cold-treated and the whole scion/rootstock tissues above/below graft junction were collected to explore mRNA movement assay. (j) The relative expression of AtVNI2. Lines on the bars show standard error of the mean. AtACTIN8 was used as reference gene for qRT-PCR calculation. Asterisk indicates significant difference at p < 0.05 by Student’s t-test. As expected, there was no expression data of AtVNI2 detected in the vni2/vni2 homografted Arabidopsis. CA cold acclimation, DW deep winter, DA cold de-acclimation, RS rootstock, SC scion, G202 ‘Gala’/‘G202’, M9 ‘Gala’/‘M9’.

Figure 5

A proposed scheme for mobile mRNA signals associated with cold tolerance in stem vasculature of grafted apple. Candidate mRNA signals that potentially mediate interactive metabolic responses between rootstock and scion are shown with the movement direction. CA cold acclimation, DW deep winter, DA cold de-acclimation.

Figure 6

Metabolomics profile using UPLC-QTOF-MS. (a,b) Principal component analysis and expression pattern of metabolites were analyzed by Partial Least-Squares Discriminant Analysis (PLS-DA). (a) rootstock; (b) scion. (c-d) Venn diagram of differentially accumulated metabolites across three winter stages (CA, DW, DA) with the criteria of |fold change|≥ 1.5 and p-value ≤ 0.05. (c) rootstock; (d) scion. (e) Venn diagram of the shared metabolites between scion and rootstock tissues in apple graft combinations pooled from all stages. (f) Expression profile of eight identified metabolites from the 14 shared scion metabolites. Within the genotype columns, left is rootstock labeled as RS and right is scion labeled as SC. (g) Flavonoid biosynthesis pathway in which shared scion metabolites were enriched and potentially considered to be involved. Detected metabolites are marked with asterisks. CA cold acclimation, DW deep winter, DA cold de-acclimation, RS rootstock, SC scion, G202 ‘Gala’/‘G202’, M9 ‘Gala’/‘M9’.

Figure 6

Metabolomics profile using UPLC-QTOF-MS. (a,b) Principal component analysis and expression pattern of metabolites were analyzed by Partial Least-Squares Discriminant Analysis (PLS-DA). (a) rootstock; (b) scion. (c-d) Venn diagram of differentially accumulated metabolites across three winter stages (CA, DW, DA) with the criteria of |fold change|≥ 1.5 and p-value ≤ 0.05. (c) rootstock; (d) scion. (e) Venn diagram of the shared metabolites between scion and rootstock tissues in apple graft combinations pooled from all stages. (f) Expression profile of eight identified metabolites from the 14 shared scion metabolites. Within the genotype columns, left is rootstock labeled as RS and right is scion labeled as SC. (g) Flavonoid biosynthesis pathway in which shared scion metabolites were enriched and potentially considered to be involved. Detected metabolites are marked with asterisks. CA cold acclimation, DW deep winter, DA cold de-acclimation, RS rootstock, SC scion, G202 ‘Gala’/‘G202’, M9 ‘Gala’/‘M9’.