RNA-Seq analysis of genes affected by Cyclophilin A/DIAGEOTROPICA (DGT) in tomato root development

Cyclophilin A/DIAGEOTROPICA (DGT) has been linked to auxin-regulated development in tomato and appears to affect multiple developmental pathways. Loss of DGT function results in a pleiotropic phenotype that is strongest in the roots, including shortened roots with no lateral branching. Here, we present an RNA-Seq dataset comparing the gene expression profiles of wildtype (‘Ailsa Craig’) and dgt tissues from three spatially separated developmental stages of the tomato root tip, with three replicates for each tissue and genotype. We also identify differentially expressed genes, provide an initial comparison of genes affected in each genotype and tissue, and provide the pipeline used to analyze the data. Further analysis of this dataset can be used to gain insight into the effects of DGT on various root developmental pathways in tomato.


Introduction
The tomato (Solanum lycopersicum) cyclophilin DIAGEO-TROPICA (DGT) has been linked to auxin-regulated development through identification of the gene affected by the diageotropica (dgt) mutation 1,2 . Tomato dgt mutants are auxinresistant and display a pleiotropic phenotype that includes slow gravitropic response, lack of lateral root initiation, altered vascular development, reduced ethylene production in response to auxin, reduced apical dominance, impaired shoot and root growth, reduced fertility, and impaired fruit growth [3][4][5][6][7][8][9][10][11] . DGT likely interacts with auxin transport and signaling in a complex manner. For instance, DGT has been shown to negatively regulate auxin efflux via PIN-FORMED (PIN) transporters by altering subcellular localization of PINs and expression of some PIN genes 9 . DGT, in turn, is downregulated by auxin at the root tip 8 , suggesting functional feedback between DGT, auxin, and PINs. Through targeted gene expression quantification with RT-PCR and northern blots, DGT has been demonstrated to affect expression of a number of other auxin-related genes in addition to PIN genes. The dgt mutation reduces auxin-induced expression of genes encoding certain 1-aminocyclopropane-1-carboxylic acid synthases (ACCSs; key ethylene biosynthesis regulatory enzymes), SMALL AUXIN UPREGULATED RNA (SAUR) genes, and several members of the auxin-regulated Aux/IAA gene family in a tissue-and developmental stage-specific manner 4,6,10-14 . However, full transcriptomic profiling of tomato dgt mutants in different developmental zones has not been performed. Given the complex role of DGT in auxin pathways, an exploration into the widespread effects of DGT on the transcriptome may provide valuable insights into its potentially extensive and multifaceted role in development.
In this study, we perform a global analysis of gene expression in dgt roots that compares the root meristem, elongation zone, and differentiation zone in wildtype ('Ailsa Craig') and dgt tomato plants. The root tip provides an excellent system for studying development-related plant gene expression because cell division, elongation, and maturation are not only temporally but also spatially separated in this growth region; this allows for anatomical dissection and analysis of specific developmental zones, such as the meristem, elongation zone, and differentiation zone. The root tip is also most appropriate for this study because the dgt phenotype is the strongest in the root tip and has been characterized morphologically in the root tip 7,8 . We have previously performed histological analyses of tomato root tips including the meristem, elongation zone, and maturation zone, and demonstrated a decrease in length and number of cells of the dgt meristem and elongation zone 8 , whereas the initiation of lateral root primordia was abolished in the dgt root maturation zone [7][8][9] . Here, we present an RNA-Seq dataset containing raw reads and abundance estimates for three replicates in each zone and genotype, the pipeline used for analysis, and an initial exploration of expressed and differentially expressed genes in each developmental region that can be used to guide future investigations.

Materials and methods
Plant material and growth conditions Three biological replicates for each tissue and genotype were performed. Seven-day old tomato (Solanum lycopersicum) wildtype (WT) and dgt 1-1 plants in the 'Ailsa Craig' background 7 were used. Seeds were sterilized in 20% commercial bleach for 30 min and rinsed four times for 10 min with sterile water. Sterilized seeds were vernalized at 4°C for two days to ensure even germination and then planted on media containing 0.2× Murashige and Skoog basal medium with vitamins (PhytoTechnology), 1% sucrose, 10 mM MES buffer pH 5.7, and 0.8% agar. Seeds were germinated in Magenta boxes (16 seeds per box) in a growth chamber at 21°C under long day (16h light, 8h dark) conditions and light intensity as in Ivanchenko et al. (2013). Root samples were dissected under DIC optics at 4x objective as in Ivanchenko et al. (2006). On average, 50-100 root portions were collected per biological replicate. The meristem was dissected between the root tip and the root portion where tissue becomes more transparent. The elongation zone was collected from the proximal meristem border and the first hair bulge, and approximately 1 cm portions were collected above the first hair bulge for the differentiation zone.

