Identification of a master transcription factor and a regulatory mechanism for desiccation tolerance in the anhydrobiotic cell line Pv11.

Water is essential for living organisms. Terrestrial organisms are incessantly exposed to the stress of losing water, desiccation stress. Avoiding the mortality caused by desiccation stress, many organisms acquired molecular mechanisms to tolerate desiccation. Larvae of the African midge, Polypedilum vanderplanki, and its embryonic cell line Pv11 tolerate desiccation stress by entering an ametabolic state, anhydrobiosis, and return to active life after rehydration. The genes related to desiccation tolerance have been comprehensively analyzed, but transcriptional regulatory mechanisms to induce these genes after desiccation or rehydration remain unclear. Here, we comprehensively analyzed the gene regulatory network in Pv11 cells and compared it with that of Drosophila melanogaster, a desiccation sensitive species. We demonstrated that nuclear transcription factor Y subunit gamma-like, which is important for drought stress tolerance in plants, and its transcriptional regulation of downstream positive feedback loops have a pivotal role in regulating various anhydrobiosis-related genes. This study provides an initial insight into the systemic mechanism of desiccation tolerance.

Introduction

Living organisms require water for active life. Losing water, desiccation stress, is lethal for most organisms. Terrestrial organisms are particularly vulnerable to desiccation stress, and have evolved mechanisms to tolerate it. Entering an ametabolic state (anhydrobiosis) enables organisms to survive even if more than 99% of body water is lost, to recover after rehydration, and to resume their normal life cycle [1]. Rotifers, tardigrades, nematodes, plants, and larvae of the African midge Polypedilum vanderplanki acquired the molecular mechanisms underlying anhydrobiosis and can avoid mortality caused by desiccation [2-6]. An embryonic cell line established from P. vanderplanki, Pv11, is desiccation tolerant and is a tool to study anhydrobiosis [7]. After pretreatment with culture medium containing 600 mM trehalose, Pv11 cells survive 12 days of desiccation (<10% relative humidity) and resume proliferation a few days after rehydration [8]. Thus, the mechanism of the desiccation tolerance of P. vanderplanki is at the cellular level.

After the onset of desiccation stress, P. vanderplanki begins to synthesize and accumulate trehalose, a nonreducing disaccharide C12H22O11 [9]. Trehalose is a compatible solute; according to the water replacement hypothesis, it protects phospholipid membranes and intracellular macromolecules [10]. Initially, trehalose accumulation is crucial for preventing lethal damage caused by desiccation. Then, genes encoding antioxidant proteins, late embryogenesis abundant (LEA) proteins, and protein-L-isoaspartate (D-aspartate) O-methyltransferases (PIMTs), which reduce the harmful effects of desiccation, are upregulated [11-15]. Antioxidant proteins reduce the harmful effects of reactive oxygen species (ROS) generated during desiccation and rehydration. LEA and PIMT proteins prevent protein denaturation caused by desiccation. The above genes are also upregulated in Pv11 cells after pretreatment with a high concentration of trehalose; thus, trehalose could play a role in the initiation of the transcription of anhydrobiosis-related genes [16]. Genes for oxygen-binding hemoglobins and aquaporins, which controls osmotic pressure on the phospholipid membrane, are upregulated [11, 17]. The genes related to the ribosome complex are downregulated, reducing energy consumption for protein translation; in desiccated Pv11 cells, no significant changes in gene expression have been detected and the expression of each gene is steady [16]. The addition of water to dried Pv11 cells significantly upregulates the genes related to DNA repair, especially to homologous recombination and nucleotide excision repair [16, 18]. These genes are needed not only to repair DNA damage that occurred during desiccation despite the activity of anhydrobiosis-related protective genes, but also to help resume proliferation.

Although various genes related to desiccation tolerance have been detected and molecular mechanisms for the anhydrobiosis of P. vanderplanki have been inferred [9-18], the transcription factors and regulatory mechanisms that drive the expression of these genes have not been clarified. Mazin et al. revealed that heat shock factor (HSF) has a substantial role in regulation of these genes [19]; the HSF-binding element is enriched in the enhancers of anhydrobiosis-related genes. In Pv11 cells, HSF knock-down significantly decreases transcription of the genes for the antioxidant protein thioredoxin (TRX), a LEA protein, and trehalose-6-phosphate phosphatase, but not that of other anhydrobiosis-related genes [19]. This knock-down does not completely suppress cell proliferation after rehydration, and some cells are able to recover. Thus, HSF is required but not sufficient to enhance the desiccation tolerance system. These data illustrate that transcription factors control anhydrobiosis and provide evidence of both regulation- and gene structure-level evolution in anhydrobiotic species.

Using Pv11 cells as a model, we conducted a transcription factor–oriented genome-wide study to clarify the anhydrobiosis transcriptional regulatory network. We showed that nuclear transcription factor Y subunit gamma-like (NF-YC), which is important for drought stress tolerance in plants [20], regulates various anhydrobiosis-related genes. We revealed that many positive feedback loops, which are important for the ON/OFF switching of the system [21], were regulated by NF-YC and could initiate the transcription of various anhydrobiosis-related genes, such as LEA, TRX, aquaporins, hemoglobins, as well as genes related to trehalose synthesis and DNA repair. We also showed that HSF is involved in the positive feedback loop regulation. To the best of our knowledge, this is the first report of the clarification of a comprehensive gene regulatory network for desiccation tolerance; we suggest that the driving mechanism to induce the anhydrobiosis-related genes consists mainly of NF-YC and many downstream positive feedback transcriptional loops. Thus, our study confirmed the role of changes in the gene regulatory network as part of changes leading to the establishment of anhydrobiosis in insects.

Results

Quantification of gene expression by RNA-seq and analysis of differentially expressed genes during trehalose pretreatment and rehydration

We first obtained samples of total RNA from Pv11 cells under control conditions (T0) and at different stages of anhydrobiosis: trehalose pretreatment for 12 h (T12), 24 h (T24), 36 h (T36), or 48 h (T48), subsequent desiccation for 10 days (R0), and rehydration of R0 cells for 3 h (R3), 12 h (R12), 24 h (R24), or 72 h (R72). Experiments were performed in biological triplicate. All extracted RNA reads were mapped to the scaffolds of the P. vanderplanki genome; the mapping rates were higher than 90% except for R3, R12, and R24 (S1 Table). Not only the nuclear genome but also the mitochondrial genome is transcribed in Pv11 cells [22]. Mapping the reads of all samples to both genomes (S1 Fig) revealed that the mapping rates to the mitochondrial genome were much higher in R3, R12, and R24. To quantify the expression of each gene in all samples, we counted the mapped reads in each gene region using the P. vanderplanki annotated genome and performed principal component analysis (PCA) (Fig 1). The variance was considerably lower within each sample than between samples. A smaller distance between two points in a PCA plot indicates a more similar pattern of overall expression between the corresponding samples. Based on this mathematical principle, this result indicated that the expression patterns were more similar between the biological replicates than between different samples. Therefore, the samples were adequate for the intended analysis.

