Transcriptomic Analysis of Oenococcus oeni SD-2a Response to Acid Shock by RNA-Seq

Oenococcus oeni can be applied to conduct malolactic fermentation (MLF), but also is the main species growing naturally in wine. Due to the high stress tolerance, it is an interesting model for investigating acid response mechanisms. In this study, the changes in the transcriptome of O.oeni SD-2a during the adaptation period have been studied. RNA-seq was introduced for the transcriptomic analysis of O. oeni samples treated with pH 4.8 and pH 3.0 at 0 and 1 h, respectively. Gene ontology (GO) and Kyoto encyclopedia of genes and genome (KEGG) were performed to compare the transcriptome data between different treatments. From GO analysis, the majority of differentially expressed genes (DEGs) (pH 3.0_1 h-VS-pH 4.8_1 h, pH 3.0_1 h-VS-pH 4.8_0 h, and pH 4.8_1 h-VS-pH 4.8_0 h) were found to be involved in the metabolic process, catalytic activity, cellular process, and binding. KEGG analysis reveals that the most functional gene categories affected by acid are membrane transport, amino acid metabolism and carbohydrate metabolism. Some genes, like the heat shock protein Hsp20, malate transporter and malate permease, were also over-expressed in response to acid stress. In addition, a considerable proportion of gene indicate a significantly different expression in this study, are novel, which needs to be investigated further. These results provide a new viewpoint and crucial resource on the acid stress response in O. oeni.


INTRODUCTION
Malolactic fermentation (MLF) is a biological process involved in winemaking, in which tarttasting dicarboxylic malic acid, naturally present in grape must, is converted to softer-tasting monocarboxylic lactic acid and carbon dioxide by decarboxylation (Spano and Massa, 2006). Through MLF, Oenococcus oeni can bring stabilization, sensory impacts, and deacidification to most red wines, so MLF and O. oeni are crucial in the process of winemaking (Wang et al., 2015). O. oeni is the main lactic acid bacteria existing in MLF. MLF and the growth of O. oeni are clearly inhibited by several of the physiochemical properties of wine (Betteridge et al., 2015). The four main stress factors in wine affecting MLF are ethanol (10-16% v/v), low pH (3.0-3.5), SO 2 (over 10 mg/L), and low temperature (can be below 12 • C) (Spano and Massa, 2006;Betteridge et al., 2015;Olguin et al., 2015;Darsonval et al., 2016). Many efforts have been put to investigate the mechanism of stress response of O.oeni (Spano and Massa, 2006;Olguin et al., 2015).
Low pH appears as a crucial parameter that limits bacterial growth in wine (Fortier et al., 2003). Currently, several studies have been launched to understand how O. oeni response under acid stress conditions, such as membrane composition and fluidity, pH homeostasis, oxidative stress response, DNA, and protein damage repair (Darsonval et al., 2016). But the mechanism of stress adaption in O. oeni still needs a further research.
The transcriptome of O. oeni has been studied and quantified via traditional approaches, like hybridization, fingerprinting, and tilling microarrays (Marques et al., 2012;Olguin et al., 2015;Margalef-Català et al., 2016). According to the transcriptomics and proteome results, the mechanism of stress response in O. oeni is believed very complicated, which involves series of proteins (GroEL, GroES, etc.), genes (dnaJ, dnaK, and hsp18, etc.) and metabolic pathways (amino acid transport and metabolism, malate, and citrate metabolism, etc.) (Margalef-Català et al., 2016). Nevertheless, there are still some disadvantages of these techniques, for example, non-specific cross-hybridization usually cause a high background level which limits the detection range, the transcripts can be detected only with high copy number, and the total coverage of the transcripts are almost unknown (Liu et al., 2015). Additionally, it is difficult and arduous to normalize methods and compare the expression data from different experiment.
RNA-seq is a revolutionary method with many advantages, like rapidness, high precision, reproducibility, and low cost (Liu et al., 2015). This technique is mainly applied to study the transcriptome differences from different treatments. The complexity, plasticity, and regulation of bacterial transcriptomes have been gradually appeared with the application of RNA-seq technology (Sorek and Cossart, 2010).
To provide genetic information on the acid response mechanisms of O. oeni, the transcriptome dataset was generated by using Illumina HiSeq TM 2500 platform. The transcriptomes of cells with and without acid stress were compared to determine the changes in the gene transcription level, as well as the functions and KEGG pathways of differentially expressed genes (DEGs) were analyzed. The RNA-seq data and the expression patterns are valuable genetic resources, that can advance knowledge on acid stress response of O.oeni or other bacteria's. Understanding the stress response mechanisms may help us to improve MLF starter robustness without using genetic engineering. Several works were done on stress mechanisms in O.oeni in wine-like medium /wine /microvinification etc. (Olguin et al., 2015;Margalef-Català et al., 2016). The originality, in this case, is the use of RNA-seq to investigate low pH response in O.oeni.

