Transcriptomic analysis of Chlorella sp. HS2 suggests the overflow of acetyl-CoA and NADPH cofactor induces high lipid accumulation and halotolerance

of Chlorella sp. suggests the overflow of acetyl-CoA and NADPH cofactor induces high lipid accumulation and halotolerance. Abstract Previously, we isolated Chlorella sp. HS2 (referred hereupon as HS2) from a local tidal rock pool and demonstrated its halotolerance and high biomass productivity under different salinity conditions. To further understand acclimation responses of this alga under high salinity stress, we performed transcriptome analysis of triplicated culture samples grown in freshwater and marine conditions at both exponential and stationary growth phases. The results indicated that the transcripts involved in photosynthesis, TCA, and Calvin cycles were downregulated, whereas the upregulation of DNA repair mechanisms and an ABCB subfamily of eukaryotic type ABC transporter was observed at high salinity condition. In addition, while key enzymes associated with glycolysis pathway and triacylglycerol (TAG) synthesis were determined to be upregulated from early growth phase, salinity stress seemed to reduce the carbohydrate content of harvested biomass from 45.6 dw% to 14.7 dw% and nearly triple the total lipid content from 26.0 dw% to 62.0 dw%. These results suggest that the reallocation of storage carbon toward lipids played a significant role in conferring the viability of this alga under high salinity stress by remediating high level of cellular stress partially resulted from ROS generated in oxygen-evolving thylakoids as observed in a direct measure of photosystem activities.


| INTRODUCTION
Microalgae exhibit a greater biomass yield than most terrestrial crops and can be grown with excess nutrients in wastewater sources, prompting its industrial utilization as a biofeedstock for the production of nutraceuticals, pharmaceuticals, cosmetics, and biofuels (Hu et al., 2008;Quinn & Davis, 2015;Smith et al., 2010;Unkefer et al., 2017;Yun, Cho, Lee, Heo, et al., 2018). However, commercial production of algal biomass is not yet considered to be economically competitive because of high energy inputs associated with biomass harvesting and downstream extraction of desirable biomolecules (Laurens et al., 2017;Stephens et al., 2010;Valizadeh Derakhshan et al., 2015). Importantly, the productivity and operational stability of algal cultivation platforms are prone to be compromised by unpredictable meteorological conditions and culture contamination (McBride et al., 2014;Wang et al., 2016;Yun et al., 2016Yun et al., , 2019, which has led to multifactorial efforts to develop robust algal "crops" under changing environments, just as in the case of conventional agriculture. Of environmental conditions that determine the productivity of biomass and desirable biomolecules from industrial crops, salinity appears on the top of the list because of high crop sensitivity to presence of high concentrations of salts in the soil or waters (Flowers et al., 1977;Peng et al., 2014;Yuge Zhang & Liang, 2006). In particular, the extensive application of chemical fertilizer facilitates accumulation of salts in agricultural fields, which in turn could lead to a positive feedback loop by necessitating an increased application of synthetic fertilizer (Yuge Zhang & Liang, 2006). Notably, industrial algal cultivation platforms require continuous provision of nutrient salts with some studies demonstrating the utilization of saline wastewater sources enriched with nitrogenous and phosphorus nutrients as growth media to drive down the costs of commercial operation of algal cultivation systems Yun et al., 2015;Zhu et al., 2013). In addition, the direct application of salinity stress for algal cultivation systems has been demonstrated as an effective abiotic inducer of high lipid accumulation and an environmental barrier inhibiting the proliferation of undesirable alien invaders in cultivation systems (Church et al., 2017;Kakarla et al., 2018;Lee et al., 2016). Kakarla et al., for instance, supplemented 60 g/L of NaCl into concentrated Chlorella cultures for 48 h and reported ca. 58% increase in algal lipid productivity, supporting the possibility of deploying high salinity stress as a promising post-treatment for the cultivation systems targeting to produce algal lipids (Kakarla et al., 2018). Moreover, while high salinity stress could act as an effective method of crop protection in reducing freshwater cyanobacterial or ciliate contaminants, it was successfully demonstrated to facilitate algal harvesting by enlarging cellular diameter and increasing algal settling rates (von Alvensleben et al., 2013;Church et al., 2017;Lee et al., 2016). Even though general osmosensitivity of algal crops has been acknowledged (Flowers et al., 1977), there is thus a great industrial incentive to exploit algal diversity and especially high tolerance of some algal species to highly saline environment (Yun et al., 2015).
With the apparent advantages of incorporating high salinity stress into the management of industrial algal cultivation platforms, bioprospecting halotolerant algal strains that exhibit high and reliable production of biomass and/or desirable biomolecules was the motivation of our previous study in which a halotolerant Chlorella sp. was isolated from a tidal rock pool (Yun et al., 2019). While the remarkable toughness of Chlorella under different physical and chemical stress and its recognition as one of a handful of successful industrial crops have been well documented (Fogg, 2001;Yun et al., 2019), this isolated Chlorella sp. HS2 (referred to hereupon as HS2) exhibited relatively high growth under a wide range of salinity conditions (i.e., 0%-7% (w/v) of supplemental NaCl) compared to reference Chlorella strains (Yun et al., 2019). Importantly, substantial shifts in the composition of fatty acid methyl ester (FAME) and the amount of carotenoid pigments under different salinity conditions led us to speculate that elucidating mechanisms behind relatively short-term (i.e., few days) algal acclimation to high salinity stress would enable maximizing the industrial potential of HS2 by guiding ongoing efforts in metabolic and process engineering (Oh et al., 2019;Rathinasabapathi, 2000;Yun et al., 2019).
In previous studies, transcriptome analysis has served as an important tool to understand intricate algal responses to changing salinity conditions. For example, Foflonker et al. challenged Picochlorum cells with high or low salinity shock and used transcriptomic and chlorophyll fluorescence analyses to elucidate salinity-tolerance mechanisms (Foflonker et al., 2016); the authors identified photoprotective mechanisms, oxidative stress response, cell wall and membrane rearrangement, nitrogen assimilation, and diverting resources from growth and PSII repair in favor of maintaining homeostasis as the main responses against a challenging environment (Foflonker et al., 2016). Moreover, Perrineau et al. compared salt-acclimated and progenitor populations of Chlamydomonas reinhardtii, and reported downregulation of acetyl-CoA, Chlorella sp. HS2, halotolerance, lipid synthesis, photosynthesis, RNA-seq genes involved in the salt stress response (most notably, glycerophospholipid signaling) and in transcription/translation in the salt-acclimated populations, suggesting gene-rich mixotrophic algal lineages could rapidly adapt to high salinity conditions (Perrineau et al., 2014). Importantly, the survey of existing literature suggested the presence of strain-specific algal responses that could be closely associated with the phenotypic characteristics of an algal strain of interest (Erdmann & Hagemann, 2001).
Herein, we report the transcriptome of HS2 grown in freshwater and marine conditions to accomplish mechanistic understanding of algal acclimation to high salinity stress. Triplicated cultures samples were first obtained at exponential and stationary growth phases in freshwater and marine growth media for RNA-seq analysis, and the proximate analysis of harvested biomass was additionally performed along with the measure of photosystem II (PSII) activity. Combined together with the results in our previous study, we were able to elucidate how vital metabolic pathways were shifted under high salinity stress, and an important role of allocating storage carbon toward the synthesis of lipids in conferring the viability of HS2 and remediating high oxidative stress under high salinity stress.