Fig 1

Principal component analysis of all samples.

For each point on the plot, the label refers to the hours of trehalose pretreatment (T) or rehydration (R). The axes show the principal components (PCs). The contribution ratio was 23.6% for PC1, 14.4% for PC2, and 8.2% for PC3; the cumulative contribution ratio was 46.2%.

We performed differentially expressed gene (DEG) analysis and detected 9,558 DEGs for trehalose pretreatment and 11,997 for rehydration (likelihood ratio test, adjusted p-value < 0.05, Benjamini-Hochberg (BH) method, S1 Data). The number of DEGs annotated with GO:003700 (DNA-binding transcription factor activity) was 5 for trehalose pretreatment, 31 for rehydration, and 76 for both (Fig 2A). The number of other DEGs, which we considered as the target genes (such as LEA genes), was 1,249 for trehalose pretreatment, 3,662 for rehydration, and 8,228 for both (Fig 2B). The number of DEGs for both trehalose and rehydration treatments was much higher than that for each individual treatment.

Fig 2

Venn diagram of differentially expressed genes in trehalose pretreatment and rehydration.

(A) Transcription factor genes, (B) other genes (DESeq2, adjusted p-value < 0.05, Benjamini-Hochberg method).

Inference of transcriptional regulatory networks

Based on time-series expression data on differentially expressed transcription factors, we inferred the transcriptional regulatory networks activated during trehalose pretreatment and rehydration (S2 and S3 Figs). We then integrated these networks including the transcription factors and estimated regulations in both treatments (Fig 3). Feed-forward loop (FFL) and feedback loop (FBL) transcriptional regulation among any three transcription factors controls the characteristic dynamics of gene expression [21, 23]. FFLs can be classified into the coherent type (the same sign of direct and indirect regulation) and incoherent type (opposite signs of direct and indirect regulation) [24], and FBLs can be classified into positive and negative [23] (Fig 4A). The estimated network contained 11 coherent FFLs, 5 incoherent FFLs, 2 positive FBLs, and 2 negative FBLs (Fig 4B). In the randomization test, only the number of coherent FFLs was significantly higher (p-value = 0.0177) than that in a random network with the same number of nodes and edges. A coherent FFL is able to cut a transient input but transmit a sustained signal to the output [24], and acts as a noise filter. The fact that the transcriptional regulatory network in Pv11 cells contained a significant number of coherent FFLs implies that this desiccation tolerance system might accurately regulate anhydrobiosis-related genes.

Fig 3

Inferred transcriptional regulatory network obtained by integration of those of trehalose pretreatment and rehydration.

Each rectangular node refers to a transcription factor gene with the ID defined in the previous study [19]. Genes differentially expressed in trehalose pretreatment (yellow), in rehydration (blue), and both treatments (green) are shown. Red arrows, positive regulation; blue arrows, negative regulation.

Fig 4

Conceptual diagram of important structures in gene regulation in the integrated network.

(A) Structure of feed-forward loops (FFLs) and feedback loops (FBLs) [21], [23]. (B) The numbers of these loops in the integrated network. * Randomization test, p-value < 0.05.

Detection of groups of coexpressed genes (modules)

We detected 27 modules (Fig 5, S2 Data) and performed gene ontology (GO) enrichment analysis to identify the modules essential for desiccation tolerance (Fisher’s exact test, p-value < 0.05, Table 1, S3 Data). Various GOs corresponding to gene functions reportedly necessary for desiccation tolerance were enriched in many modules (S1 Text). To annotate the modules containing anhydrobiosis-related genes, we counted the genes for LEA, TRX, PIMT, hemoglobin, aquaporin, and genes related to synthesis of trehalose and DNA repair in each module (S2 Table, S4 Data). Based on these results, we defined modules that included the maximum number of genes of each family previously reported as anhydrobiosis-related.

Fig 5

Heatmap of the time-series Z-scores for all genes clustered in modules as detected by weighted gene coexpression network analysis.

Each module was assigned an easily recognizable color. Functional annotations (GOs) for each module are shown in Table 1, S2 Table and S3 Data.

Table 1

Number of genes in each module and representative enriched GOs (GOs with the lowest p-value in Fisher’s exact test).

Module	Number of Genes	Representative GO ID	GO Term	p-value
Black	630	GO:0046034	ATP metabolic process	8.69 × 10−6
Blue	1222	GO:0006260	DNA replication	1.20 × 10−12
Brown	1045	GO:0005488	binding	1.34 × 10−13
Cyan	371	GO:0016310	phosphorylation	2.96 × 10−5
Dark green	255	GO:0098800	inner mitochondrial membrane protein complex	1.11 × 10−12
Dark grey	225	GO:1901565	organonitrogen compound catabolic process	4.98 × 10−4
Dark orange	213	GO:0042254	ribosome biogenesis	3.62 × 10−6
Dark red	265	GO:0004719	protein-L-isoaspartate (D-aspartate) O-methyltransferase activity	8.38 × 10−12
Dark turquoise	227	GO:0009055	electron transfer activity	1.56 × 10−4
Green	814	GO:0008417	fucosyltransferase activity	2.17 × 10−4
Green yellow	437	GO:0042592	homeostatic process	5.11 × 10−4
Grey	2	No Hits
Grey60	296	GO:0016667	oxidoreductase activity, acting on a sulfur group of donors	2.14 × 10−4
Light cyan	298	GO:0044444	cytoplasmic part	3.03 × 10−8
Light green	294	GO:0043227	membrane-bounded organelle	1.72 × 10−8
Light yellow	289	GO:0003735	structural constituent of ribosome	5.71 × 10−98
Magenta	552	GO:0140098	catalytic activity, acting on RNA	5.40 × 10−10
Midnight blue	340	GO:0016149	translation release factor activity, codon specific	5.04 × 10−4
Orange	217	GO:0008191	metalloendopeptidase inhibitor activity	9.04 × 10−5
Pink	578	GO:0005488	binding	1.35 × 10−11
Purple	517	GO:0004645	phosphorylase activity	1.71 × 10−3
Red	760	GO:0000502	proteasome complex	5.17 × 10−24
		GO:1905369	endopeptidase complex	5.17 × 10−24
Royal blue	279	GO:0036459	thiol-dependent ubiquitinyl hydrolase activity	2.40 × 10−4
		GO:0101005	ubiquitinyl hydrolase activity	2.40 × 10−4
		GO:0019783	ubiquitin-like protein-specific protease activity	2.40 × 10−4
Salmon	379	GO:0000786	nucleosome	3.56 × 10−13
		GO:0032993	protein-DNA complex	3.56 × 10−13