Bacterial Strain
The MLF starter used in this study is O. oeni SD-2a, which shows strong abilities to survive in stress conditions, and more active than commercial type strain (Viniflora R Oenos) in MLF ability. It was isolated from Chinese wines regions (Shandong province) and stored in College of Enology, Northwest A&F University (Liu, 2002;Wang et al., 2003;Zhang, 2008;Li et al., 2016). Many studies have been done on the commercial application of O. oeni SD-2a. The strain O. oeni SD-2a has obtained patent protection (02123444.2).

Growth Conditions
O. oeni SD-2a was cultured at 28 • C in a flask containing FMATB broth medium at pH 4.8 (glucose 5 g/L, D, L-malate 5 g/L, yeast extract 5 g/L, peptone 10 g/L, MgSO 4 •7H 2 O 0.2 g/L, MnSO 4 •4H 2 O 0.05 g/L, Cysteine/HCl 0.5 g/L, and tomato juice 250 mL) (Li et al., 2009). When cultures reached the midexponential phase (OD600 nm ≈ 1) they were mixed and divided into six equal parts. Then cells were harvested by centrifugation (12, 000 rpm for 1 min at 25 • C). Immediately, they were washed by FMATB broth medium at pH 3.0 and pH 4.8 (control) into a same sterile flask, respectively. The possible effect of centrifugation and time were evaluated using the control assay with pH 4.8. All assays were performed in triplicate using independent cultures and incubated at 28 • C. Samples were taken at time zero just before acid shock, and then at one hour with or without acid shock (Margalef-Català et al., 2016).

RNA Extraction
Cells were harvested and kept by following the protocol of (Margalef-Català et al., 2016). Total RNA was extracted by using the RNAprep pure Cell/Bacteria Kit (Tiangen, Beijing, China) following the manufacturer's instructions. To determine the concentration of RNA, the absorbance at 260 nm was measured using a BioDrop µLITE Spectrophotometer (Tamar Laboratory Supllies LTD., Cambridge, England) (Margalef-Català et al., 2016). The RNA integrity number (RIN) and 28S:18S ratio were also measured, total RNA samples with RIN > 7.0 and a 28S:18S ratio > 1.8 were used in subsequent experiments (Miller et al., 2009).

cDNA Library Construction and Sequencing
Sequence libraries were generated and sequenced by CapitalBio Technology (Beijing, China). The triplicate samples of all assays were constructed an independent library, and do the following sequencing and analysis. The NEB Next Ultra RNA Library Prep Kit for Illumina (NEB) was used to construct the libraries for sequencing. NEB Next Poly(A) mRNA Magnetic Isolation Module (NEB) kit was used to enrich the poly(A) tailed mRNA molecules from 1 µg total RNA. The mRNA was fragmented into ∼200 base pair pieces. The first-strand cDNA was synthesized from the mRNA fragments reverse transcriptase and random hexamer primers, and then the second-strand cDNA was synthesized using DNA polymerase I and RNaseH. The end of the cDNA fragment was subjected to an end repair process that included the addition of a single "A" base, followed by ligation of the adapters. Products were purified and enriched by polymerase chain reaction (PCR) to amplify the library DNA. The final libraries were quantified using KAPA Library Quantification kit (KAPA Biosystems, South Africa) and an Agilent 2100 Bioanalyzer. After quantitative reverse transcription-polymerase chain reaction (RT-qPCR) validation, libraries were subjected to paired-end sequencing with pair end 150-base pair reading length on an Illumina HiSeq sequencer (Illumina) (Kwon et al., 2016).