| Strain selection and cultivation
conditions HS2 was previously isolated from a local tidal rock pool, and its high tolerance to a wide range of salinity conditions was acknowledged (Yun et al., 2019). While the results of HS2 cultivation in 1-L cylindrical PBRs were reported in our previous study (Yun et al., 2019), both autotrophic cultures grown in freshwater inorganic medium and in marine inorganic growth medium supplemented with 3% (w/v) sea salt were subjected to transcriptome analysis. These triplicated cultures were grown under pre-determined optimal light and temperature conditions with continuous supplementation of 5% CO 2 at 0.2 vvm and agitation at 120 rpm.

| PSII activity measurement and proximate analysis
While pigment and FAME composition of harvested HS2 biomass in both freshwater and marine conditions were reported previously, photoautotrophically grown cells in exponential and stationary growth phases were subjected to measurements of the photosynthetic parameters in vivo using Multi-Color-PAM (Heinz Walz, Germany) (Shin et al., 2017). After adapting cells under dark condition for 20 min, the light response curves of the relative electron transport rate (rETR), the quantum yields of non-photochemical quenching (Y(NPQ)), and non-regulated excess energy dissipation (Y(NO)) were measured in biological triplicates while increasing the actinic light intensities of 440 nm LEDs with a step width of 2 min (Shin et al., 2017). In addition, proximate analysis of the biomass harvested at stationary growth phase was performed in biological triplicates to further elucidate metabolic shifts in HS2 under high salinity stress. The lipid content of harvested biomass was first analyzed by extracting total lipids from freeze-dried biomass with chloroformmethanol (2:1 (v/v)) following a slightly modified version of Bligh and Dyer's method (Bligh & Dyer, 1959). Samplesolvent mixtures were then transferred into a separatory funnel and shaken for 30 min and the lipid fraction was separated from the separatory funnel; the solvent was evaporated using a rotary evaporator and the weight of the crude lipid obtained from each sample was measured using an analytical balance following Yun, Cho, ). In addition, the protein content was determined using the method of Lowry using ca. 2 mg (dry weight) of the cell pellet resuspended in 0.5 ml of 1 M NaOH and boiled for 5 min (Illman et al., 2000;Lowry et al., 1951); the carbohydrate content was measured using the phenol sulfuric acid method of Dubois et al. using ca. 0.5 mg (dry weight) of the cell pellet resuspended in 1 ml of water (Dubois et al., 1956;Illman et al., 2000). Finally, the ash content was analyzed gravimetrically after exposing dry biomass to 500°C in a muffle furnace for 8 h (Kent et al., 2015).