Tan	381	GO:0140097	catalytic activity, acting on DNA	9.47 × 10−5
Turquoise	1366	GO:0005488	binding	1.43 × 10−9
Yellow	887	GO:0044281	small molecule metabolic process	1.90 × 10−4

Integration of transcriptional regulatory network and modules

We performed motif enrichment analysis to detect the transcription factor–binding elements enriched in the upstream regions of the genes in each module (CLOVER [25], randomization test, p-value < 0.05; S5 Data). We also tested Granger causality between each transcription factor and the genes in the module at each time point (Granger causality test, adjusted p-value < 0.05, BH method; S5 Data). Finally, we extracted the regulatory relationships satisfying both criteria (Fig 6, Table 2, S5 Data). The results indicated that HSF (PvG018882) regulates all modules regardless of whether the regulation is direct or indirect. Thus, HSF was confirmed as a regulator of anhydrobiosis-related genes, consistent with the previous data [19].

Fig 6

Inferred gene regulatory network with the integration of the regulatory relationships of transcription factors and modules.

Rectangular nodes, transcription factors; circular nodes, modules; arrows, inferred regulatory relationships. The colors of rectangular nodes and arrows between transcription factors are the same as in Fig 3. The colors of the circular nodes are the same as in Fig 5. Regulatory relationships between transcription factors and modules are shown in grey (CLOVER, p-value < 0.05 and Granger causality test, adjusted p-value < 0.05, Benjamini-Hochberg method).

Table 2

Direct regulation of modules by transcription factors and dependency on heat shock factor (HSF, PvG018882) for each module.

Module	Direct Regulators	HSF(PvG018882) Dependency
Black	PvG029178; PvG034773; PvG029323; PvG031138	TRUE
Blue	PvG029323	TRUE
Brown	PvG029178; PvG006426; PvG000336; PvG013747	TRUE
Cyan	PvG030844; PvG021633; PvG033380; PvG004799;	TRUE
Dark green	PvG006426; PvG030844; PvG029323; PvG025794; PvG005038	TRUE
Dark grey	PvG006426; PvG029323	TRUE
Dark orange	PvG034773; PvG003458; PvG005038	TRUE
Dark red	PvG030844; PvG004799; PvG005038	TRUE
Dark turquoise	PvG013738; PvG034773; PvG003006; PvG029323	TRUE
Green	PvG029178	TRUE
Green yellow	PvG029178; PvG034773; PvG013747; PvG005038	TRUE
Grey	PvG027886; PvG021123;	TRUE
Grey60	PvG027886; PvG030844; PvG021633; PvG029323; PvG025794; PvG005530	TRUE
Light cyan	PvG030844; PvG021633; PvG003851; PvG034773; PvG029323	TRUE
Light green	PvG022603; PvG029323; PvG004799	TRUE
Light yellow	PvG034773; PvG004799; PvG013747; PvG003458	TRUE
Magenta	PvG021123; PvG000336; PvG004799	TRUE
Midnight blue	PvG025913; PvG018882; PvG017947; PvG033843; PvG029323	TRUE
Orange	PvG000336; PvG031138;	TRUE
Pink	PvG030844; PvG034773; PvG022743; PvG013747	TRUE
Purple	PvG027886; PvG029178; PvG021633; PvG013738; PvG031138	TRUE
Red	PvG022603; PvG030844; PvG033843; PvG022743; PvG013747; PvG025794; PvG005038	TRUE
Royal blue	PvG025913; PvG003854;	TRUE
Salmon	PvG033843; PvG022743; PvG005038	TRUE
Tan	PvG029178	TRUE
Turquoise	PvG025913; PvG027886; PvG029178; PvG030844; PvG012246; PvG000336; PvG003006; PvG004799; PvG025794; PvG003458; PvG005038	TRUE
Yellow	PvG025913; PvG030844; PvG022743; PvG013747; PvG005038	TRUE

Comparative analysis of gene regulatory network between Pv11 cells and Drosophila melanogaster, a desiccation-sensitive species

We expected the estimated gene regulatory network to include a regulatory mechanism specific to desiccation tolerance in anhydrobiotes. We performed a comparative analysis between the estimated network and the transcriptional regulatory network of D. melanogaster, which is a desiccation-sensitive species [26] (S4 Fig), yet the results were too complex to understand which structure of gene regulation is crucial for desiccation tolerance. Therefore, we contracted the network, mainly by extracting the transcription factors that have no homologues in the reported transcriptional regulatory network of D. melanogaster [26] and that directly regulated the modules (Fig 7). The contracted network was consistent with the previous reports [19] that HSF (PvG018882) is an important regulator of desiccation tolerance in P. vanderplanki. In the contracted network, PvG003116 was the most upstream transcription factor; it was homologous to nuclear transcription factor Y subunit gamma-like (NF-YC, accession number XP_029711307.1, blastp e-value = 1.0e-75). Most of the previously reported anhydrobiosis-related genes were regulated by PvG003116. The contracted network contained 3 FFLs and an FFL that included 3 FFLs, which we called TripleFFL (Fig 8). Four positive FBLs were downstream of NF-YC: between TripleFFL and FFL3; between FFL2 and FFL3; among TripleFFL, FFL2, and FFL3; and among HSF, FFL2, and TripleFFL. The modules that included the genes for LEA, TRX, aquaporin, hemoglobin, and the genes related to trehalose synthesis and DNA repair were regulated by NF-YC and these FBLs.

Fig 7

Contracted network obtained from the comparative analysis of the transcriptional regulatory networks between D. melanogaster and Pv11 cells.

Rectangular nodes, transcription factors; circular nodes, modules; arrows, inferred regulatory relationships. Transcription factors with sequence similarity to the transcription factors in the D. melanogaster transcriptional regulatory network are shown in green (blastp, e-value < 1.0e-15); other transcription factors are shown in red. Regulatory relationships confirmed in both Pv11 cells and D. melanogaster are shown as green arrows, and those specific for Pv11 cells are shown as red arrows. Transcription factors that are components of feed-forward loops (FFLs) are shown in blocks FFL1, FFL2, FFL3, and TripleFFL.

Fig 8

Final contracted gene regulatory network for desiccation tolerance in Pv11 cells.

The network includes four positive FBLs: between TripleFFL and FFL3; between FFL2 and FFL3; among TripleFFL, FFL2, and FFL3; and among HSF, TripleFFL, and FFL2. Modules that include the anhydrobiosis-related genes for LEA, thioredoxin (TRX), hemoglobin (HB), aquaporin (AQP), and the genes related to the synthesis of trehalose and DNA repair are regulated by NF-YC homologue through the 4 positive FBLs.

Discussion

