Illumina RNA sequencing and transcriptome assembly
For samples from experiment I, we extracted total RNA following the TRIZOL RNA extraction protocol (Chomczynski, 1993) and checked RNA quality with a Bioanalyzer 2100 (Agilent, Japan). RNA samples with RIN>7 were sent to the DNA Sequencing Section at the Okinawa Institute of Science and Technology (OIST) for library construction and Illumina Hiseq 4000 150 × 150-bp paired-end RNA sequencing. To remove adapters and low-quality sequences, r aw reads from Illumina RNA sequencing were trimmed with Trimmomatic v0.36 (Bolger, Lohse, & Usadel, 2014) and were checked for quality using FastQC (Andrews, 2010). To remove RNA reads from Symbiodinium and other coral-associated microbiomes, we first mapped trimmed reads (~877 million read pairs from all sequenced libraries) to a reference Pocillopora damicornis genome assembly (GenBank assembly accession: GCA_003704095.1) (Cunning et al., 2018) using Tophat v2.1.1 (Kim et al., 2013). Successfully aligned RNA reads (~168 million read pairs; successful alignment rate: 19.1%) were then subjected to de novo assembly using Trinity v2.8.4 (Grabherr et al., 2011). RNA read quantification was conducted using RSEM v1.3.2 (Li & Dewey, 2011) and reconstructed transcripts shorter than 200 nt or with low coverage (<5 transcripts per million) were discarded. BUSCO analysis was performed using the metazoa_odb9 dataset to examine completeness of the transcriptome assembly (Simão, Waterhouse, Ioannidis, Kriventseva, & Zdobnov, 2015).
Functional annotation, differential expression (DE), and gene ontology (GO) enrichment analyses
We performed functional annotation of the transcriptome assembly by searching the SwissProt eukaryotic protein database (downloaded on 3rd Jul 2019) with blastx, using a threshold of E=10-5. Although principle component analysis (PCA) of the annotated transcripts showed some variation among biological replicates (PC2: 14.9%), most of the variation could be attributed to the hyperosmotic treatment (PC1: 57.9%; Fig. 1). We therefore conducted DE analysis on annotated transcripts for pairwise comparisons of S43/S35 (representing response at 12 h) and S46/S35 (response at 24 h), designed to cancel out intra-individual variation and to identify transcriptomic response induced specifically by the treatment. The DE analysis was performed using edgeR (Robinson, McCarthy, & Smyth, 2010) with the criteria of >2-fold absolute change and >5 CPM (counts per million) for at least 3 libraries in a given pairwise comparison. GO enrichment analyses were then performed separately for upregulated and downregulated DEGs at 12 h (S43/S35) and 24 h (S46/S35), using DAVID bioinformatics resources v6.8 (Huang, Sherman, & Lempicki, 2008, 2009) to examine the participation of specific cellular processes and signaling pathways in polyp bail-out.