| RNA extraction, library construction, and Illumina sequencing
Each of salt-stressed and control PBR cultures was harvested during exponential and stationary growth phases by centrifugation at 4105 g for 10 min. Total RNA was then extracted using the Trizol reagent (Invitrogen, Carlsbad, CA, USA), according to manufacturer's instructions. Subsequently, the RNA samples were treated with DNase I for 30 min at 37°C to remove genomic DNA contamination, and the quantity and integrity of the total RNA were verified using an Agilent 2100 bioanalyzer. The cDNA libraries were developed according to manufacturer's instructions (Illumina, Inc., San Diego, CA, USA), and sequenced on the Illumina HiSeq 2000 platform at Seeders Co. (Daejeon, Korea) (Liu et al., 2017). In addition, RNA-Seq paired-end libraries were prepared using the Illumina TruSeq RNA Sample Preparation Kit v2 (catalog #RS-122-2001, Illumina, San Diego, CA). Starting with total RNA, mRNA was first purified using poly (A) selection or rRNA depletion, then RNA was chemically fragmented and converted into single-stranded cDNA using random hexamer priming; the second strand was generated next to create double-stranded cDNA. Library construction began with generation of blunt-end cDNA fragments from ds-cDNA. Thereafter, A-base was added to the blunt-end in order to make them ready for ligation of sequencing adapters. After the size selection of ligates, the ligated cDNA fragments which contained adapter sequences were enhanced via PCR using adapter specific primers. The library was quantified with KAPA library quantification kit (Kapa biosystems KK4854) following the manufacturer's instructions. Each library was loaded on Illumina Hiseq2000 platform, and the desired average sequencing depth was met while performing high-throughput sequencing.

| De novo assembly and analysis
De novo assembly was performed using Trinity 2.8.5 (Grabherr et al., 2011) using raw 100 bp paired-end reads. Assembly quality assessment was carried out with BUSCO 3.0.2 (Simão et al., 2015), for which the chlorophyte database of OrthoDB 10 was employed as datasets at an e-value cutoff of 1e−5 (Kriventseva et al., 2018); high-quality reads were mapped onto genome sequences using Bowtie2 2.3.5. Thereafter, the quantification of the number of reads (i.e., counts mapped per transcripts) was performed following alignment and abundance estimation of each Trinity script using RSEM 1.3.2 and Bowtie 1.2.2, respectively (Langmead et al., 2009;Li & Dewey, 2011). Transcripts with no count across all sampling points were removed. The matrix of counts for unigenes (i.e., a collection of expressed sequences that are aligned or located to the same position on genome) was used for downstream analyses.

| DEG analysis and functional annotation
Prior to functional annotation, differential expression analysis (DEA) was performed first to avoid determining the most relevant transcript for each unigene based on unnecessary assumptions at the early stage. In addition, given that quantitative asymmetry between up-and downregulated unigenes was strong, SVCD 0.1.0, which does not assume the lack-of-variation between up-and downregulated unigene counts (Evans et al., 2017;Roca et al., 2017), was used for normalization of unigenes. The mean of raw counts greater than the first quartile (i.e., 5.9 raw counts) as recommended (Roca et al., 2017) was used during normalization. To determine DEGs, we used DESeq2 1.20.0, and the DEGs between exponential and stationary growth phases were based on the adjusted p-values (i.e., DEGs were determined as unigenes with adjusted p-value < 0.01).
Functional annotation of DEGs was subsequently performed using Swiss-Prot, Pfam, and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases. First, following Trinotate 3.2.0's recommendation, we predicted transcript coding regions that could be assigned to putative proteins using TransDecoder 5.5.0 (Haas et al., 2013). Thereafter, homologies were identified using in parallel BLASTp from BLAST+ 2.9.0; to identify pfam domains, hmmscan from HMMER 3.2.1 was used (Camacho et al., 2009;Eddy, 2011). BLASTp and hmmscan were run twice from the predicted proteins. SignalP 5.0b (http://www.cbs.dtu.dk/servi ces/Signa lP/) was used to determined eukaryotic signal peptides within transcripts. We also used BLASTx to find homologues, which allows to identify sequence similarities within all six reading frames of the transcript. All BLAST runs were performed against the Swiss-Prot database through DIAMOND 0.8.36 (Buchfink et al., 2015) with an e-value cutoff of 1e−10. Then, KEGG cross-references associated with BLASTx or BLASTp hits were retrieved to assign each BLAST hit with a KEGG Orthology number (KO). Transcripts without a BLASTx or BLASTp hit were excluded, and a pair of transcript and coding region was removed when the KOs of corresponding transcript and coding regions were not identical. In addition, when one gene had multiple KOs, the mean of average e-values was computed and the KO with the lowest mean was selected as the most relevant KO. Metabolic pathway maps were constructed using KEGG mapper based on the organism-specific search results of Chlorella variabilis (cvr) and biological objects for each KO were determined using KEGG BRITE. Enrichment was performed by implementing GSEAPreranked from Gene Set Enrichment Analysis with the conda package GSEApy 0.9.15 (Mootha et al., 2003;Subramanian et al., 2005). A term was considered to be significantly enriched when its false discovery rate (FDR) was lower than 0.25. All data generated from our transcriptome analysis are available at the NCBI GEO repository: GSE146789 at https://www.ncbi.nlm.nih.gov/geo/query/ acc. cgi?acc=GSE14 6789.

| Phenotypic shifts of HS2 under high salinity stress
Shifts in growth, FAME, and pigment composition of HS2 during autotrophy in freshwater (i.e., 0% (w/w) of supplemental sea salt) and marine (i.e., 3% (w/w) of supplemental sea salt) media were reported in the previous study (Yun et al., 2019). Briefly, the results indicated a nearly 10-fold decrease in the maximum cell density of the autotrophic PBRs in marine medium at stationary growth phase, whereas only a twofold decrease in the average dry cell weight (DCW) was observed (Yun et al., 2019) ( Figure S1). As microscopic observation revealed, a non-proportional decrease in DCW of HS2 under high salinity stress corresponded to roughly 50% increase in cellular diameter or 3.4-fold increase in cellular volume. While previous study also reported substantial decreases in the amount of algal pigments and relative amount of polyunsaturated fatty acids under high salinity stress (Yun et al., 2019), TEM images of harvested algal cell suggested the formation of large lipid droplets under high salinity stress ( Figure 1): Indeed, proximate analysis of harvested biomass indicated a significant increase in lipid content from 25.0 dw% to 62.0 dw% under high salinity stress, contrasting a nearly three-fold decrease in the amount of carbohydrate (Figures 1 and 2).
While relatively high amounts of carotenoid pigments (i.e., β-carotene and lutein) under high salinity stress observed in the previous study suggested their possible contribution to the protection of photosynthetic machinery (Talebi et al., 2013;Yun et al., 2019), the measures of relative electron transport rate (rETR), the quantum yields of non-photochemical quenching (Y(NPQ)) and non-regulated excess energy dissipation (Y(NO)) using multi-color-PAM indicated that rETR was reduced early during the exponential growth phase under high salinity stress and was recovered at later stationary growth phase. Although differences in Y(NPQ) and Y(NO) were not observed, respectively, at exponential and stationary phases, a significant difference in Y(NPQ) was observed during stationary phase only at high light intensities and Y(NO) of salt-shocked culture was significantly greater than that of control across all light intensities during exponential growth phase (Figure 3).

| Summary of de novo assembly
To determine differential transcriptomic regulation of HS2 under freshwater and marine conditions, RNA-seq was performed using Illumina Hiseq 2000 platform, followed by de novo RNA-seq assembly and mapping of data to the newly assembled and processed transcriptome. Alignment statistics from Trinity and Bowtie2 2.3.5 mapping results were summarized in Table S1. Overall, 57,640 unigenes were obtained out of 290 million raw reads, and the assessment of assembly quality indicated 89% of complete BUSCOs following the removal of 4870 unigenes with 0 count in any of the treatments.

| Functional annotation of differentially expressed genes
To elucidate differentially expressed genes (DEGs), read normalization was first performed using SVCD normalization following standard DEGseq2 statistical test; a total of 9117 DEGs were subsequently obtained from 52,770 unigenes corresponding to 39,469 transcripts. While 3573 DEGs were commonly observed across all conditions, 2334 and 3120 DEGs were distinctively observed at exponential and stationary phases, respectively ( Figure 4). Overall, global observation of transcriptome changes indicated general transcriptional downregulation under high salinity stress, highlighting substantial metabolic constraints and subsequent biochemical shifts that presumably facilitated the survival of algal cells under high salinity stress. It should be also noted that a substantial difference in terms of the overall DEG expression was observed between exponential and stationary growth phases, with more transcriptional shifts toward downregulation during stationary growth phase. Finally, KO annotation of DEGs yielded 2795 DEGs (i.e., 31% of all DEGs) with 1982 unique consensus KOs, and these DEGs represented one third of genes of Chlorella variabilis NC64A's genome (Eckardt, 2010).

| Functional enrichment of differentially expressed genes
Enrichment analysis was performed with the first and second elements of functional hierarchies of KEGG BRITE. While the terms with a p-value lower than 0.05 and a false discovery rate (FDR) equal to or lower than 0.25 were considered to be enriched, the results indicated high enrichment of ribosomal proteins ( Figure 5). In addition, papain family of intramolecular chaperones and heparan sulfate/heparin glycosaminoglycan binding proteins were enriched. Notably, even though FDR values below the cutoff were not observed, many enriched terms with a p-value lower than 0.05 were related to protein processing and membrane trafficking.

| KEGG pathway analysis
To elucidate metabolic pathways associated with the acclimation of HS2 to high salinity stress, we mapped DEGs to 120 reference KEGG pathways; pathways enriched with 20 or more DEGs were summarized in Table S2. 3.5.1 | Genes involved in cell cycle and DNA replication Upon exposure to high salinity stress, the growth of HS2 seemed to be inhibited with an apparent enlargement of cellular biovolume (see 3.1). Correspondingly, most unigenes homologous to genes identified to be involved in cell cycle were downregulated (Table 1). Additionally, DNA replication seemed to be downregulated as well, although Mcm4 of MCM complex (helicase) and DNA polymerase delta subunit 1 [EC: 2.7.7.7] were upregulated (Appendix S1), suggesting the inhibition of DNA replication under high salinity stress. Likewise, most of the unigenes associated with RNA degradation seemed to be downregulated under high salinity stress (Table 1), except CNOT3 (Appendix S1). Furthermore, most genes associated with RNA transport seemed to be downregulated under high salinity stress; and genes associated with aminoacyl-tRNA biosynthesis were downregulated, except glutaminyl-tRNA synthetase [EC: 6.1.1.18] and cysteinyl-tRNA synthetase [EC: 6.1.1.16] ( Table 1 and Appendix S1). Although these results generally supported the impairment of both DNA and RNA processing under high salinity stress, it should be emphasized that a number of unigenes associated with repair mechanisms (i.e., nucleotide excision repair, base excision repair, mismatch repair) seemed to be upregulated (Appendix S1).

| Genes involved in protein processing, MAPK signaling pathway, and ABC transporters
While salinity stress is known to substantially influence the processing and function of protein (Erdmann & Hagemann, 2001;Perrineau et al., 2014), the results indicated the downregulation of enzymes associate with protein processing in endoplasmic reticulum, except mannosyl-oligosaccharide alpha-1,3-glucosidase [EC:3.2.1.207] (GIcII), protein disulfide-isomerase A6 [EC: 5.3.4.1], and protein transport protein SEC24 (Table 1 and Appendix S1). Moreover, most of the ribosomal proteins were downregulated under high salinity stress: of 89 unigenes enriched on KEGG mapper's ribosome pathway, only S9, S16, and S26e of ribosomal proteins seemed to be upregulated at the exponential or stationary growth phases. In addition, while mitogen-activated protein kinase (MAPK) signaling cascades are widely recognized for their role in stress response and signal transduction in eukaryotes (Yang et al., 2018), most of the genes associated with MAPK signaling pathway seemed to be downregulated, except P-type Cu + transporter (RAN1) ( Table 1). Although enriched unigenes indicated that all of the genes associated with protein export were also downregulated under high salinity stress, 3 protein subunits associated with the PA700 (base) of proteasome seemed to be upregulated along with an ABCB subfamily of ABC transporters (i.e., ATM) under high salinity stress (Appendix S1).

Calvin cycle
There was a clear trend that all of the genes associated with PSII and PSI were downregulated from exponential phase under high salinity stress, corroborating with the measure of PSII activity that indicated a significant reduction in rETR during early growth phase. It should be, however, noted that these genes seemed to be less-downregulated or reverse its downregulation at later stationary growth phase (Table 1 and Appendix S1). Notably, there were more than threefold downregulation of transcripts (based on log 2 fold change) associated with PSI-D, -E, -F, -H, -K, and -O subunits and PSII Psb27 protein during exponential growth phase under high salinity stress; however, most of the transcripts associated with these subunits were upregulated during stationary growth phase, except those associated with PSI-K and PSII Psb27, which exhibited the downregulation with less than an absolute log 2 fold change of 1.0 (Appendix S1). Similarly, all of the proteins associated with light-harvesting complex (LHC) of HS2 seemed to be downregulated initially under high salinity stress at transcriptional level, whereas Lhcb2 and Lhcb4 were upregulated at the later growth phase.
While these results suggested an early compromise in photosynthesis, it should be pointed out that most of the enriched genes involved in carbon fixation via Calvin cycle were downregulated as well (Table 1 and Appendix S1). However, the upregulation of alanine transaminase [EC: 2.6.1.2] was observed under high salinity stress and no differential expression in RUBISCO [EC: 4.1.1.39] was observed. In addition, although the results of our transcriptome analysis did not indicate differential expression of ferredoxin-NADP + reductase, an enzyme that catalyzes the reaction generating NADPH in PSI (Medina & Gómez-Moreno, 2004), malate dehydrogenase (oxaloacetate-decarboxylating) (NADP + ) [EC: 1.1.1.40], the third-class malic enzyme that catalyzes the oxidative decarboxylation of malate to pyruvate by the reduction of NADP + into NADPH was upregulated during the exponential growth phase (Spaans et al., 2015). Furthermore, our transcriptome analysis suggested that glucose-6-phosphate 1-dehydrogenase [EC: 1.1.1.49], one of the key enzymes involved in the generation of NADPH during the oxidative phase of pentose phosphate pathway, was substantially upregulated (Spaans et al., 2015). It is thus likely that these enzymes associated with central carbon metabolism played a significant role in enhancing NADPH supply upon the induction of high salinity stress.3.5.4. Genes associated with glycolysis and TCA cycle. High salinity stress seemed to induce the upregulation of important genes associated with the conversion of glucose to acetyl-CoA (Table 1 and Appendix S1). In particular, pyruvate dehydrogenase E1 component alpha subunit [EC: 1.2.4.1], which is involved in the first step of converting pyruvate to acetyl-CoA was upregulated along with pyruvate decarboxylase [EC: 4.1.1.1]. Moreover, phosphoglucomutase [EC: 5.4.2.2], the enzyme involved in the first step of glycolysis, was upregulated. On the contrary, our results clearly indicated the downregulation of TCA cycle under high salinity stress: most unigenes corresponded to the known genes on TCA cycle were downregulated, suggesting the inhibition of cellular respiration (Table 1 and Appendix S1). In particular, 3 transcripts associated with citrate synthase [EC: 2.3.3.1], which mediates the first step of TCA cycle of converting acetyl-CoA to citrate, were substantially downregulated during both growth phases; and a transcript associated with isocitrate dehydrogenase [EC: 1.1.1.42], which catalyzes the rate-limiting step of the oxidative decarboxylation of isocitrate to α-ketoglutarate, was downregulated during exponential growth phase (Bellou & Aggelis, 2013). Collectively, these results suggested that acetyl-CoA became more available for other cellular metabolisms, including lipid synthesis, under high salinity stress (Bellou & Aggelis, 2013).

| Genes associated with fatty acid and TAG accumulation
Although the genes involved in the synthesis of fatty acids at upstream were downregulated, fatty acyl-ACP thioesterase A [EC: 3.1.2.14] and acyl-desaturase [EC: 1.14.19.2] were upregulated. Provided that the combined amount of C16:1, C18:0, and C18:1 was increased under high salinity stress (Yun et al., 2019), it is especially notable that these two upregulated genes are directly associated with the synthesis of these groups of mono-saturated or saturated fatty acids. Moreover, while the genes enriched on KEGG mapper indicated that fatty acid elongation and the biosynthesis of unsaturated fatty acids were not upregulated, survey of the fatty acid degradation pathway indicated the inhibition of fatty acid degradation under high salinity stress (Table 1 and Appendix S1). Most notably, transcripts associated with acyl-CoA dehydrogenase [EC:1.3.8.7], enoyl-CoA hydratase [EC:4.2.1.17], and acyl-CoA oxidase [EC:1.3.3.6] were substantially downregulated during both exponential and stationary growth phases. Given that these enzymes facilitate fatty acid β-oxidation in mitochondria or in peroxisome (Gross, 1989;Kong et al., 2017), the results suggested their role in decreasing fatty acid turnover rate and in possibly preserving fatty acids under high salinity stress.
As the upregulation of lipid synthetic pathway in marine medium was postulated based on the increased lipid content in harvested biomass (see 3.1), the transcriptome analysis also identified that genes essential for the synthesis of triacylglycerol (TAG) were upregulated: both phosphatidate phosphate [EC: 3.1.3.4] and diacylglycerol O-acyltransferase 2 [EC: 2.3.1.20] that both are involved in the conversion of 1,2,-Diacyl-sn-glycerol 3-phosphate to 1,2,-Diacyl-snglycerol and in the generation of TAG from 1,2,-Diacyl-snglycerol seemed to be substantially upregulated under high salinity stress during early growth phase.

| Genes associated with carotenoid synthesis
Of 5 unigenes enriched on KEGG mapper's carotenoid biosynthesis pathway, all of the genes were downregulated, including a gene involved in the conversion of alpha-carotene to lutein (i.e., carotenoid epsilon hydroxylase [EC: 1.14.14.158]) (Appendix S1). In addition, two genes associated with the conversion of phytoene to lycopene, an important intermediate for the synthesis of other carotenoids, were downregulated (i.e., zeta-carotene isomerase [EC: 5.2.1.12] and zeta-carotene desaturase [EC: 1.3.5.6]) (Appendix S1). Interestingly, both relative and absolute amounts of lutein were increased under high salinity stress (Yun et al., 2019); these results suggest the provision of far-upstream precursors could have played an important role in lutein synthesis.

| DISCUSSION
Given that high salinity stress strongly influences the viability and biochemical composition of algal crops and thus the economic feasibility of entire algal biorefinery (Kakarla et al., 2018;Laurens et al., 2017;Oh et al., 2019), this study was set out to elucidate transcriptional responses that give rise to the salt tolerance of highly productive HS2. While genetic engineering approaches have been extensively explored with an aim of obtaining robust algal crops, the results clearly indicated that halotolerant HS2 undergoes systematic acclimation responses against high salinity stress, identifying potential target pathways of interest for further genetic modifications or process optimization efforts (Ajjawi et al., 2017;Oh et al., 2019;Qiao et al., 2017). Of these acclimation responses, our results particularly identified a significant role of allocating available carbon toward the synthesis of algal lipids.
These results support a preferential role of lipid as a carbon and energy reserve under growth-inhibiting stress in HS2. Being an energy-dense biomolecule, previous studies indeed identified the role of lipids as a reserve facilitating cellular survival and growth upon the alleviation of growth-inhibiting stress conditions (Juergens et al., 2016). Similarly, our results indicated the upregulation of enzymes associated with glycolysis and the accumulation of lipid throughout entire growth stages: These results clearly suggest that a "push" of the acetyl-CoA precursor from glycolysis toward lipid synthesis is a major driver of lipid accumulation. Accordingly, the shift in the allocation of storage carbon resulted in an increase in algal lipids and a corresponding decrease in carbohydrates from the harvested biomass. In addition, KEGG pathway analysis of carotenoid synthesis pathway and TCA cycle suggested that these competing pathways for the "pulling" of acetyl-CoA precursor were downregulated, thereby positively contributing to the redirection of acetyl-CoA toward glycerolipid synthesis ( Figure 6).
Recent studies, however, further revealed that lipid droplets are essential and dynamical components of the cellular stress response in terms of maintaining energy and redox homeostasis (Jarc & Petan, 2019), suggesting another important metabolic function of algal lipids besides simple storage reserve. In particular, the accumulation of TAG and/ or starch could prevent cellular damage by utilizing excess photosynthetic energy and/or carbon inputs as postulated in the overflow hypothesis (OH) (Juergens et al., 2016;Neijssel & Tempest, 1975;Tan & Lee, 2016). Provided that Y(NO) of PSII represents non-regulated losses of excitation energy and thus indirectly indicate the relative amount of reactive oxygen species (ROS) (GmbH, 2012;Klughammer & Schreiber, 2008), our results suggested a strong reduction of PSII accepters and photodamage via formation of ROS during early growth phase under high salinity stress, which seemed to be subsequently resolved at stationary phase with no substantial compromise in non-photochemical quenching (NPQ). In addition, while our results indicated no differential expression of D1 protein of HS2 under high salinity stress, overall downregulation of protein processing, including subunits of the proteasome, under high salinity stress hints at a decrease in D1 protein turnover in PSII (Andersson & Aro, 2001;Erdmann & Hagemann, 2001), which likely further contributes to the increased oxidative stress due to the inhibition of the recovery of damaged PSII and could elicit cellular remediative responses, including lipid synthesis (Zhang et al., 2000).
Importantly, the synthesis of glycerolipid necessitates NADPH as a cofactor (Tan & Lee, 2016): being an electron donor, NADPH is synthesized along with ATP during the light reaction of photosynthesis and has been acknowledged for its role as an oxidative stress mediator (Valderrama et al., 2006). It should be, however, noted that there was no substantial upregulation of ferredoxin-NADP + F I G U R E 6 Simplified scheme of carbon and energy flows in Chlorella sp. HS2 for putative early responses against high salinity stress. Red and blue dashed arrows, respectively, indicate upregulation and downregulation of a given conversion or response based on transcriptome or phenotypic analyses. Glycerate-3p, Glycerate-3-phosphate; NADPH, Nicotinamide adenine dinucleotide phosphate; ROS, Reactive Oxygen Species; TAG, Triacylglycerol reductase in photosystems based on our transcriptome analysis. Nonetheless, the upregulation of glucose-6-phosphate 1-dehydrogenase and malate dehydrogenase (oxaloacetate-decarboxylating) (NADP + ) suggests that these enzymes coupled to central carbon metabolism likely made a substantial contribution to an increased NADPH pool in HS2 under high salinity stress (Spaans et al., 2015). Furthermore, given that KEGG pathway analysis suggested the downregulation of Calvin cycle under high salinity stress, the excess NADPH not utilized in carbon fixation was likely to be also directed to the high accumulation of fatty acid and/or glycerolipid, which in turn could play an important role in remediating excess oxidative stress in PSII ( Figure 6).
In addition to carbon allocation to lipid accumulation, common cellular responses under high salinity stress involve the upregulation of anti-oxidative enzymes, including catalase, superoxide dismutase (SOD), and glutathione reductase (GR) as well as the upregulation of DNA repair mechanisms and ABC transporters (Fu et al., 2014;Huang et al., 2006;Valderrama et al., 2006). Although substantial upregulation of anti-oxidative enzymes was not observed at least at transcriptional level, the degree to which each mitigation response contributes to the overall acclimation of HS2 under high salinity stress across different growth stages remains to be elucidated. Importantly, the results also indicated upregulation of P-type Cu + transporter (RAN1) on MAPK signaling pathway in HS2-the activity of RAN1 was determined to be positively correlated with plant cold resistance; overexpression of RAN1 was further reported to increase abiotic stress tolerance in Arabidopsis thaliana (Xu & Cai, 2014;Xu et al., 2016;Yang et al., 2018). Moreover, the increased relative proportion of saturated and mono-saturated fatty acids in HS2 under high salinity stress corresponded to the upregulation of enzymes involved in the synthesis of palmitoleate (C16:1), stearate (C18:0), and oleate (C18:1n9c) (Guo et al., 2019). Hence, the putative remediation of oxidative stress under growth-inhibiting high salinity condition could concurrently involve signal transduction and a shift in membrane fluidity (Guo et al., 2019), in addition to directing acetyl-CoA precursor and excess cofactor toward lipid synthesis.
While the orchestration of each of elucidated responses likely conferred the relatively high salt tolerance of HS2, lack of some of common algal responses under high salinity stress could offer potential targets along with the identified responses when aiming to further enhance the robustness of HS2 as an industrial algal crop. First, violaxanthin deepoxidase (VDE) and zeaxanthin epoxidase (ZEP) that are, respectively, involved in the synthesis of zeaxanthin and violaxanthin were not differentially expressed in HS2 under high salinity stress. Zeaxanthin, however, is known to be associated with several types of photoprotection events of the PSII reaction center (Dall'Osto et al., 2012); therefore, VDE upregulation has been acknowledged as one of common algal responses under high oxidative stress . Given that the relative amount of carotenoid pigments in HS2 was increased under high salinity stress (Yun et al., 2019), enhancing the content of zeaxanthin by either upregulating VDE or downregulating ZEP may further enhance the halotolerance of HS2. Furthermore, although NPQ was not changed under high salinity stress, the elevation of NPQ has been denoted as one of the common algal responses under stress conditions (Cui et al., 2017). It would be, therefore, interesting to modulate the NPQ activity of HS2 as part of an effort to confer a greater halotolerance or induce a higher lipid productivity. As an example of the latter, reducing the expression levels of peripheral light-harvesting antenna proteins in PSII was demonstrated to decrease NPQ of Chlorella vulgaris, thereby improving biomass productivity by funneling more photosynthetic energy toward the electron transport chain (Shin et al., 2016). A similar approach can be adapted to direct more light energy toward the electron transport chain and/ or to possibly increase the available NADPH pool, although cautions should be taken to avoid the possibility of antagonistic interactions between competing metabolic pathways.