RNA-Seq: Data Analysis
The genome of O. oeni SD-2a was used as reference (unpublished). The sequencing quality were assessed with FastQC (Version 0.11.5) and then low quality data were filtered using NGSQC (v0.4).The clean reads were then aligned to the reference genome using HISAT2 (Johns Hopkins University, USA) with default parameters (Liu et al., 2015).
The processed reads from each sample were aligned using HISAT (Johns Hopkins University, USA) against the corresponding O. oeni SD-2a reference genome. The gene expression analyses were performed with Cuffquant and Cuffnorm (Cufflinks 2.2.1).
Cuffdiff was used to analyze the DEGs between samples. The standardization method of Cuffdiff is geometric, with the per-condition and pooled as the discrete model (Trapnell et al., 2013). Thousands of independent statistical hypothesis testing were conducted on DEGs, separately. Then a p-value was obtained, which was corrected by FDR method. And Corrected P-value (q-value) was calculated by correcting using BH method. p-value or q-value were used to conduct significance analysis. Parameters for classifying significantly DEGs are ≥2fold differences (|log 2 FC|≥1, FC: the fold change of expressions) in the transcript abundance and q < 0.05 (Parreira et al., 2016).
By searching the ENSEMBL, NCBI, Uniprot, GO, and KEGG databases, the BLAST (Basic Local Alignment Search Tool) alignment was performed to determine the functional annotation of DEGs. The best matches were selected to annotate the DEGs. Finally, DEGs were subjected to GO functional analysis and KEGG, utilizing default parameters, to annotate the DEGs' major GO, and KEGG categories (Liu et al., 2014;Parreira et al., 2016).