RNA isolation and RNA-Seq
Tissue samples were collected in Plant RNA Reagent (Life Technologies) on ice, and total RNA was prepared using the RNeasy Mini kit (Qiagen) according to manufacturer's recommendations. The RNA pellets were dried in 1.7 mL centrifuge tubes and solubilized in 178 µL 1X RNA Secure Reagent (Ambion) preheated at 65°C, and incubated at 65°C for 10 min, mixing by pipetting a few times. Then 20 µL 10X DNase I buffer and 2 µL RNase-free DNase I (Ambion) were added to each tube, and tubes incubated for 10 min at 37°C. 700 µL RLT buffer was added (to which 7 µL 2-Mercaptoethanol was freshly added) to each sample, which were then mixed by vortexing. 500 µL ethanol was then added and samples were mixed again by vortexing. Each sample (2 X 700 µl) was applied to an RNeasy Mini spin column from the RNeasy kit and RNA cleanup performed following the manufacturer's instructions. Each sample was eluted with 30 µL nuclease-free ultrapure H 2 O, and the RNA concentrations were measure using a NanoDrop 1000 spectrophotometer (Thermo Fisher).
Illumina sequencing of RNA RNA-Seq library preparation and sequencing were performed at the Oregon State University Center for Genome Research and Biocomputing. Libraries were prepared using the TruSeq RNA Library Prep Kit v2 (Illumina) and sequenced as single-end 51 bp reads on the Illumina HiSeq 2000 using a total of two lanes.

Gene alignment and expression level analysis
Sequencing produced 5-8 FASTQ files per replicate, which were merged into a single FASTQ file per replicate. Reference transcripts were extracted and preprocessed from the NCBI Heinz 1706 genome sequence 15 using RSEM's 16 (RNA-Seq by Expectation Maximization, version 1.3.1) rsem-prepare-reference function. Using RSEM, raw sequence reads were then aligned to the reference transcript sequences and estimated transcripts per million (TPM) values were calculated using default parameters. Resulting gene-level estimates from biological replicates were merged into a single input matrix and EBSeq 17 (version 1.26.0) was then used to test for differential expression between dgt and WT for each root-tip zone. Raw RNA-Seq reads, abundance estimates, and differential expression analysis are publicly available in the Sequence Read Archive (SRA) and Gene Expression Omnibus (GEO), see Data availability.

Data visualization
As a preliminary analysis to guide future investigations, a basic comparison of presence/absence calls and differential gene expression is presented. Expression analysis was performed using TPM and differential expression analysis was performed using the posterior probability that the gene is differentially expressed (PPDE) and posterior fold change (postFC). Genes with an average of TPM > 2 across biological replicates were compared between zones and genotypes ( Figure 1; data summarized in Data File 1 18 ). Following EBSeq analysis, genes were filtered for PPDE = 1 to identify those which were differentially expressed.
Of the differentially expressed genes, those with postFC (dgt over WT) > 2 were considered upregulated in dgt and those with postFC < 0.5 were considered downregulated in dgt. Upregulated and downregulated genes were compared between zones and genotypes (Figure 2; data summarized in Data File 2 19 ). Principal component analysis (PCA) of the data was performed using TPM values of each gene and was calculated using the scikit-learn 20 PCA function with default parameters (Figure 3).

Conclusion
It is clear from the outcomes in Figure 1 and Figure 2 that reduced function of DGT has sweeping effects on the transcriptome in all three developmental zones examined in this study, supporting the concept that DGT very likely plays important and potentially   complex roles in multiple developmental pathways. Additionally, PCA using TPM values for each sample demonstrated that while there was some variance within the replicate pools, replicates from different genotypes were distinctly separate from each other (Figure 3). Further functional genomics studies are needed to narrow down the most likely direct interactions with DGT, leading to the identification of specific functional roles. We hope that this dataset will be of value to the community in future studies in this area.

Aurélien Bailly
Department of Plant and Microbial Biology, University of Zurich, Zurich, Switzerland Ivanchenko and colleagues describe here a transcriptomics data set profiling gene expression in root tips of the tomato mutant dgt and matching wildtype. As for several plant immunophilins, loss of function of DIAGETROPICA triggers dramatic structural and behavioural phenotypes in tomato. The link between immunophilins and hormonal control has been extensively recognised; however, the detailed molecular mechanisms behind them remain elusive.
The interest in this RNAseq dataset lies within the dissection of the root into the three standard spatio-temporal zones of development. Gaining insight into the effect of the dgt mutation in the transcriptomes of the meristematic, elongation and mature zones of the root would potentially lead to the isolation of long-sought components of root growth and development machineries.
The work presented here is technically sound and used state of the art transcripts-reading methods. Data deposition was correctly performed and gene lists produced after alignments and differential expression analysis through regular scripts are clear and readily exploitable. Although the level of variation between replicates, as mentioned by authors, do vary in the different root zones samples, I trust that they separate enough in PCA to yield interesting outcomes. In addition, the relatively reduced number of differentially expressed genes filtered in this analysis will allow for a rapid identification of proteins associated with the dgt phenotype, perhaps outside of the classical pool of auxin-related genes.
I would express only two regrets: 1) the DGT protein was shown to hold a conserved auxin signalling function down to moss, but it is not mentioned here; and 2) other immunophilins present overlapping roles with DGT, a paradigm worth mentioning somewhere in the introduction.

Is the rationale for creating the dataset(s) clearly described? Yes
Are the protocols appropriate and is the work technically sound?