5.1 Plant material
This study was performed at La Platina Research Center that belongs to the Instituto de Investigaciones Agropecuarias (INIA), located in Santiago, Chile (33°34′20′′S; 70°37′32′′W; 630 m.a.s.l.), during the 2015-2016, 2016-2017 and 2017-2018 seasons. The plant material (plants of 5-10 years old grown on their own roots, conducted with the Spanish trellis system and managed under standard irrigation, fertilization and pest management programs) of two seedless table grape varieties, cv. Thompson Seedless and L23 (‘Line #23’: F1 from cv. Ruby Seedless x cv. Centennial Seedless). Six individual plants from each genotype were used for qPCR experiments; four biological replicates were used for RNAseq experiments. Each plant had treated (GA3) and untreated (Control) clusters to reduce any bias related to individual variation.
Timing of samples was based on the Einchorn-Lorenz phenology system adapted for Vitis vinifera L. by Coombe, 1995 [42]. Stage(s) evaluated in each experiment are indicated in each figure.
Thompson Seedless (also known as ‘Sultanina’) is a free, ancestral variety, supposedly native to Iran; was introduced into Chile by the beginning of XX century from California. L23 is a segregant belonging to INIA’s breeding program. These parental varieties were both originally released by UC Davis, CA, but are currently free. As ‘Thompson seedless’, they were also introduced into Chile many years ago and kept through vegetative propagation.
5.2 Growth regulator treatment
Clusters of cv. Thompson Seedless and L23 were subjected to the following growth regulator treatments: (1) Thompson Seedless were sprayed five times: once with 10 ppm gibberellic acid (GA3, Pro-Gibb 40% Valent Biosciences Chile S.A., Santiago, Chile) during pre-bloom for bunch elongation; 15 ppm at full bloom for thinning and 200 ppm for berry enlargement (three doses of 40, 100 and 60 ppm with seven-day intervals) when berries reached 4~5 mm diameter. (2) L23 were sprayed with 120 ppm GA3 (3 x 40 ppm, intervals of seven days) for berry enlargement when berries reached 4~5 mm diameter. Dosage for cv. Thompson Seedless and L23 was defined in previous assays, optimizing berry enlargement at the lowest possible value of berry drop, as has been described in [11].
5.3 Physiological parameters
The parameters measured were berry and pedicel diameter (mm), soluble solids (% w/w g sucrose per 100 g solution), titratable acidity (g L-1 tartaric acid) and firmness (g mm-1) of 10 bunches collected from different plants. For each variable, 30 healthy and homogenous berries with their cap stems were randomly sampled per cluster. Diameter was measured with a digital caliper. Soluble solid content of fruit was measured with a refractometer (ATC-1E, Atago, Tokyo, Japan). Firmness of berries (considering skin and flesh) was determined using a firmness tester (Firmtech II, BioWorks, Wamego, KS). Titratable acidity was measured in the pooled juice of 10 berries per cluster by titration with 0.1 N NaOH (pH 8.2).
Random samples of pedicels from clusters were pooled (300 ~ 500 pedicels). Plant material was immediately frozen with liquid nitrogen and stored at -80 ºC for molecular analysis.
5.4 Berry drop
Clusters were harvested based on the average soluble solids content (18° Brix). For postharvest storage experiments, 10 clusters per treatment (control and GA3-treated) were stored after harvest for 15, 30 and 45 days at 0 ºC and then berry drop was measured. To determine percentage of berry drop (shattering), 10 bunches were weighed before and after being shaken for 30 s; weight of detached berries was recorded. Berry drop percentage was calculated as BD% = [(DB + SB)/TB] ⁄ 100, where DB represents the weight of dropped berries during storage time at 0 °C, SB represents the weight of berries dropped after shaking of bunches and TB represents the total weight of the bunch.
5.5 Determination of lignin
For lignin estimation, quantitative determination by the acetyl-bromide method [21] was followed with slight modifications. Pedicels were dried at 105 °C for 16 hours in an oven and cooled in a vacuum desiccator until the next step. Dry samples (0.3 g) were homogenized in 50 mM potassium phosphate buffer (5 mL, pH 7.0) using an Ultra Turrax T25 homogenizer (IKA-Werke GmbH & Co., Staufen, Germany). The solution was centrifuged (1400xg, 5 min) and the pellet was washed by successive stirring and centrifugation as follows: twice with phosphate buffer (pH 7.0; 7 mL), three times with 1% (v/v) Triton X-100 in pH 50mM potassium phosphate 7.0 buffer (7 mL), and six times with acetone (5 mL). All supernatants were monitored by spectrophotometric measurements at 280 nm to ensure no contamination with protein and UV-absorbing materials in downstream steps. The pellet was dried in an oven (60 °C, 24 hours) and kept cool in a vacuum desiccator until further reaction with acetyl bromide. Samples and a blank solution were prepared as follows: dry matter obtained from the previous step (considered as a protein-free cell wall fraction) was placed in a screw cap 15 mL centrifuge tube (20 mg per sample) containing 0.5 mL 25% acetyl bromide (v/v in glacial acetic acid) and incubated at 70 °C for 30 min. After digestion, samples were cooled in an ice bath and mixed with 0.9 mL of 2 M NaOH and 0.1 mL of 5 M hydroxylamine-HCl. Samples and blank solution were solubilized in 8 mL of glacial acetic acid and absorbance of supernatant was measured at 280 nm immediately after centrifugation (1400xg, 5 mins). Dilution in glacial acetic acid was performed to avoid signal saturation. The standard curve was constructed using alkali lignin (Aldrich 37, 096-7) with a reported absorptivity value of 23.03 g-1 L cm-1 (Additional file 1: Fig. S1B). Results were expressed as mg lignin g-1 cell wall as appears in the cited protocol [21].
5.6 Dry Matter
For the determination of changes in pedicel dry matter accumulation elicited by GA3 treatments, samples from different conditions at harvest were measured in an analytical balance. Samples taken from bunches when berries reached an average of 18 °Brix were pooled into a group of 30~60 pedicels and weighed for fresh weight determination (at least 0.3 g fresh weight to avoid technical error). Then replicates were dried at 105 °C in an oven for 16 hours and dry weight was measured [21]. Percentage of dry matter and weight of pedicels was determined, considering the number of pedicels weighed in each batch.
5.7 RNA isolation
Total RNA was isolated from 1.0 g frozen tissue using a modified hot borate method [43], following all the indications listed in the cited protocol with the only exception that sample tissue was reduced to a third of the recommended amount and polyvinylpyrrolidone was doubled to avoid phenolic compound contamination before RNA isolation. The quantity and quality of RNA were assessed with Qubit® 2.0 fluorometer (InvitrogenTM by Life Technologies, Singapore). Spectrophotometric determination of A260/280 and A260/230 ratios and electrophoresis in 1.2% formaldehyde-agarose gels verified the quality and integrity of extracted RNA.
5.8 Library synthesis and sequencing
For the RNAseq experiment, pooled pedicel samples of individual plants were obtained seven days after the last GA3 application (corresponding to EL31 stage according to the Einhoch-Lorenz phenology system [42], plus seven days—defined in this study as 7DAT). Two levels were used in the treatment factor, control and GA; and two levels were considered in the genotype factor, cv. Thompson Seedless and L23. Four biological replicates were used for each condition, summing to 16 libraries which were synthesized and sequenced.
Prior to library synthesis, total RNA was assessed by fragment analyzer PROSize® 2.0 version 1.3.1.1 (Advanced Analytical Technologies, Inc., Ames, IA, USA). Mean RNA quality in all samples was 7.35 (SD: 0.56), confirming the integrity of extractions (refer to Table SIV for individual scores). Then 2.5 μg aliquots were used to isolate poly(A) mRNA for preparation of libraries using TruSeq RNA Sample Prep Kit v2, following the manufacturer’s instructions described in the TruSeq RNA Sample Preparation v2 Guide, Part #15026495 Rev. F (Illumina, Inc.). Libraries were sequenced using the HiSeq 4000 platform (Illumina). Libraries were sequenced as paired-end data (2 x 100 bp).
5.9 Pre-processing of reads
Low-quality reads (Phred score Q<25, nucleotides with undefined base assignment N>1 and read length<50) were removed using the wrapper tool Trim Galore! [44]. Contamination with Illumina adapters was handled with the same tool, removing sequences matching the adapters.
5.10 Mapping and matrix count extraction
Alignment of clean reads to the V. vinifera L. PN40024 12X reference genome [45] was performed using STAR software, considering paired-end data default parameters [46]. Uniquely mapped reads were kept for generation of a count matrix using HTseq [47]. The gene model was extracted from the CRIBI database. Genes reported were based on Genoscope structural annotation. Summarized mapping results are given in Table S4.
Data from the count matrix of uniquely mapped reads were visualized to confirm major effects expected from the experimental design. Assessment using multidimensional scaling and hierarchical clustering based on counts per million reads showed no anomalies; differences were given mainly by genotype and treatment factors, with low variability of replicates within groups (for experimental details please see Fig. S10).
5.11 Differential expression analysis
Count matrices were analyzed with EdgeR package [48] under R 3.4.4 software. Correction for composition bias in each sample was handled with the CalcNormFactor option. Then limma-voom transformation of data with was performed with the limma package [25]. Since the number of DEGs was relatively high, an additional filter of log2-fold-change was established. To keep a low false discovery rate, the method described in [22] was followed, considering a prior log2-fold-change of 1. For differential expression analysis considered an adjusted p-value (FDR) < 0.05.
5.12 Gene set enrichment analysis and Revigo
Gene set enrichment analysis was conducted on the Agrigo platform v2.0 [49] using hypergeometric distribution of data; p-values were adjusted by the Benjamini-Yeukiteli method and significant terms (FDR<0.05) were analyzed with the Revigo algorithm to reduce redundancy of terms based on semantic relatedness. The cutoff value used was C: 0.7 [23].
5.13 MapMan analysis
Mapping of differentially expressed genes as determined in 2.8 was performed using MapMan software (X4 version). Annotation of V. vinifera L. was downloaded from BioMart of the Phytozome v12.1 resource page; amino acid sequences were fed into the Mercator 4 online tool, creating an up-to-date mapping file. The Wilcoxon test with adjusted p-value by the Benjamini-Yeukiteli method was used to detect significant differences among BIN categories provided by MapMan [24].
5.14 cDNA synthesis
Synthesis of cDNA was performed by reverse transcription with 2 µg of total RNA as template using MMLV-RT reverse transcriptase and oligo dT primers (Promega, Madison, WI), according to standard procedures. Concentration of cDNA was assessed by measuring the absorbance at 260 nm.
5.15 Real-time qPCR assays
Transcript abundance was analyzed by qPCR with the LightCycler® 96 system from Roche (Roche Diagnostics, Mannheim, Germany), using SYBR Green to measure the amplified RNA-derived DNA products as described in [50]. Gene-specific primers (Table IV) were designed using the online software Primer3plus [51] and were synthesized by IDT (Integrated DNA Technologies). The qPCR assays were performed on six biological samples in duplicate.
The reference gene was GSVIVG01011810001 – probable fructose-bisphosphate aldolase 3 chloroplastic (AFLC3). This gene was identified based on RNAseq results; all non-differentially expressed genes were filtered and ranked according to the lowest variance and coefficient of variation, considering the 16 libraries. Five genes were selected, and primers were synthetized (Additional file 13: Table S5). AFLC3 gene values were used in downstream calculations related to relative expression.
5.16 Experimental design and statistical analyses
A fully randomized experimental design was used. The two main factors considered were Genotype factor (two levels: cv. Thompson Seedless and L23) and Treatment factor (two levels: Control and GA3).
Data were previously analyzed with Levene test to verify homogeneity of variance assumption (p>0.05). Then, Analysis of Variance test (ANOVA) and Tukey’s Honest Significance Difference test (Tukey HSD) were conducted to find significant differences between conditions (p<0.05). These analyses were applied to find significant differences on Fig. 1B (analysis was performed on each season, separately. Since that ‘season’ effect could not be fitted to a linear model), Fig. 2, Fig.7 and Table 1. The number of observations (n) can be found on the corresponding description of each figure.
To extract differential expressed genes from RNAseq experiments, data was contrasted against a negative bimomial data distribution through EdgeR package [48], obtained p-values were applied with the Benjamini-Hochberg method to control false discovery ratio (FDR). As was explained before, since the number of differentially expressed genes (DEGs) was high. The method described by [22] on limma package [25] was applied to consider additional log2 fold-change (FC) of the treated over control samples between each genotype with a prior log2 FC = 1, p-values were corrected according to Benjamini-Hochberg method and a final FDR<0.05 was considered to identify DEGs. Please refer to Additional file 3 to see results regarding this set of statistical analyses.
For gene set enrichment analysis on gene ontology (GO) terms, hypergeometric distribution of data was assumed, and calculated p-values were applied with Benjamini-Yekuiteli method to control false discovery at a rate of 0.05 (FDR<0.05). From here, revigo analysis was conducted [23] (with distance between GO terms based on relative similiraty index). Cutoff value was considered 0.7 in cv. Thompson Seedless and 0.5 in L23.
For mapman [24] results, Wilcoxon test was conducted on assigned categories to detect significant differences p-values were corrected by Benjamini-Hochberg method (FDR<0.05).
Hierarchal clustering method on heatmap figure (Figure 6) was accomplished by calculating Euclidean distance matrix among genes (rows) and among replicates (columns), the quantitative variable was the fragment per kilobase of exon per million reads (FPKM). This quantity was normalized to Z-score for the purpose of comparing the abundance on each cell generated by this method.