In our new RNA-seq datasets, the mapping rates of the R3, R12, and R24 samples were considerably lower for the nuclear genome and were much higher for the mitochondrial genome than those of other samples (S1 Fig). These samples were obtained soon after the onset of rehydration after complete desiccation. Considering that mitochondria generate energy, we infer that the mechanism of recovery of Pv11 cells includes mitochondrial activation during rehydration. Currently, there are no reports of the active transcription of the mitochondrial genome just after rehydration. Subsequent detailed analysis of the process to activate the mitochondrial transcription may contribute to the discovery of the active recovery mechanism after rehydration.

In this study, the result of GO enrichment analysis for each module was not used for any prior knowledge to extract the finally contracted gene regulatory network (Fig 8). This result was used only for the verification of whether GOs related to the functions previously revealed as important for desiccation tolerance are enriched in modules or not (S1 Text). This result showed that almost all of the important GOs for desiccation tolerance was enriched in any of the modules. However, there are various GOs that were not focused on this verification (S3 Data), and the relationships between these GOs and desiccation tolerance was not experimentally clarified yet. The subsequent detailed analysis of these GOs could lead to the discovery of novel molecular biological mechanisms of desiccation tolerance.

Our analysis indicates that NF-YC is the most upstream transcription factor regulating various anhydrobiosis-related genes. It is one of the three subunits of nuclear factor Y (NF-Y) [27], the other two being NF-YA and NF-YB. NF-YC binds to the CCAAT motif through the interaction with NF-YA and NF-YB and regulates the transcription of various genes [28]. Transgenic Arabidopsis thaliana constitutively expressing NF-YB was alive even after 8 days of severe drought stress, whereas wild-type plants were completely wilted [29]. The overexpression of NF-YC in rice enables it to stay turgid after 10 days of water withholding, whereas the wild-type rice showed extreme wilting [20]. These studies indicate that NF-Y is essential for drought stress tolerance in plants, but there have been no reports of the relationship between NF-Y and desiccation tolerance in animals. To our knowledge, this is the first report implicating NF-Y in the extreme desiccation tolerance of an anhydrobiotic animal. We examined NF-YC expression at the organismal level by re-analyzing our previously reported dataset [14]. The expression increased after desiccation and rehydration in P. vanderplanki larvae (24 h of desiccation/control ratio = 1.22; 3 h of rehydration/48 h of desiccation ratio = 7.22). Polypedilum nubifer is closely related to P. vanderplanki but is desiccation sensitive; it has an NF-YC homologue, Pn.05756 (blastn, e-value = 9e-32), and the increase in its expression was similar to that in P. vanderplanki (24 h of desiccation/control = 1.21) [14]. In D. melanogaster, nuclear factor Y-box C (FBtr0070974) is homologous to P. vanderplanki NF-YC (tblastn, e-value = 1.35e-38), but there are no reports about its expression under desiccation stress. We infer that the NF-YC gene is conserved in insects and contributes to mild desiccation tolerance; in P. vanderplanki, NF-YC has evolved to strengthen desiccation tolerance.

Why is P. nubifer not desiccation tolerant despite an increase in the expression of NF-YC similar to that of P. vanderplanki? Our estimated network could provide an answer to this question: the transcription factors connecting NF-YC to the modules that includes anhydrobiosis-related genes are absent in P. nubifer. These transcription factors include (1) PvG02163 in FFL2, which regulates FFL1 and modules containing genes for LEA, aquaporin, and hemoglobin, and trehalose-related genes; (2) PvG024515 in FFL2, which regulates FFL1 and TRX and DNA repair modules; and (3) PvG029323 in TripleFFL (Fig 7), which regulates modules containing genes for LEA, aquaporin, and hemoglobin, and trehalose- and DNA repair–related modules (blastn, e-value < 1.0e-5) [14]. This result implies that P. nubifer is not desiccation tolerant because of the absence not only of anhydrobiosis-related genes, such as LEA genes, but also of critical transcription factor genes regulating anhydrobiosis-related genes.

Our estimated gene regulatory network includes several positive FBLs from NF-YC to modules that include anhydrobiosis-related genes (Fig 8). A positive FBL of transcriptional regulation controls the activation (ON) or inactivation (OFF) of the system depending on the strength of the input signal [21]. Four positive FBLs in the estimated gene regulatory network imply that the transcription of anhydrobiosis-related genes is activated by trehalose pretreatment but not in the normal hydrated state (Fig 9).

Fig 9

Conceptual diagram of the inferred ON/OFF mechanism of desiccation tolerance in Pv11 cells.

The inferred gene regulatory network for desiccation tolerance contains four positive FBLs (Fig 8). A positive FBL of transcriptional regulation can either activate (ON) or inactivate (OFF) the system depending on the strength of the input signal [21]. This property and our results suggest that the activation or inactivation of anhydrobiosis-related genes may be controlled by these positive FBLs.

Our estimated gene regulatory network is sufficient to explain various previous data on the expression of anhydrobiosis-related genes. However, the essential and required structure of gene regulation in the estimated network has not been tested yet. A recent study used RNAi-mediated knock-down in Pv11 cells and achieved silencing of a specific gene [30]. This gene silencing technique would be useful to test which transcription factor in our estimated network, such as NF-YC or a component of a positive FBL, is essential for desiccation tolerance; a significant decrease in the survival rate of Pv11 cells would be expected after gene silencing of such a transcription factor. By inferring de novo a gene regulatory network, we discovered the putative induction mechanism of desiccation tolerance. The estimated network also gave us the insight into the acquisition of desiccation tolerance in evolution and may allow us to control desiccation tolerance in the future.

Materials and methods

RNA sample preparation, Illumina HiSeq-1500 sequencing, and data preprocessing

Pv11 cells were maintained in IPL-41 medium supplemented with 10% FBS (T0) and were incubated in 600 mM trehalose solution containing 10% of supplemented IPL-41 medium for 12, 24, 36, or 48 h (T12, T24, T36, T48). The cells were desiccated in a desiccator for 10 days (R0; relative humidity, <10%). After desiccation, they were rehydrated with IPL-41 containing 10% FBS for 3, 12, 24, or 72 h (R3, R12, R24, R72). Samples were prepared in biological triplicate. Total RNA was extracted with a TruSeq RNA sample prep kit v2 (Illumina, San Diego, CA, USA) and purified with a NucleoSpin RNA kit (Macherey-Nagel, Düren, Germany). RNA quality was assessed with a Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA), and RNA was quantified with a Qubit 2.0 fluorometer (Thermo Fisher, Waltham, MA, USA).