Validation of RNA-Seq Data by RT-qPCR
To validate the RNA-seq data, RT-qPCR was introduced. Several genes were selected for the validation. Some genes were selected due to their involvement in stress response according to previous studies (Beltramo et al., 2006;Olguin et al., 2009Olguin et al., , 2015, and others were randomly selected ( Table 1). The RNA samples used are same as used in RNA-seq analysis. The primers were selected and analyzed by the Primer Premier Software (version 5.0). In this work, five genes (ldhD, dpoIII, dnaG, gyrA, and gyrB) were evaluated as internal controls for RT-qPCR, using the primers described in Table 1 (Desroche et al., 2005;Costantini et al., 2011;Margalef-Català et al., 2016). The five internal controls were calculated on their geometric mean for the normalization of RT-qPCR data (Sumby et al., 2012). The Real Time PCR System iQ5 (Bio-Rad) was used for the amplification of RT-qPCR. The threshold value used in this study was automatically determined by the instrument. Results were analyzed using the 2 − CT method, and the amount of target RNA was adjusted to the geometric mean of the five internal controls as previously described (Livak and Schmittgen, 2001).

RESULTS AND DISCUSSION
To better understand the stress response and regulation mechanism of O. oeni, functional analysis based on comparative transcriptomics was used in this study. The genes most affected by acid shock were mainly studied in this paper. mRNA from the control (pH 4.8) at t = 0 h and t = 1 h, and from acid treated samples at t = 1 h were used to conduct transcriptional analysis. The RNA-seq data using in this article have been submitted to Sequence Read Archive (SRA) database with an accession number of SRP105332.
The results obtained from the RNA-seq were validated by RT-qPCR with the same RNA samples, and 11 genes were selected in this section ( Table 1). For all the 11 genes tested, the RT-qPCR data have a general accordance with RNA-seq data (Figure 1 and Supplementary Figure 1). Of the 11 genes in different groups, most were clearly correlated using both techniques. Indicating no significant changes through this technique, although some genes display low correlated in the group VS2 (pH 4.8_1 h-VS-pH 4.8_0 h). Overall, the correlation between RT-qPCR and RNA-seq is good, suggesting that the RNA-seq data are valid.

Global Analysis of Functions Affected during Acclimation after Acid Shock
In order to identify the biological processes influenced by acid shock, transcriptomic data were grouped by functional categories. pH 4.8_0 h and pH 4.8_1 h, as the reference conditions, were used to normalize data. Under the control conditions, the expression level of some genes was decreased, probably due to the influence of centrifugation (data not shown). However, acid shock is the biggest influencing factor in gene expression. Table 4 shows some DEGs from each functional category after acid shock at pH 3.0 (t = 1 h). Genes within a wide range of functional classes were influenced by acid shock.
A total of 955 DEGs were detected by the RNA-seq. It is significantly higher than those identified by Margalef-Català et al. (2016). But in the three separate comparison groups, the numbers of DEGs were almost the same or less than that in Margalef-Català et al. (2016). Of these, as in Figure 2, compared to pH 4.8_0 h, 235 genes decreased their expression 1 h after acid shock and 406 genes increased their expressions. Compare to pH 4.8_1 h, 158 genes decreased in their expression after 1 h acid shock and 249 genes increased in their expression. Compared to the research of Margalef, apart the techniques, the media (WLM in the case of Margalef), strain (PSU-1 in the case of Margalef) were different as well. These differences could be corrected by the setting of control groups. The samples of pH 4.8_0 h and pH 4.8_1 h were set as control groups in this study, while the samples at 0 h were control groups in the case of Margalef. The DEGs number of comparisons VS1 (pH 3.0_1 h-VS-pH 4.8_0 h) was almost the same as Margalef, which was much higher than VS3 (pH 3.0_1 h-VS-pH 4.8_1 h). With the group of pH 4.8_0 h as control, the result may overlooked the genes differentially   expressed due to the time changing, which is normal in the process of bacterial growth. Therefore, the group of pH 4.8_1 h was set as the primary control in this study.
All the specific DEGs numbers changed by different conditions were shown in a Venn diagram (Figure 3). Ninetynine genes were identified to be expressed in significant difference within all comparisons. The analysis of comparisons VS1 and VS2 had the same transcription patterns with 375 genes, while only 99 genes showed some differences in comparison VS3, which suggests that with different control samples, the transcription patterns are also exist differences.
A hierarchical heat map (Figure 4) was adopted to show the global DEGs patterns occurring in the experimental conditions. The expression profiles under different growth conditions were shown in this map, obviously. In this study, the key factor influencing the cluster patterns of DEGs was growth condition. Thus, the DEGs patterns were similar within parallel samples, showing good correlation between parallel samples.
All the DEGs were also shown in the format of scatter diagram and volcano plot (Additional Supplementary Figures 2, 3).

Functional Analysis and Classification of DEGs
To better understand the transcriptome of O. oeni SD-2a, the function of predicted genes was classified by GO and KEGG.
GO enrichment was used to identify the putative function of all the DEGs in every group, which can provides DEGs a statistical support in GO terms. In general, the enrichment  analyses of DEGs showed that VS1, VS2, and VS3 were mainly belong to one category: biological processes ( Table 2, Additional Supplementary Figures 4-6). Among them, the majority of DEGs of all groups (VS1, VS2, and VS3) were found to be involved in the metabolic process (GO:0008152), catalytic activity(GO:0003824), cellular process(GO:0009987), and binding(GO:0005488). Among the KEGG enrichment of DEGs, some groups seemed to be less affected by acid (like transport and catabolism, metabolism of terpenoids, and polyketides, biosynthesis of other secondary metabolites, transcription and signal transduction mechanisms), while others were more sensitive (like amino acid metabolism, carbohydrate metabolism, membrane transport, and energy metabolism) ( Table 3, Additional Supplementary Figures 7-9).

Highly Induced/Suppressed Genes
Eighty-eight genes were identified in the comparison VS3, with the limits of (1) q < 0.05 and (2) log2FC≥2 (highly induced) or log2FC≤ −2 (highly suppressed) (Supplementary Table 1). Among them, compared to the report of Margalef-Català et al. (2016), some genes showed different changes in their expressions. For example, some genes which were found decrease or increase their expressions in this study were not mentioned in Margalef 's report, and others showed opposite expression changes (Supplementary Tables 2, 3). In the Supplementary Table 2, there were some genes shown the same expression patterns between the data from group VS1 and Margalef-Català et al. (2016) but the group VS3 was crosscurrent (orf00275 and orf00218). This again demonstrates the importance of control group (pH 4.8_1 h).

Main Metabolisms Modified by Acid Shock
This study was designed to identify genes differentially expressed in pH 3.0 and pH 4.8, by using the RNA-seq technique. The comparisons were carried out between pH 4.8_0 h, pH 4.8_1 h, and pH 3.0_1 h. Table 4 showed the transcriptomic analysis of the relative expression of genes between time 0 and 1 h after the Xenobiotics biodegradation and metabolism 10 8 5

VS1, is the comparison of pH 3.0_1 h-VS-pH 4.8_0 h; VS2, is the comparison of pH 4.8_1 h-VS-pH 4.8_0 h; VS3, is the comparison of pH 3.0_1 h-VS-pH 4.8_1 h.
acid shock. The table showed a selection of the most inhibited or promoted genes with known functions. Figure 5 showed the hierarchical clustering in heat map format of DEGs shown in Table 4. From the RNA-seq data, the glycosyltransferase genes related to the carbon source in the medium were overexpressed, while the expression levels of other irrelevant genes were decreased. And the expression of transport proteins, the membrane-like ion transport proteins, amino acid transporter, etc., were significantly increased. These changes are related to the responses of O. oeni SD-2a to acid shock.

Malate and Citrate Metabolism
One of the strategies that microorganism defense the acid stress is to decrease the internal high concentration proton, as well-known, O. oeni can do this by MLF. The way of malate transport into cells is through malate permease (mleP), which was up-regulated in this study. Among the MLF, oxidative decarboxylation is an important process. One of the enzymes that catalyzes such a reaction is 3-isopropylmalate dehydrogenase (IPMDH), a member of the β-hydroxyacid oxidative decarboxylase family, to which malate dehydrogenase (decarboxylating) also belong (Pallo et al., 2014). In our study, the 3-isopropylmalate dehydrogenase (leuB) and malate dehydrogenase (maeA) genes were over-express compares to pH 4.8_0 h, which can offset the influence of low pH in some ways.  Besides, 2-isopropylmalate synthase (leuA) and malate transporter gene were also over-expressed. The observed transcriptional activation of maeA, mleP, and malate transporter under acid conditions are in accordance with previous studies about wine-related conditions (Augagneur et al., 2007;Costantini et al., 2015;Margalef-Català et al., 2016). But, the expression of maeA and mleP at pH 3.0_1 h did not have significant differences compared to pH 4.8_1 h, which have not reported before. Citrate lyase is a key enzyme of citrate fermentation, the prosthetic group of citrate lyase is catalyzed by CitG and CitX in Escherichia coli. These two genes are part of the citrate lyase gene cluster, citCDEFXG (Schneider et al., 2002). The expression of the citrate lyase operon were induced in this study, which has been previously reported over-expressed under low pH and multi-stress conditions (Olguin et al., 2009;Margalef-Català et al., 2016). The over-expression of genes related to malate transporter and citrate consumption indicated that the consume of L-malate and citrate were associated with acid stress response, and may be as an alternative energy source to sugar metabolism, just as described by Margalef-Català et al. (2016).
Significant changes were also observed within genes involved in diacetyl utilization. Diacetyl is the main aromatic compound associated to MLF and is derived from citrate consumption. Diacetyl reductase and acetoin reductase showed transcriptional inhibition. The expression patterns of these two genes at 1h after  acid shock (pH 3.0) were in accordance with previous studies (Margalef-Català et al., 2016). Since there is no subsequent time monitoring, the expression changes are not clear.

Amino Acid Transport and Metabolism
As nutrition and flavoring ingredients, amino acids play a key role in the quality of wine. They are the precursors of higher alcohols, esters, and aromatic thiols, which are the flavor active compounds of wine. During wine fermentation, their biosynthetic, and catabolic pathways also play a central role in the biosynthesis and releasing of aroma (Holt et al., 2012). Chorismate synthase (CS) was over-expressed in VS1 and VS3 comparison group. It catalyzes the biosynthesis of chorismate by 5-enolpyruvylshikimate 3-phosphate. It is the seventh enzyme in the shikimate pathway (SP), and in the biosynthesis of numerous aromatic compounds, the product of this reaction is the last common precursor in bacteria (Macheroux et al., 1999). This reaction catalyzed by CS can release two moles carbonyl, which can combine with free H + , decrease the concentration of H + . The fifth enzyme of the SP, shikimate kinase (SK), was also over-expressed after acid shock 1 h (Vianna and de Azevedo, 2012). The SP is important for the synthesis of some aromatic amino acids, such as phenylalanine, tyrosine, tryptophan, and other functional aromatic compounds, which will participate in the signaling, electron transport, UV protection, and wound response (Macheroux et al., 1999). These changes can help O. oeni to synthesize aromatic compounds and defense the damages caused by acid shock.
Argininosuccinate synthase (ASS) is involved in the biosynthesis of arginine together with argininosuccinate lyase (ASL), and ASS is the rate-limiting enzyme for arginine biosynthesis (Locke et al., 2016). They were over-expressed after acid shock 1 h, which means the up-regulation of arginine synthetase. Arginine can stimulates the expression of some stress-responsive genes, such as ftsH and omrA, and it also can increase the cell number of O.oeni at low pH (Arena and de Nadra, 2005;Bourdineaud, 2006). The up-regulation of ASS gene can promote the reproduction of O. oeni under stress conditions, for accumulating cell number to resist stress and start MLF.
The enzyme proline iminopeptidase, which releases proline from the N-terminus of small peptides, was over-expressed at 1 h after acid shock. Since peptides account for the largest proportion of total nitrogen in wine (Margalef-Català et al., 2016), it is important for O.oeni to utilize them under low pH. Due to the inhibition of proline permease by nitrogen metabolic by-products, during MLF the consume of proline is very few. But the existence of proline can improve the growth of O.oeni (Lv, 2012).
The 4-aminobutyrate aminotransferase gene, which transforms gamma aminobutyric acid (GABA) into succinate semialdehyde and L-glutamate, was 3-fold over-expressed after 1 h adaption to acid shock, as reported by Margalef-Català et al. (2016). In bacteria like Corynebacterium glutamicum and E. coli, GABA can be utilized as the form of carbon and/or nitrogen source, but its assimilation in O. oeni is not clear yet (Bartsch et al., 1990;Zhao et al., 2012;Margalef-Català et al., 2016). It is worth a further study.
The cystathionine beta-lyase (CBL) gene was over-expressed after acid shock. CBL is involved in the biosynthesis of methionine. CBL catalyzes the conversion of cystathionine into homocysteine in an α, β-elimination reaction, which will convert to methionine in a later step. The CBL activity plays an important role in aromatic thiol release (Holt et al., 2012). It can improve the flavor and quality of wine during MLF.
Acetolactate synthase (ALS) is the first rate-limiting enzyme for branched-chain amino acid biosynthesis, like valine, leucine, and isoleucine. It converts 2 mol of pyruvate to acetolactate, using thiamine diphosphate (ThDP) as a cofactor (Duggleby and Pang, 2000;Pang et al., 2002;Zheng et al., 2015). Acetolactate can participate in the synthesis pathway of diacetyl (2, 3 -butyl ketone) and its derivatives, which are the main flavor compounds generated by MLF. The up-regulation of ALS gene could help O. oeni accomplish the MLF and increase the abundance of aroma compound in wine.

Stress Response
To alleviate the challenge of reduction in internal pH caused by high concentration proton, the bacteria cytoplasm will be alkalization. Among the efflux systems of harmful-compounds and cell detoxification, one of the important parts is ABC transporters (Leverrier et al., 2004). In this study, there were 121 DEGs detected by RNA-seq related to ABC transporters. Meanwhile, as an important molecular marker of stress response in O. oeni, the protein Hsp20 was also over-expressed in this essay. ATPase activity has been associated to MLF. Fortier et al. (2003) described the increase of F 0 F 1 -ATPase β subunit mRNA in response to low pH (Fortier et al., 2003). In this work several genes codifying other ATPase subunits (β, γ, and ε) were upregulated after the acid shock at pH 3.0_1 h (Table 4). However, the F 0 F 1 ATP synthase subunit α was up-regulated under acid conditions compare to pH 4.8_0 h, but did not have significant differences compare to pH 4.8_1 h. This result is opposite to the report of O.oeni PSU-1 under wine like medium (1 h), and similar with the situation after 6 h inoculation (Margalef-Català et al., 2016), meanwhile it indicate that when cells are exposing at low pH, the ATPase activity is increased more quickly than wine like medium, and agrees with the role of this enzyme in the regulation of the cytoplasmic pH and in the acid stress response of O. oeni.
The D-alanyl-D-alanine carboxypeptidase (dacC) gene related to cell envelope biogenesis was over-expressed, and it was 6-fold over-expressed at 1 h, and also over-expressed in transcriptomic analysis by Costantini et al. (2015) and Margalef-Català et al. (2016) after adaption with ethanol and WLM. This result is consistent with earlier reports on the barrier and homeostasis functions of cell membranes in the stress response of O.oeni, which is well-known (Grandvalet et al., 2008). But previous studies showed that, several genes related to cell wall biosynthesis were significant differentially expressed (Margalef-Català et al., 2016), which point out the relevance of some genes involved in cell wall protection against stress challenges. This point was also confirmed by our study. The expression level of phosphoglycerol transferase gene was significant up-regulated at 1h after acid shock. The dacC gene is involved in the pathway of lipoteichoic acid biosynthesis, and is a part of cell wall biogenesis. The role of cell wall in the stress response of O.oeni is worth a further study.

CONCLUSIONS
This is the first transcriptome study using RNA-seq on O.oeni under different conditions. The RNA-seq study is useful to identify the metabolisms mostly altered due to low pH conditions. Our results revealed the relevance of carbohydrate metabolism, amino acid metabolism and membrane transport as key metabolisms involved in the adaptation of O.oeni SD-2a to acid stress. From GO analysis, the majority of DEGs of all groups (VS1, VS2, and VS3) were found to be involved in the metabolic process, catalytic activity, cellular process and binding. In addition, a considerable proportion of genes are novel, which have a significantly differently expression in this study. These results provide a new viewpoint and crucial resource on the acid stress response in O. oeni.

AUTHOR CONTRIBUTIONS
TW, JS, and HW conceived the idea of the work. LL, HZ and SP designed the experiments and performed the experiments. LL, YL, HL, and HW analyzed the data and wrote the manuscript.

FUNDING
This work was supported by National Natural Science Foundation of China (Grant No. 31471708). This work was also financially supported by Shaanxi special finance for agriculture "Construction of technological system for Shaanxi vitis industry-2016."

ACKNOWLEDGMENTS
The authors would like to thank reviewers for their comments and suggestions which greatly improved the original version of the article.

SUPPLEMENTARY MATERIAL
The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fmicb. 2017.01586/full#supplementary-material Supplementary Figure 1 | Validation of RNA-seq data using RT-qPCR. Eleven representative genes were chosen to validate the RNA-Seq data by RT-qPCR. The black bars represent mean values of log2-transformed fold change obtained from three biological replicates of RT-qPCR with error bars stand for standard deviations. And the red bars represent RNA-Seq data. (A-E) Represent gene dnaG, dpoIII, gyrA, gyrB, and ldhD as internal controls, respectively. The number 1-3 represent group pH 3.0_1 h-VS-pH 4.8_0 h, pH 4.8_1 h-VS-pH 4.8_0 h, and pH 3.0_1 h-VS-pH 4.8_1 h, respectively. Frontiers in Microbiology | www.frontiersin.org