First-strand cDNAs were synthesized and sequenced on an Illumina HiSeq 1500 with a length of 51 nucleotides using a single-end sequencing protocol in the Rapid mode at the Functional Genomics Facility of the National Institute for Basic Biology (Okazaki, Japan). Read quality was checked with FastQC [31]. Based on the results of quality check, the reads were preprocessed by cutting the adaptor region and the first 20 bases from the 3’-end if these bases represented a poly(A) tail; FaQCs [32] was used for preprocessing. Reads with a total length less than 14 bp and reads with more than 2 ‘N’ bases were removed. The remaining reads were mapped to the scaffold of the P. vanderplanki genome (midgeBase: http://150.26.71.110/midgebase/) using HISAT2 [33]. The mapping rate was calculated as the ratio of mapped reads to total remaining reads. The mapped reads were counted using the annotated P. vanderplanki genome [19] with the qCount function in the QuasR package [34].

Principal component analysis

The read count data of each sample were transformed to reads per kilobase of exon per million mapped reads (RPKM) and used for PCA with the prcomp function in R [35]. The PC1–3 scores were calculated for each sample and plotted using the scatterplot3D package [36] in three-dimensional space.

Analysis of differentially expressed genes

DEG analysis was performed using raw count data of each sample with the nbinomLRT function in the DESeq2 package [37]. In the log-likelihood ratio test, the following hypothesis testing was performed for each gene by comparing the mean expression of regarding gene in the overall time-series of trehalose treatment or rehydration, and the expression of regarding gene in each time point.

H0: The mean expression of a gene is the same in overall time points (Null hypothesis).

H1: The mean expression of a gene at least at a one-time point was different from its mean expression in overall time points (Alternative hypothesis).

DEGs were detected when the false discovery rate of the one-sided test calculated by the BH method [38] was less than 0.05. DEGs were classified as those for transcription factors if annotated with GO:0003700 (sequence-specific DNA binding transcription factor activity) [39], or as others if they were not. The numbers of DEGs for transcription factors and other DEGs were visualized as a Venn diagram using the venn.diagram function in the VennDiagram package [40].

Inference of transcriptional regulatory network based on time-series expression data

The regulatory relationships between transcription factors detected as DEGs were inferred using time-series RPKM of each transcription factor in trehalose pretreatment or rehydration using BTNET, which currently has the best performance [41]. Gradient boosting in BTNET was used. The number of boosting trees was set to 100 and the time lag of BTNET was set to 1. Determination of the threshold of the obtained score for each regulation inferred by BTNET was based on the fact that the transcriptional regulatory networks of Escherichia coli and Saccharomyces cerevisiae have scale-free topology [42]. The scale-free property of the network was defined as where k denotes the degree (the number of regulations of the transcription factor), γ is the coefficient of scale-free topology, and C is the normalization coefficient. The above equation can be log-log transformed as follows:

Based on the transformed equation, the scale-free property of the network constructed by adopting the regulation whose score was more than score a was evaluated by linear regression between degree and degree distribution using the lm function in R. At this step, γ, the p-value of one-sided F-test, and the R2 value were calculated for a from 0 to 1 per estimated score (S5 and S6 Figs). Finally, the minimum score satisfying the condition that γ was negative and the p-value of the F-test was less than 0.05 was determined as the threshold of the estimated score by BTNET. The positive or negative regulations were classified based on the sign of the value of cross-correlation with a time lag of 1 using the ccf function in R, i.e., the regulation was considered positive when its value was more than 0 and negative when it was less than 0. The estimated transcriptional regulatory networks for trehalose pretreatment and rehydration were integrated to extract the transcription factors and regulation present at least in one of the estimated networks. The transcriptional regulatory networks were visualized using yEd [43]. FFLs and FBLs were extracted using the graph.get.subisomorphism.vf2 function in the igraph package [44]. The extracted FFLs were classified into coherent or incoherent when the signs of direct and indirect regulation were the same or opposite, respectively. The extracted FBLs were classified into positive or negative when the number of negative regulations was even or odd, respectively. Whether these numbers were significantly higher than that of a random network was tested using a randomization test as follows. A random network whose number of nodes and edges was exactly the same as the number of transcription factors and regulations in the integrated network was generated on the basis of the Erdös–Renyi model [45] using the erdos.renyi.game function in the igraph package. The sign of each regulation in the random network was based on a binomial distribution under the condition that the probability of success in the binomial distribution was equal to the ratio of positive regulations to that of all regulations in the integrated network. The numbers of coherent and incoherent FFLs and positive and negative FBLs in the random network were calculated as for the integrated network. These calculations were iterated 10,000 times and the networks with a higher number of each type of FFL and FBL than that of integrated network were counted (S7 Fig). Finally, the p-value was calculated for each FFL and FBL by dividing the counts for this structure by 10,000; p-values less than 0.05 were considered significant.

Detection of modules and gene ontology enrichment analysis

One transcription factor binds regulatory regions of various genes and induces the coordinated expression (coexpression) of these genes, called “gene batteries” [46], which we refer to as modules. To detect such modules, we performed weighted gene coexpression network analysis (WGCNA) [47]. For all DEGs other than those for transcription factors, a pair-wise correlation matrix was computed using time-series RPKM of each gene in trehalose pretreatment and rehydration. Using the pickSoftThreshold function in the WGCNA package, we chose the power of 7 to satisfy the scale-free topology criterion so that the signed R2 value of the fitting of the scale-free equation was more than 0.75 (S8 Fig). Then the adjacency matrix for each pair of genes was obtained by raising the correlation matrix to the power of 7. This adjacency matrix was transformed into a robust measure of network interconnectedness (topological overlap measure, TOM [48]). The dissimilarity of each pair of genes was calculated by subtracting 1 from each element of the TOM-based adjacency matrix. This dissimilarity matrix was used as an input for average linkage hierarchical clustering using hclust in R (S9 Fig).

Modules were detected with the cutreeDynamic function in the WGCNA package. This function requires setting the minimum number of genes per module, minClusterSize. We determined this value as the number that increased the pseudo-F [49], a clustering metric, from 1 to 200; the minClusterSize value based on this criterion was 165 (S10 Fig). Finally, 27 clusters were detected as modules. Each module was assigned an easily recognizable color. To clarify what functional genes belong to each module, we performed GO enrichment analysis using the GSEAGOHyperGParams function in the GOstats package [50], and GOs with p-values < 0.05 were detected as enriched in each module (S3 Data).

Prediction of the regulation of coexpressed genes in a module by a transcription factor

Transcription factor–binding motif enrichment analysis was performed with CLOVER [25]. The information on binding motif for each transcription factor of P. vanderplanki is not yet available; therefore, we collected information for D. melanogaster for each transcription factor from the CISBP [51], FlyFactorSurvey [52], and JASPER2018 [53] databases using motifDB package [54]; we collected 134, 291, and 140 binding motifs, respectively. Then we extracted amino acid sequences of the corresponding transcription factors and searched for sequence similarities against the whole genome of P. vanderplanki using bidirectional search in blastp (e-value < 1.0e-5). We collected 72 transcription factors of P. vanderplanki homologous to D. melanogaster proteins, which were linked to the binding motif information. Using these binding motifs and the upstream regions of genes in each module as an input to CLOVER, we performed motif enrichment analysis. The upstream regions were designated as 100, 500, 1,000, 2,000, 3,000, 4,000, 5,000, and 10,000 bp, and the scaffold of P. vanderplanki was used as a control sequence. A transcription factor was considered to bind when the p-value of one-sided test in CLOVER was less than 0.05 for at least one upstream region, but whether this binding induces gene transcription in the module was not tested. To take this into consideration, we tested Granger causality [55] between the transcription factor and the genes in the module in the time-series expression data for trehalose pretreatment or rehydration. In this test, we transformed all RPKM values for each DEG into Z-score and created a vector autoregressive (VAR) model, in which the time-series Z-score of a transcription factor was set as the cause and the mean time-series Z-score of genes in the module was set as the effect. The VAR model was created with the VAR function in the vars package [56], and the time lag of the model was optimized with the VARSelect function in the same package for each candidate transcription factor–module pair. We tested whether the variance from the above model was significantly smaller than that of the null model using one-sided F-test with the causality function in the vars package. The Granger causality was detected when the adjusted p-value of the F-test by the BH method was less than 0.05. The regulatory relationship was detected if it passed both tests, motif enrichment analysis and Granger causality test.

Comparative analysis of transcriptional regulatory networks between Pv11 cells and D. melanogaster and contraction of the estimated network

To reveal the part of the regulatory network specific to desiccation tolerance, we compared the transcriptional regulatory networks (nodes and edges) between Pv11 cells and D. melanogaster. That of D. melanogaster was obtained from modEncode (http://intermine.modencode.org/release-33/flyRegulatoryNetwork.do, network was shown in Fig 7A in [26]). To compare the nodes, we examined sequence similarity between P. vanderplanki and D. melanogaster transcription factors included in the modEncode network (blastp, e-value < 1.0e-15). To compare the edges, the regulatory relationships were judged to exist in P. vanderplanki when either direct or indirect regulation from the source to target transcription factors of D. melanogaster homologous to those of P. vanderplanki existed in the modEncode network. To extract the structure important for desiccation tolerance, we contracted the obtained network. We extracted the transcription factors that had no homology with the transcription factors in modEncode and that directly regulated the modules, and selected those that satisfied any of the following criteria: 1) the most upstream transcription factors, i.e., those with no indegree edge, 2) FFL components, and 3) HSF and transcription factors to connect the transcription factors extracted based on criteria 1) and 2). We selected the regulations that were judged to exist only in P. vanderplanki. Finally, we presented each FFL in the contracted network as one block and eliminated the transcription factors between HSF and FFLs.

Results of mapping for each sample.

(XLSX)

Click here for additional data file.

Modules with the maximum number of anhydrobiosis-related genes.

(XLSX)

Click here for additional data file.

DEGs in trehalose pretreatment and rehydration.

(XLSX)

Click here for additional data file.

List of P. vanderplanki genes in each module.

(XLSX)

Click here for additional data file.

GO terms enriched in each module.

(XLSX)

Click here for additional data file.

BLASTN results for the known anhydrobiosis-related genes annotated in a previous study [14] vs. genes of P. vanderplanki annotated in a previous study [19], and the module for each gene.

(XLSX)

Click here for additional data file.

Transcription factor–binding elements in P. vanderplanki genes included in modules detected by motif enrichment analysis (CLOVER) and Granger causality test.

(XLSX)

Click here for additional data file.

This document contains detailed information about the results of GO enrichment analysis for the module that includes anhydrobiosis-related genes reported previously.

(PDF)

Click here for additional data file.

Mapping rate to the chromosomal and mitochondrial genomes of P. vanderplanki for each sample.

Three biological replicates per sample are shown. The mapping rate to the mitochondrial genome increased sharply just after rehydration.

(PDF)

Click here for additional data file.

Transcriptional regulatory network inferred for trehalose pretreatment.

Each rectangular node refers to a transcription factor; transcription factor IDs are as defined in [19]; arrows show the inferred regulatory relationships (red, positive; blue, negative). The node color was based on the value of the mean Z-score in each sample calculated at each time point as reads per kilobase of exon per million mapped reads (RPKM) for each gene: red, high Z-score; blue, low Z-score.

(PDF)

Click here for additional data file.

Transcriptional regulatory network for rehydration.

Designations are as in S2 Fig.

(PDF)

Click here for additional data file.

Comparison of the transcriptional regulatory networks between D. melanogaster and Pv11 cells.

Rectangular nodes, transcription factors; circular nodes, modules. Transcription factors with sequence similarity to those in the D. melanogaster transcriptional regulatory network are shown in green (blastp, e-value < 1.0e-15); other transcription factors are shown in red. Regulatory relationships confirmed in both P. vanderplanki and D. melanogaster are shown as green arrows, and those specific for P. vanderplanki are shown as red arrows. Regulatory relationships between transcription factors and modules are shown as grey arrows (CLOVER, p-value < 0.05 and Granger causality test, adjusted p-value < 0.05, Benjamini-Hochberg method). Thus, the candidate regulatory network specific to desiccation tolerance is represented by red nodes and arrows.

(PDF)

Click here for additional data file.

Determination of the threshold for the score of inferred regulation by BTNET in trehalose pretreatment.

(upper left) The value of γ, which is the property of scale-free topology. (upper right) p-value for γ in F-test. (lower left) R2 values obtained by fitting linear relationships between degree (the number of regulations of the transcription factor) and degree distribution. Red lines show the minimum value of the threshold (0.160963) at which γ was negative and p-value was less than 0.05. (lower right) A plot of the log-log transformed degree distribution based on threshold = 0.160963. The line shows fitting result for the scale-free equation.

(PDF)

Click here for additional data file.

Determination of the threshold for the score of inferred regulation by BTNET in rehydration.

For panel descriptions and designations, see S5 Fig. The threshold value was 0.151133.

(PDF)

Click here for additional data file.

Results of the randomization test of the significant number of FFLs and FBLs in the inferred transcriptional regulatory network obtained by integration of those of trehalose pretreatment and rehydration.

Histograms show the number of networks generated by the Erdös–Renyi model (#) categorized against the number of specific structures in the network. Red columns show the numbers of generated networks obtaining more than that of the inferred transcriptional regulatory network obtained by integration of those of trehalose pretreatment and rehydration. FFL, feed-forward loop; FBL, feedback loop. Rates of the red columns: coherent FFLs, 0.0177; incoherent FFLs, 0.5331; positive FBLs, 0.3087; negative FBLs, 0.3009. Thus, only the number of coherent FFLs was significant in the integrated network at the significance level of α = 0.05.

(PDF)

Click here for additional data file.

Signed R2 value calculated by fitting the scale-free equation against the power value in weighted gene coexpression network analysis (WGCNA).

7 of power satisfied that the signed R2 was more than 0.75 (red line).

(PDF)

Click here for additional data file.

Dendrogram of dissimilarity matrix calculated by topological overlap measurements (TOM) subtracted by 1 in WGCNA analysis.

TOM was calculated from the Pearson correlation matrix of time-series Z-score for gene pairs raised by the estimated power of 7.

(PDF)

Click here for additional data file.

Change of the clustering metric, pseudo-F, for each minModuleSize of a dynamic tree cut from the dendrogram.

Pseudo-F increased from minModuleSize 1 to 165 and drastically decreased after 166; therefore, the minModuleSize was determined as 165.

(PDF)

Click here for additional data file.

9 Jan 2020

PONE-D-19-27695

Identification of a master transcription factor and a regulatory mechanism for desiccation tolerance in the anhydrobiotic cell line Pv11

PLOS ONE

Dear Dr. Funahashi,

Thank you for submitting your manuscript to PLOS ONE. After careful consideration, we feel that it has merit but does not fully meet PLOS ONE’s publication criteria as it currently stands. Therefore, we invite you to submit a revised version of the manuscript that addresses the points raised during the review process.

The manuscript has been evaluated by two reviewers, and their comments are available below. The reviewers have raised a few issues that they feel a revised submission should address. In particular, they feel that there are some parts of the data analyses that could be clarified. Please revise the manuscript to carefully address the concerns raised.

We would appreciate receiving your revised manuscript by Feb 23 2020 11:59PM. When you are ready to submit your revision, log on to https://www.editorialmanager.com/pone/ and select the 'Submissions Needing Revision' folder to locate your manuscript file.

If you would like to make changes to your financial disclosure, please include your updated statement in your cover letter.

To enhance the reproducibility of your results, we recommend that if applicable you deposit your laboratory protocols in protocols.io, where a protocol can be assigned its own identifier (DOI) such that it can be cited independently in the future. For instructions see: http://journals.plos.org/plosone/s/submission-guidelines#loc-laboratory-protocols

Please include the following items when submitting your revised manuscript:

A rebuttal letter that responds to each point raised by the academic editor and reviewer(s). This letter should be uploaded as separate file and labeled 'Response to Reviewers'.

A marked-up copy of your manuscript that highlights changes made to the original version. This file should be uploaded as separate file and labeled 'Revised Manuscript with Track Changes'.

An unmarked version of your revised paper without tracked changes. This file should be uploaded as separate file and labeled 'Manuscript'.

Please note while forming your response, if your article is accepted, you may have the opportunity to make the peer review history publicly available. The record will include editor decision letters (with reviews) and your responses to reviewer comments. If eligible, we will contact you to opt in or out.

We look forward to receiving your revised manuscript.

Kind regards,

Hanna Landenmark

Associate Editor

PLOS ONE

Journal Requirements:

When submitting your revision, we need you to address these additional requirements:

1. Please ensure that your manuscript meets PLOS ONE's style requirements, including those for file naming. The PLOS ONE style templates can be found at http://www.plosone.org/attachments/PLOSOne_formatting_sample_main_body.pdf and http://www.plosone.org/attachments/PLOSOne_formatting_sample_title_authors_affiliations.pdf

2. Please provide details of the origin and vendor information for the cell line used.

3. We note that you have indicated that data from this study are available upon request. PLOS only allows data to be available upon request if there are legal or ethical restrictions on sharing data publicly. For more information on unacceptable data access restrictions, please see http://journals.plos.org/plosone/s/data-availability#loc-unacceptable-data-access-restrictions.

In your revised cover letter, please address the following prompts:

a) If there are ethical or legal restrictions on sharing a de-identified data set, please explain them in detail (e.g., data contain potentially sensitive information, data are owned by a third-party organization, etc.) and who has imposed them (e.g., an ethics committee). Please also provide contact information for a data access committee, ethics committee, or other institutional body to which data requests may be sent.

b) If there are no restrictions, please upload the minimal anonymized data set necessary to replicate your study findings as either Supporting Information files or to a stable, public repository and provide us with the relevant URLs, DOIs, or accession numbers. For a list of acceptable repositories, please see http://journals.plos.org/plosone/s/data-availability#loc-recommended-repositories.

We will update your Data Availability statement on your behalf to reflect the information you provide.

[Note: HTML markup is below. Please do not edit.]

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

Reviewer #2: Yes

**********

2. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: Yes

Reviewer #2: Yes

**********

3. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

Reviewer #2: Yes

**********

4. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

**********

5. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: The anhydrobiosis is an outstanding phenomenon in animal and plant. Nevertheless, many efforts were invested to analyze the genetic basis of desiccation tolerance, transcriptional regulatory mechanisms underlying induction of desiccation-associated genes are still poorly studied. The authors thoroughly analyzed gene regulatory network in embryonic cell line Pv11 of African midge, Polypedilum vanderplanki, and compared the obtained results with that of Drosophila melanogaster. The main goals are clearly enunciated. Research design and methodologies in addressing the raised questions are very well thought out, appropriate and fitting. A huge amount of strong and convincing data and deep and wide discussion of desiccation tolerance mechanisms are presented in the paper. The experimental procedures are described so well and scrupulous that no questions remain unanswered. The style of writing makes it easily understood by a diverse audience of scientists. Data presentation looks attractive, elegant and convincing. Figures and tables are clear and sufficiently illustrative. The suggested approach and obtained results will certainly forward understanding of desiccation tolerance mechanisms.

I strongly support publication.

Reviewer #2: This study uses a cell line from a highly desiccation tolerant insect to investigate gene regulatory networks underlying anhydrobiosis. Nuclear transcription factor Y subunit gamma-like (NF-YC) is identified as a key upstream regulator of feedback loops regulating anhydrobiosis-related genes. NF-YC is known to be important in plant drought tolerance; this study is an important addition to the list of common mechanisms (e.g. LEA proteins, antioxidant processes) responsible for dehydration-tolerance in plants and animals.

As a physiologist, I must state at the start that I am not very familiar with systems biology and have only recently been collaborating with a systems biologist. I outlined the general analysis trajectory of this study to her, and she felt (without knowing details) that the approach seemed reasonable. I will note, except for the section on Inference of transcriptional regulatory network…… (page 16/28), I found the Materials and methods to follow a clear logical sequence.

Experimental methods for RNA isolation and sequencing were standard. The differential expression analysis was a bit unclear. Fold-change differences are provided for each treatment. Differences compared to what? Untreated cells? Are these the data used to generate S2_Fig and S3_Fig?

The GO enrichment analysis (S1_Text) was restricted to genes previously identified as being involved in anhydrobiosis. This seems to be biased to confirm previous work. There are many unbiased GO enrichment methods published – why not look for novel GO categories that could lead to new understanding?

I like the comparison to the Drosophila melanogaster transcriptional regulatory network, but I was left wondering: which D. melanogaster network? Citation 26 includes one in the full text and 4 in the supplemental material.

The Discussion was concise and clearly placed the results of this work in the greater context of anhydrobiosis. In general, I found the entire manuscript to be well written and well edited.

**********

6. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: No

[NOTE: If reviewer comments were submitted as an attachment file, they will be attached to this email and accessible via the submission site. Please log into your account, locate the manuscript record, and check for the action link "View Attachments". If this link does not appear, there are no attachment files to be viewed.]

While revising your submission, please upload your figure files to the Preflight Analysis and Conversion Engine (PACE) digital diagnostic tool, https://pacev2.apexcovantage.com/. PACE helps ensure that figures meet PLOS requirements. To use PACE, you must first register as a user. Registration is free. Then, login and navigate to the UPLOAD tab, where you will find detailed instructions on how to use the tool. If you encounter any issues or have any questions when using PACE, please email us at figures@plos.org. Please note that Supporting Information files do not need this step.

15 Jan 2020

Reply to Reviewer #1

Thank you so much for your positive comments.

Reply to Reviewer #2

Thank you for your precious comments. We hope that our revised manuscript will fulfill your request.

Minor comments

1. Experimental methods for RNA isolation and sequencing were standard. The differential expression analysis was a bit unclear. Fold-change differences are provided for each treatment. Differences compared to what? Untreated cells? Are these the data used to generate S2_Fig and S3_Fig?

Response:

You are correct that the explanation about DEG analysis was unclear in terms of the sample comparison. In this study, the detection of DEGs was performed by comparing the mean expression of a gene in the overall time-series of trehalose treatment or rehydration, and its expression in each time point. That means a gene is detected as DEG when the mean expression of the gene at least at a one-time point was different from the mean expression in overall time points. Therefore, we did not compare the expression of a gene with that from untreated cells. In S2_Fig and S3_Fig, there are nodes (transcription factors) that were detected as DEG based on the above DEG analysis. Clarifying these points, we have added the sentence to our manuscript in materials and methods "Analysis of differentially expressed genes'' (pp. 16, line 264 to 271).

2. The GO enrichment analysis (S1_Text) was restricted to genes previously identified as being involved in anhydrobiosis. This seems to be biased to confirm previous work. There are many unbiased GO enrichment methods published – why not look for novel GO categories that could lead to new understanding?

Response:

Firstly, we declare that the result of GO enrichment analysis for each module was not used for any prior knowledge to extract the finally contracted gene regulatory network (Fig. 8). This result was used only for the verification of whether GOs related to the functions previously revealed as important for desiccation tolerance are enriched in modules or not (S1 Text). As you pointed out, there are various GOs that were not focused on this verification (S3 Data). The reason why we did not analyze it in detail is that we mainly focused on finding not the mechanism regulated by effector genes but the transcriptional regulatory mechanism in this study. However, your point is undoubtedly right in terms of uncovering the mechanism regulated by effector genes. Therefore, we have added the sentence that "The subsequent detailed analysis of these GOs could lead to the discovery of novel molecular biological mechanisms of desiccation tolerance. '' and aforementioned scope of this research in Discussion (pp. 12, line 167 to 176).

3. I like the comparison to the Drosophila melanogaster transcriptional regulatory network, but I was left wondering: which D. melanogaster network? Citation 26 includes one in the full text and 4 in the supplemental material.

Response:

Your point is correct. The explanation of the resource about the transcriptional regulatory network of D.melanogaster was not enough. We used Fig. 7A in the full text of citation 26. We downloaded it from modEncode database, and the URL is http://intermine.modencode.org/release-33/flyRegulatoryNetwork.do. The main reason to use it was that this network was the most comprehensive transcriptional regulatory network compared with four networks in supplementary materials of citation 26. We have added this information in materials and methods "Comparative analysis of transcriptional regulatory networks between Pv11 cells and D. melanogaster and contraction of the estimated network'' (pp. 20, line 386 to 387).

25 Feb 2020

Identification of a master transcription factor and a regulatory mechanism for desiccation tolerance in the anhydrobiotic cell line Pv11

PONE-D-19-27695R1

Dear Dr. Funahashi,

We are pleased to inform you that your manuscript has been judged scientifically suitable for publication and will be formally accepted for publication once it complies with all outstanding technical requirements.

Within one week, you will receive an e-mail containing information on the amendments required prior to publication. When all required modifications have been addressed, you will receive a formal acceptance letter and your manuscript will proceed to our production department and be scheduled for publication.

Shortly after the formal acceptance letter is sent, an invoice for payment will follow. To ensure an efficient production and billing process, please log into Editorial Manager at https://www.editorialmanager.com/pone/, click the "Update My Information" link at the top of the page, and update your user information. If you have any billing related questions, please contact our Author Billing department directly at authorbilling@plos.org.

If your institution or institutions have a press office, please notify them about your upcoming paper to enable them to help maximize its impact. If they will be preparing press materials for this manuscript, you must inform our press team as soon as possible and no later than 48 hours after receiving the formal acceptance. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information, please contact onepress@plos.org.

With kind regards,

Axel Imhof

Academic Editor

PLOS ONE

Additional Editor Comments (optional):

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. If the authors have adequately addressed your comments raised in a previous round of review and you feel that this manuscript is now acceptable for publication, you may indicate that here to bypass the “Comments to the Author” section, enter your conflict of interest statement in the “Confidential to Editor” section, and submit your "Accept" recommendation.

Reviewer #1: All comments have been addressed

Reviewer #2: All comments have been addressed

**********

2. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

Reviewer #2: Yes

**********

3. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: Yes

Reviewer #2: Yes

**********

4. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

Reviewer #2: Yes

**********

5. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

**********

6. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: (No Response)

Reviewer #2: (No Response)

**********

7. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: No

4 Mar 2020

PONE-D-19-27695R1

Identification of a master transcription factor and a regulatory mechanism for desiccation tolerance in the anhydrobiotic cell line Pv11

Dear Dr. Funahashi:

I am pleased to inform you that your manuscript has been deemed suitable for publication in PLOS ONE. Congratulations! Your manuscript is now with our production department.

If your institution or institutions have a press office, please notify them about your upcoming paper at this point, to enable them to help maximize its impact. If they will be preparing press materials for this manuscript, please inform our press team within the next 48 hours. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information please contact onepress@plos.org.

For any other questions or concerns, please email plosone@plos.org.

Thank you for submitting your work to PLOS ONE.

With kind regards,

PLOS ONE Editorial Office Staff

on behalf of

Dr. Axel Imhof

Academic Editor

PLOS ONE