Skip to main content
Download PDF

Gene expression associated with white syndromes in a reef building coral, Acropora hyacinthus

BMC Genomics volume 16, Article number: 371 (2015) Cite this article

Abstract

Background

Corals are capable of launching diverse immune defenses at the site of direct contact with pathogens, but the molecular mechanisms of this activity and the colony-wide effects of such stressors remain poorly understood. Here we compared gene expression profiles in eight healthy Acropora hyacinthus colonies against eight colonies exhibiting tissue loss commonly associated with white syndromes, all collected from a natural reef environment near Palau. Two types of tissues were sampled from diseased corals: visibly affected and apparently healthy.

Results

Tag-based RNA-Seq followed by weighted gene co-expression network analysis identified groups of co-regulated differentially expressed genes between all health states (disease lesion, apparently healthy tissues of diseased colonies, and fully healthy). Differences between healthy and diseased tissues indicate activation of several innate immunity and tissue repair pathways accompanied by reduced calcification and the switch towards metabolic reliance on stored lipids. Unaffected parts of diseased colonies, although displaying a trend towards these changes, were not significantly different from fully healthy samples. Still, network analysis identified a group of genes, suggestive of altered immunity state, that were specifically up-regulated in unaffected parts of diseased colonies.

Conclusions

Similarity of fully healthy samples to apparently healthy parts of diseased colonies indicates that systemic effects of white syndromes on A. hyacinthus are weak, which implies that the coral colony is largely able to sustain its physiological performance despite disease. The genes specifically up-regulated in unaffected parts of diseased colonies, instead of being the consequence of disease, might be related to the originally higher susceptibility of these colonies to naturally occurring white syndromes.

Background

Increasing rates of disease have contributed greatly to global coral population declines over the last few decades [1,2]. The broadly defined “white syndromes” in Indo-Pacific regions, characterized in the field by tissue loss resulting in exposure of the coral skeleton, have been attributed to Vibrio spp. [3], a genus of bacteria involved in several coral diseases [4-8]. Other reports find no evidence of pathogenic bacteria in diseased corals [9] and instead blame stress-triggered programmed cell death for the manifestation of symptoms [10]. These conflicting conclusions, drawn mostly from culturing assays and histological observations, are further confounded by insufficient knowledge of the cnidarian immune response.

Corals, like all invertebrates, rely entirely on innate immunity for protection from invading pathogens. Features of innate immunity in corals include physical barriers [11], molecular pattern recognition [12], secretion of antimicrobial macromolecules [13], and cellular responses (e.g., phagocytosis) [14-16]. Recent efforts to characterize those features of immunity using various cnidarian genome and transcriptome sequence databases have identified putative components of coral stress management and immune response pathways by homology with better-studied organisms. The Acropora digitifera genome project revealed striking differences in innate immunity complexity in corals compared to a closely related cnidarian, Nematostella vectensis [17]. Whereas the N. vectensis genome encodes only a single Toll/Toll-like receptor (TLR), the A. digitifera genome included at least four TLRs, along with other related immune signaling molecules. Miller et al. reported the presence of TLR signaling components, including adaptor proteins that link that cascade with other signaling events, in the expressed sequence tag library of another acroporid coral, A. millepora [12]. Together these elements suggest an ability of corals to respond to pathogen-associated molecular patterns (PAMPs) via TLR recognition and integrate that signal to cellular responses such as inflammation and apoptosis. Toll/TLR signaling can activate NF-κB transcription factor that, upon nuclear localization, up-regulates transcription of immune response genes. In corals, the identities of those response genes and the roles they play remain unclear. Some suggested immune response genes include lectins, complement c3 and apextrin, proteins involved in non-self recognition, aggregation and cell lysis respectively [18]. Metabolism and calcification genes demonstrated differential expression in addition to immunity genes in A. millepora challenged with bacterial and viral immunogens, providing a more comprehensive picture of cellular events during an acute infection [19]. Global RNA-sequencing of A. cervicornis displaying signs of White Band Disease (WBD) revealed that disease significantly affected the expression of genes involved in immune processes and apoptosis [20]. The up-regulation of phagocytic cell surface receptors and reactive oxygen species (ROS) producing enzymes suggested that the phagocytosis and degradation of damaged cells drives the WBD response in corals.

These coral sequencing projects and experimental immune challenges have provided conclusive evidence that corals are capable of launching defensive responses upon direct contact with pathogens. A coral’s ability to communicate the recognition of that pathogen along the colony, however, is less understood. Coral polyps utilize a gastrovascular system lined with flagellated gastrodermal cells to transport organic products and zooxanthellae within the colony [21]. These channels are used to allocate energetic resources to areas that need them most, such as fast-growing branch tips [22-24] and wounded regions [25]. Radiolabeled carbon accumulation experiments have shown that corals preferentially direct energetic resources towards physically damaged regions [25] but away from disease-induced lesions [26]. These findings suggest that healthy coral tissues might possess means to detect and respond to an advancing disease lesion, but it is still unclear what the physiological consequences of this action might be.

Here we examine the gene expression profiles of A. hyacinthus displaying white syndromes (Figure 1) to determine the molecular consequences of the diseased condition. White syndromes advance along a colony in a way such that a distinct lesion forms between affected and unaffected tissues. Tissues ahead of the lesion are presumably healthy, while tissues at the lesion boundary are actively sloughing cells in response to infection. We compared gene expression profiles among three health states: affected tissues (diseased, “D”), apparently healthy tissues from diseased colonies (“ahead of the lesion”, “AL”), and tissues from completely unaffected colonies (healthy, “H”). Comparing healthy regions of diseased colonies to completely disease-free individuals provided an opportunity to look for expression patterns that might indicate a colony-wide systemic effect of infection and/or disease susceptibility. We used tag-based RNA-Seq [27] followed by weighted gene correlation network analysis [28] to achieve systems-level insight into molecular responses to chronic disease in corals.

Figure 1
figure 1

White syndrome in A.hyacinthus sampled in this study. (A) Close-up of the lesion area. (B) Position of sampled locations in diseased colonies: “D” – diseased, “AL” – ahead of the lesion. Tissues were also sampled from completely healthy individuals (“H”, not pictured). Photo credit: Carly Kenkel.

Full size image

Results

Differential gene expression between health states

Sequencing yielded an average of 6,367,219 reads per sample. An average of 19.5% of these remained after filtering and of these, an average 31.45% mapped to the transcriptome. A total of 44,701 isogroups (clusters of contigs representing the same gene, from here on referred to as “genes”) were detected. Reads were converted to unique transcript counts by removing PCR duplicates, yielding an average of 156,650 counts per sample (Additional file 1). A generalized linear model with contrasts between all three tissues detected differentially expressed genes between health states (Additional file 2). The disease-healthy contrast yielded 646 DEGs passing a Benjamini-Hochberg FDR cutoff of 10%. The disease-AL contrast yielded 333 DEGs passing an FDR cutoff of 10%. No genes passed the 10% FDR cutoff for the healthy-AL contrast. Between all contrasts, a total of 757 genes passed the FDR cutoff of 10%.

Principal coordinate analysis of the variance-stabilized data for all genes revealed expression differences mainly between disease (D) and the other two health states (Figure 2A). PCoA using the 3827 isogroups with an unadjusted p-value of less than 0.05 for any contrast revealed more differences between health states, but a significant overlap in expression of healthy and AL tissues remains (Figure 2B).

Figure 2
figure 2

Principal coordinate analysis clusters samples by health state. Samples cluster by the presence of disease symptoms (D vs. AL and H) when all genes are included in the PCoA (A). Differences between health states become more evident when PCoA is performed on DEGs (unadjusted p-value < 0.05) only (B).

Full size image

Gene Ontology (GO) Enrichment

Functional enrichments between all three contrasts allow a general examination of the molecular functions and biological processes being differentially regulated between health states. The enriched groups of both the disease-healthy and disease-AL contrasts were largely identical (Figure 3 and Additional file 3). Ribosomal proteins, oxidative stress responses, and translation factor activity were up-regulated in diseased tissues compared to both AL and healthy tissues. Likewise, receptor activity, regulation of biological quality, and extracellular matrix components (collagens) were down-regulated in diseased tissues compared to both healthier states. No GO terms were significantly enriched (FDR 10%) for the healthy-AL contrast.

Figure 3
figure 3

Gene ontology categories enriched by genes up-regulated (red) or down-regulated (blue) in diseased compared to fully healthy samples, summarized by molecular function (MF), biological process (BP), and cellular component (CC). The size of the font indicates the significance of the term as indicated by the inset key. The fraction preceding the GO term indicates the number of genes annotated with the term that pass an unadjusted p-value threshold of 0.05. The trees indicate sharing of genes among GO categories (the categories with no branch length between them are subsets of each other).

Full size image

Gene expression analysis by contrast

Gene expression heatmaps were constructed to show the relative expression patterns of the top most significant DEGs for each contrast (Figure 4A and Additional files 4, 5). Complete lists of all annotated DEGs with log fold changes for each contrast can be found in Additional files 6, 7 and 8.

Figure 4
figure 4

Expression of DEGs significant for disease-healthy contrast among health states. (A) Heatmap for top DEGs (FDR = 0.01). Rows are genes, columns are samples ordered as in the bottom panel: ahead-of-lesion (AL), healthy (H), and diseased (D). The color scale is in log2 (fold change relative to the gene’s mean). The tree is a hierarchical clustering of genes based on Pearson’s correlation of their expression across samples. (B) Principal coordinate analysis of all DEGs at 10% FDR for disease-healthy contrast.

Full size image

Diseased vs. healthy

Genes found to be up-regulated (Benjamini-Hochberg FDR < 0.01) in diseased tissues compared to healthy corals include key members of the oxidative stress response in corals (e.g., catalases and peroxidases) and pentose phosphate metabolism (transketolase, transaldolase, and 6-phosphogluconate dehydrogenase). Both proteinases (astacin and cathepsin L) and protease inhibitors (alpha-macroglobulin and serine proteinase inhibitor Ku-type) were up-regulated in diseased tissues. Two of the genes annotated as C-type lectin, a carbohydrate-binding protein, and malate synthase, a key enzyme of the glyoxylate cycle, were also up-regulated in symptomatic tissues. Down-regulated genes (Benjamini-Hochberg FDR < 0.01) include those encoding extracellular matrix constituents (collagens, heparin sulfate proteoglycans) and carbonic anhydrase, a key enzyme in coral skeletal deposition. Red fluorescent protein was also down-regulated in diseased tissues, a hallmark of the coral stress response [29-33].

AL vs. diseased and healthy

The expression differences between diseased and AL tissues within a colony paralleled the expression differences between diseased tissues and healthy corals. At the same significance threshold (Benjamini-Hochberg FDR < 0.01), almost the exact same top-candidate genes were identified (Additional files 4, 5). No differentially expressed genes passed the Benjamini-Hochberg 10% FDR when comparing the AL tissues and healthy corals.

To better quantify the behavior of disease-responsive genes in AL samples, the expression of genes passing a 10% FDR cutoff for the disease-healthy comparison was studied using principal coordinate analysis. Tukey’s test based on the first principal coordinate values revealed highly significant differences between H and D as well as between AL and D (P < 0.001 in both cases). Although AL samples appeared to be intermediate between H and D samples (Figure 4B), their scores along the first principal coordinate axis were not significantly different from healthy samples (P = 0.11).

Correlation between gene network modules and health states

A total of 6737 DEGs with unadjusted p-value < 0.1 were input into WGCNA for network analysis. A sample network was constructed to identify outlying samples with a standardized connectivity score of less than −2.5. One sample (diseased individual “4”) was identified as an outlier and removed from subsequent analysis (Additional file 9). Twelve unique modules, assigned arbitrarily color labels, remained after merging highly correlated modules. Of these twelve modules, eight were highly correlated to a single coral individual and one (grey) is reserved to contain genes that do not fall into any co-expression module. The remaining three modules were highly correlated with the health states (Additional file 10). Since we assembled these modules using a signed network, the sign of the correlation is equivalent to the direction of expression change with respect to the trait. For example, a module that is significantly negatively correlated to diseased corals contains genes that are down-regulated in that state.

The eigengene of the dark green module (1155 genes) was strongly correlated with diseased-healthy contrast (Pearson’s R 2 = 0.83, Pcor = 1e-6, Additional file 10). The genes within this module are up-regulated in diseased tissues and down-regulated in healthy tissues, while tending to be down-regulated in AL samples. Conversely, the turquoise module (669 genes) was up-regulated in healthy tissues and down-regulated in diseased tissues (Pearson’s R 2 = −0.83, Pcor = 9e-7, Additional file 10). Expression of these two modules in AL samples demonstrated similar direction of change as in healthy tissues, although the change was not statistically significant. Notably, one module was identified (green, 661 genes) that was significantly up-regulated in AL (Pearson’s R 2 = 0.64, Pcor = 0.001) and (to a lesser extent) down-regulated in diseased tissues (Pearson’s R 2 = −0.44, Pcor = 0.04), while remaining unchanged in healthy tissues (Additional file 10).

Within module gene expression analysis

Hierarchically clustered gene expression heatmaps were constructed to show the relative expression patterns of the genes within each module that best represent the module and show significant correlations to the health state, based on module membership and gene significance values (Figure 5 and Additional file 11).

Figure 5
figure 5

Gene expression heatmaps of annotated DEGs with a module membership and gene significance score greater than 0.6. Rows are genes, columns are samples ordered as in the bottom panel: ahead-of-lesion (AL), healthy (H), and diseased (D). The color scale is in log2 (fold change relative to the gene’s mean). The trees are hierarchical clustering of genes based on Pearson’s correlation of their expression across samples. The color block of the trees indicates the module to which these genes belong (dark green, turquoise, and green).

Full size image

Discussion

The major pattern of variation in gene expression was between asymptomatic (healthy and AL) and diseased (D) corals (Figure 2). There were no statistically significant differences in gene expression between healthy colonies (H) and asymptomatic parts of diseased colonies (AL), suggesting that these states are physiologically similar. The genes that are differentially regulated between diseased and healthy corals show a subtle trend towards disease-like gene expression in asymptomatic tissues of diseased colonies (Figure 4). This trend is, however, not statistically significant, indicating that white syndromes have little effect on the physiology of the unaffected portion of A. hyacinthus colony. Gene co-expression network analysis revealed groups of genes co-regulated with respect to each of the three states, including a group of genes specifically up-regulated in AL samples.

Diseased tissues up regulate immune response elements

Innate immunity provides immediate protection against non-self and responds to physical injury. Three general steps are involved in an innate immune response: detection, defense activation, and effector responses to neutralize the threat. Tissues sampled from the lesion of disease progression in corals exhibiting white syndromes have enhanced expression of genes involved in each of these three immune response phases. C-type lectins act as pattern recognition receptors to activate pathogen elimination through phagocytosis in invertebrates [34]. Cnidarian genomes encode c-type lectin genes with highly variable substrate regions, leading to hypotheses that these proteins recognize a large variety of pathogens [35]. In A. millepora, mannose-binding C-type lectins have been shown to respond immediately following an immune challenge (only 45 minutes after lipopolysaccharide injection in [36]), but show no significant response at later time points [37]. The up-regulation of C-type lectins in tissues at the lesion may suggest that these tissues have very recently become infected. The second phase of an immune response prepares targets for elimination via antimicrobial peptide synthesis and immune cell activation. While this experiment did not discover any differentially regulated antimicrobial peptides, we do detect the activation of immune activating proteins C4, alpha-macroglobulin and CD109 in diseased tissues. The lectin pathway of immune activation is triggered by lectins binding a pathogen-associated molecule and results in the activation complement component factor C4 and C3 [38]. These complement factors, along with alpha-macroglobulin and CD109, tag pathogens and secreted proteases for elimination. In the final phase of an innate immune response, foreign organisms are engulfed and destroyed by phagocytic immune cells. Lysosomes within these phagocytic cells contain proteins capable of degrading engulfed material via the production of reactive oxygen species (ROS) or proteolytic enzymes, such as cathepsins. The up-regulation of cathepsin L in diseased tissue may be a consequence of such phagocytic activity. The up-regulation of immune-related transcripts in diseased corals is consistent with previous studies of both naturally occurring disease and experimental pathogen challenges [19,20,39]. Just like the rest of genes exhibiting H-D difference, these responses are confined to the symptomatic regions of the coral. One possible explanation of this fact is that the immunity-related gene expression changes are elicited by direct contact with a pathogen rather than a systemic signal throughout the colony.

Switch to lipid-based metabolism in diseased tissues

Transcripts involving lipid and carbohydrate metabolism (triacylglycerol lipase, phosphoenolpyruvate carboxykinase), including glyoxylate cycle metabolism (malate synthase), were up-regulated in diseased tissues compared to asymptomatic tissues. The differential regulation of these metabolic genes suggests that diseased corals may be utilizing stored energy reserves more than healthy corals. Fatty acids derived from stored lipids are oxidized by the beta-oxidation pathway and release acetyl-CoA to enter the citric acid cycle. Both carbons in one molecule of acetyl-CoA are consumed during the decarboxylation steps of the citric acid cycle and energy is released. The glyoxylate cycle is an alternative route through the citric acid cycle that allows organisms to thrive on two-carbon sources by catalyzing the conversion of acetyl-CoA to malate and succinate via a glyoxylate intermediate, bypassing the decarboxylation steps in the citric acid cycle [40]. These four carbon compounds contribute to the energetic requirements of the cell and serve as building blocks for cellular components, fulfilling all the same necessities of the citric acid cycle without the need to replenish oxaloacetate from the diet. One of the most significantly up-regulated transcripts in diseased coral tissues was malate synthase, one of the key enzymes of the glyoxylate cycle along with isocitrate lyase. Glyoxylate cycle enzymes are rare throughout the animal kingdom, but bioinformatic analyses suggest they exist in cnidarians [41]. As additional support for a potential role of glyoxylate cycle metabolism in coral stress responses, glyoxylate cycle transcripts were up-regulated in A. palmata larvae subjected to thermal stress [42]. In higher plants, glyoxylate enzymes are active when the cell is switching from photosynthetic production of sugars to scavenging pathways from stored and structural lipids, as in starvation and/or senescence [43]. In corals, this metabolic shift might indicate a decline in shared energy reserves with zooxanthellae, presumably due to stress-induced symbiont loss.

Oxidative stress response genes are up-regulated in diseased tissues

Reactive oxygen species (ROS) are produced as a consequence of fatty acid oxidation. The up-regulation of antioxidants that protect the cell from these harmful byproducts corals could be a consequence of increased fatty acid metabolism. This explanation coincides well with the observed up-regulation of lipid metabolism and antioxidant (catalase, peroxidase) transcripts in diseased corals. The production of ROS is also a fundamental element of the innate immune response. While ROS are capable of neutralizing phagocytized pathogens, the harm they cause to the host must be countered if an organism is to withstand its own immune response. Catalases and peroxidases capable of hydrolyzing harmful peroxides provide a mechanism of such self-protection. The up-regulation of oxidative stress response genes is well characterized in corals experiencing thermal stress [44], physical stress [45], and infectious disease [20].

Matrix metalloproteinases are up-regulated in diseased tissues

Stony corals are subject to many potential sources of physical injury such as predators [46], boring organisms [47] and storms [48,49]. The tissue regeneration mechanisms employed by corals that have sustained a physical injury are common to wound-healing processes across metazoans [50]. One of these steps involves restructuring of the extracellular matrix to encourage tissue regeneration. Matrix metalloproteinases (MMPs) are a group of enzymes capable of such activities and have been shown to play a direct role in wound repair in Hydra [51]. In addition, MMPs act on pro-inflammatory cytokines to direct inflammation due to wounding and innate immune responses to pathogens (reviewed in [52]). The up-regulation of MMPs in response parasitic protists in a gorgonian coral suggests that these proteins are active in the immune response of cnidarians [53]. Astacin and gelatinase have matrix metalloproteinase activities and were up-regulated specifically in affected coral tissues. Additionally, a protease inhibitor alpha-macroglobulin was up-regulated, which is a vital component of the innate immune response that inactivates bacterial secreted proteases, thus compromising their virulence [54].

Calcification genes are down-regulated in diseased tissues

Calcification rates in reef-building corals are sensitive to several environmental variables such as light, pH, and temperature [1,55-57]. While this experiment did not directly measure coral calcification, the identification of DEGs with functions in biomineralization suggests that disease negatively impacts coral skeletal deposition. Both general extracellular matrix structural components and coral-specific calcification functions were differentially regulated in diseased tissues compared to asymptomatic tissues. Coral biomineralization is directed by an extracellular skeletal organic matrix comprised of secreted glycosylated proteins [58]. These proteins include collagens and negatively charged macromolecules (like chondroitin sulfate proteoglycans) that bind calcium ions to aid in crystal formation [59]. Several collagens and a protein with high similarity to nematogalectin, a collagen family protein that forms a major structural component in Hydra nematocyst tubules [60], were down-regulated in diseased tissues. The down-regulation of these genes in diseased tissues suggests a weakening of the coral skeletal organic matrix and thus a diminished capacity for biomineral deposition. The potential impact of disease on coral skeletal growth is most clearly revealed by the down-regulation of carbonic anhydrase, an enzyme that plays a fundamental role in mediating bicarbonate supplies for calcification in scleractinian corals [61-63].

AL-specific gene expression: a systemic response to disease or factors contributing to disease susceptibility?

Genes that are specifically regulated in apparently healthy tissues of diseased colonies could represent systemic (i.e., colony-wide) response to disease. However, there is another possible interpretation: since expression of these genes is correlated with natural appearance of disease, it might signify disease susceptibility rather than disease effects. A recent study employing similar a sampling scheme to investigate transcriptomic effects of yellow band disease (YBD) in Orcibella faveolata [39] found that expression in asymptomatic regions of diseased colonies was intermediate between completely healthy corals and diseased tissue, which fits well with the systemic response interpretation. In our study the AL expression was most similar to the healthy state (Figure 2 and 4A), although genes differentially regulated between diseased and healthy states demonstrated a non-significant trend towards intermediate expression in AL samples (Figure 4B). The difference between the two studies could potentially be explained by unequal levels of colony integration between Orbicella and Acropora ([64] and references therein), which could affect the extent of the systemic signaling and/or spread of the pathogen throughout the colony. The co-expression network analysis revealed a sizeable (661 genes) module that was up-regulated in the AL state compared to D and H states (Additional file 10). Among the genes most strongly associated with this module were the genes coding for the immunity-related Tolloid-like protein and the hypoxia inducible factor prolyl 4-hydroylase (HIF-P4H, Figure 5). Up-regulation of HIF-P4H suggests that healthy tissues of diseased colonies might be experiencing hypoxic conditions [65,66]. Notably, HIF-P4H has also been shown to modulate immune responses by modifying the kinase responsible for releasing NF-κB from its inhibitor [67]. Up-regulation of these genes in healthy parts of diseased colonies might therefore be a sign of altered immunity state potentially explaining higher disease susceptibility of the affected colonies in nature.

Conclusions

Our gene expression analysis identified several immune, repair, and metabolic molecular pathways expressed in coral regions affected with white syndromes. In contrast to Orbicella faveolata, A. hyacinthus does not show pronounced propagation of these responses to regions of the colony not visibly affected by disease, suggesting that the effect of chronic white syndromes on colony-wide A. hyacinthus physiology is small. Instead, asymptomatic regions of diseased colonies show gene expression signatures potentially related to higher disease susceptibility of the affected coral individuals. Further studies of natural disease-associated gene expression will contribute towards the development of diagnostic tools to predict and manage coral disease outbreaks.

Methods

Sampling

Coral fragments from 16 colonies of A. hyacinthus were sampled in the spring of 2011 along the eastern coast of Palau (7° 18.738′ N, 134° 30.423′ E) and immediately stored in RNAlater (Ambion). Eight of these colonies were visibly affected with white syndromes (Figure 1). Colonies exhibiting diffuse tissue loss along a lesion of apparently healthy tissue directly adjacent to exposed white skeleton in accordance with [68] were sampled. Colonies displayed no obvious signs of predation. The remaining eight colonies were completely symptom-free (designated healthy, H). From the eight affected colonies, coral fragments were sampled from both the lesion interface between diseased and healthy tissues (diseased, D) and areas well ahead of the lesion (AL, Figure 1B). AL coral fragments were sampled from approximately midway between the lesion boundary and the edge of the colony. AL tissues are unlikely to be in direct contact with pathogen since previous studies have demonstrated declines in pathogens only ~1 cm in advance of a white syndromes lesion [69].

Transcriptome Assembly and Annotation

The A. hyacinthus transcriptome has been generated from 5-day aposymbiotic larvae as described previously [70]. It was annotated based on two resources: the proteome of the starlet sea anemone Nematostella vectensis [71] and in-depth annotations of the Acropora digitifera proteome [72]. Based on manual verifications of a subset of A. digitifera annotations, they were pre-filtered to include only protein sequences longer than 60 amino acids with the annotation assigned based on the listed e-value = 1e-20 or better. The GO, KEGG (“Kyoto Encyclopedia of Genes and Genomes”), KOG (“euKaryotic orthologous groups”), and gene name annotations were transferred to an A. hyacinthus contig if the contig matched one or both of these two resources with e-value = 1e-4 or better in blastx [73]. The GO and KOG annotations assigned to genes that were denoted FOG (“fuzzy orthologous group”, [74]) in the N. vectensis data were removed, since such genes encode proteins with common domains and cannot be functionally annotated based on homology alone. The annotated A. hyacinthus transcriptome has been released for unrestricted use prior to this publication, http://www.bio.utexas.edu/research/matz_lab/matzlab/Data.html.

Tag-based RNA-Seq

Libraries were prepared following [27] and sequenced using Applied Biosystems SOLiD v.3 platform. Read trimming, quality filtering, mapping, and conversion to per-gene counts was performed as described previously [27] with one modification: reads mapping to the same starting coordinate in the reference and aligning with 100% identity along the full length of the shorter read were discarded as potential PCR duplicates. The current step-by-step library preparation protocol as well as bioinformatic pipeline are available at https://github.com/z0on/tag-based_RNAseq; note however that the current version of the tag-based RNA-seq method utilizes advanced procedure for PCR duplicates removal based on degenerate tags incorporated during cDNA synthesis, and assumes sequencing on the Illumina HiSeq instrument. R scripts and input files can be accessed on Dryad: doi:NN.NNNN/dryad.NNNN. Raw sequence data can be found on NCBI-SRA, Accession: PRJNANNNNNN.

Identification of Differentially Expressed Genes (DEGs)

All statistical analyses were performed using R3.1.1 [75]. DEGs were identified using a generalized linear model implemented by the R package DESeq2 [76]. No outlying samples were detected by the arrayQualityMetrics package [77]. DESeq2 performed automatic independent filtering to remove lowly abundant transcripts and maximize the rate of DEG discovery post multiple testing correction at an alpha of 0.1. P-values for significance of contrasts between all three health states were generated based on Wald statistics and were adjusted for multiple testing using the false discovery rate method [78]. The contrasts resulted in tables including adjusted and unadjusted p-values and log2 fold changes that were used in downstream analyses.

Gene coexpression network analysis

A weighted gene correlation network analysis (WGCNA, [28]) was used to identify groups of co-regulated genes in an unsupervised way. Genes with an unadjusted p-value < 0.1 for any of the three contrasts as determined by the generalized linear model testing for the effect of health state were input into WGCNA. A sample network was constructed to identify outlying samples with a standardized connectivity score of less than −2.5 [79]. A signed gene co-expression network was constructed with a soft threshold power of 24. Groups of co-regulated genes (modules) correlated with each other with the Pearson correlation coefficient 0.42 or better were merged. The eigengenes of the resulting modules (the first principal component of the expression matrix corresponding to the genes included in the module) were correlated with health states (H, AL, or D).

Assessing the robustness of the analysis

Low quality of the SOLiD sequencing resulted in low number of reads mapped, raising concerns about the reliability of the data. In tag-based RNA-seq, unlike standard RNA-seq, every count represents an observation of an independent transcript. Thus low counts could still provide sufficient quantitative information about transcript abundances. In addition, high level of biological replication (n = 8 per group) in our experiment should have compensated for the low counts within each replicate to a certain degree. To confirm that low counts did not result in inflated false discovery rate, we have simulated a series of count datasets based on the empirical per-gene total counts and coefficients of variation across samples as well as empirical sample size factors, which included no effect of experimental conditions. Analysis of these simulated datasets recovered nearly identical sample size factors and highly similar dispersion estimates as in real data (Additional file 12). When these simulated datasets were analyzed with DESeq2 using the same models as real data, at most four genes passed the 10% Benjamini-Hochberg false discovery rate (FDR) cutoff for each contrast. Compared to 646 genes passing the FDR 10% cutoff for the disease-healthy comparison and 333 genes passing the same cutoff for the disease-AL comparison in the real dataset this is much less than 10%, indicating that the real data analysis was conservative. The simulation-based p-value cutoff achieving the empirical 10% FDR (Additional file 13A-C) would have yielded 1.05-1.5X more DEGs for these comparisons than the Benjamini-Hochberg procedure. Notably, the Benjamini-Hochberg correction did not yield any DEGs for the healthy-AL comparison, and accordingly, the DEG discovery rate in this comparison was even slightly lower than the simulation-based false discovery rate (Additional file 13C), indicating that for this comparison DESeq2 analysis did not provide sufficient power. To keep the analysis conservative, we chose to report DESeq2-based DEGs discovered using the Benjamini-Hochberg procedure. The procedure for inferring empirical FDR based on simulations described above has been implemented in the R package empiricalFDR.DESeq2, hosted within the Comprehensive R Archive Network (CRAN).

WGCNA analysis provides one additional confirmation that the observed gene expression differences are driven by biological factors rather than stochasticity. WGCNA constructs gene co-expression modules from their correlation pattern across samples without using the information about how the samples are distributed among experimental conditions. The fact that post hoc the module eigengenes correlate strongly with coral condition (Additional file 10) indicates that the gene expression patterns in the data truly reflect the biological processes related to disease. To verify this, we shuffled the condition designations among samples and indeed observed that the correlations with co-expression modules disappeared (Additional file 14).

Principal coordinate analysis

Principal coordinate analysis to visualize clustering of gene expression between health states was performed using the “adegenet” package [80] using variance stabilized data for all genes and subsequently with only candidate differentially expressed genes (unadjusted p-value < 0.05), based on Manhattan distances which correspond to the sum of absolute log-fold changes across all genes. Effects of the three health states (“D”, “H”, and “AL”) were calculated using the multivariate analysis of variance function “adonis” of the R package “vegan” [81]. Tukey’s tests between specific health states based on the values of the first principal coordinate were performed using function TukeyHSD in R.

Functional Summaries

Gene ontology enrichment analysis was performed using adaptive clustering of GO categories and Mann–Whitney U tests [82] based on ranking of signed log p-values (GO-MWU, https://github.com/z0on/GO_MWU). Gene expression heatmaps with hierarchical clustering of expression profiles were created with the pheatmap package in R [83].

Abbreviations

ALk:

Ahead-of-the-lesion

DEG:

Differentially expressed genes

FDR:

False discovery rate

FOG:

Fuzzy orthologous groups

GO:

Gene ontology

HIF-P4H:

Hypoxia-inducible factor prolyl 4-hydroxylase

KEGG:

Kyoto encyclopedia of genes and genomes

KOG:

EuKaryotic orthologous groups

MMP:

Matrix metalloproteinase

NF-κB:

Nuclear factor kappa-light-chain-enhancer of activated B cells

PAMP:

Pathogen-associated molecular pattern

ROS:

Reactive oxygen species

TLR:

Toll-like receptor

WBD:

White band disease

WGCNA:

Weighted gene correlation network analysis

YBD:

Yellow band disease

References

  1. Hoegh-Gulberg O: Climate change, coral bleaching and the future of the world’s coral reefs. Mar Freshw Resour 1999;50:839–866.

  2. Porter JW, Dustan P, Jaap WC, Patterson KL, Kosmynin V, Meier OW, Patterson ME, Parsons M: Patterns of spread of coral disease in the Florida Keys. Hydrobiologia 2001;460:1–24.

  3. Sussman M, Willis BL, Victor S, Bourne DG. Coral pathogens identified for White Syndrome (WS) epizootics in the Indo-Pacific. PLoS One. 2008;3:e2393.

    Article  PubMed Central  PubMed  Google Scholar 

  4. Kushmaro A, Banin E, Loya Y, Stackebrandt E, Rosenberg E. Vibrio shiloi sp. nov., the causative agent of bleaching of the coral Oculina patagonica. Int J Syst Evol Microbiol. 2001;51:1383–8.

    CAS  PubMed  Google Scholar 

  5. Ben-Haim Y. Vibrio coralliilyticus sp. nov., a temperature-dependent pathogen of the coral Pocillopora damicornis. Int J Syst Evol Microbiol. 2003;53:309–15.

    Article  CAS  PubMed  Google Scholar 

  6. Cervino JM, Thompson FL, Gomez-Gil B, Lorence E, Goreau TJ, Hayes RL, et al. The Vibrio core group induces yellow band disease in Caribbean and Indo-Pacific reef-building corals. J Appl Microbiol. 2008;105:1658–71.

    Article  CAS  PubMed  Google Scholar 

  7. Vezzulli L, Previati M, Pruzzo C, Marchese A, Bourne DG, Cerrano C. Vibrio infections triggering mass mortality events in a warming Mediterranean Sea. Environ Microbiol. 2010;12:2007–19.

    Article  CAS  PubMed  Google Scholar 

  8. Ushijima B, Videau P, Burger A, Shore-Maggio A, Runyon CM, Sudek M, et al. Vibrio coralliilyticus strain OCN008 is an etiological agent of acute Montipora white syndrome. Appl Environ Microbiol. 2014;80:2102–9.

    Article  PubMed Central  PubMed  Google Scholar 

  9. Roff G, Kvennefors ECE, Fine M, Ortiz J, Davy JE, Hoegh-Guldberg O. The ecology of “Acroporid white syndrome”, a coral disease from the southern Great Barrier Reef. PLoS One. 2011;6:e26829.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  10. Ainsworth TD, Kvennefors EC, Blackall LL, Fine M, Hoegh-Guldberg O. Disease and cell death in white syndrome of Acroporid corals on the Great Barrier Reef. Mar Biol. 2007;151:19–29.

    Article  Google Scholar 

  11. Brown B, Bythell J. Perspectives on mucus secretion in reef corals. Mar Ecol Prog Ser. 2005;296:291–309.

    Article  CAS  Google Scholar 

  12. Miller DJ, Hemmrich G, Ball EE, Hayward DC, Khalturin K, Funayama N, et al. The innate immune repertoire in cnidaria–ancestral complexity and stochastic gene loss. Genome Biol. 2007;8:R59.

    Article  PubMed Central  PubMed  Google Scholar 

  13. Vidal-Dupiol J, Ladrière O, Destoumieux-Garzón D, Sautière P-E, Meistertzheim A-L, Tambutté E, et al. Innate immune responses of a scleractinian coral to vibriosis. J Biol Chem. 2011;286:22688–98.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  14. Mydlarz LD, Harvell CD. Peroxidase activity and inducibility in the sea fan coral exposed to a fungal pathogen. Comp Biochem Physiol Part A. 2007;146:54–62.

    Article  Google Scholar 

  15. Mydlarz LD, Holthouse SF, Peters EC, Harvell CD. Cellular responses in sea fan corals: granular amoebocytes react to pathogen and climate stressors. PLoS One. 2008;3:e1811.

    Article  PubMed Central  PubMed  Google Scholar 

  16. Olano CT, Bigger CH. Phagocytic activities of the gorgonian coral Swiftia exserta. J Invertebr Pathol. 2000;76:176–84.

    Article  CAS  PubMed  Google Scholar 

  17. Shinzato C, Shoguchi E, Kawashima T, Hamada M, Hisata K, Tanaka M, et al. Using the Acropora digitifera genome to understand coral responses to environmental change. Nature. 2011;476:320–3.

    Article  CAS  PubMed  Google Scholar 

  18. Puill-Stephan E, Seneca FO, Miller DJ, Van Oppen MJH, Willis BL. Expression of putative immune response genes during early ontogeny in the coral Acropora millepora. PLoS One. 2012;7:e39099.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  19. Weiss Y, Forêt S, Hayward DC, Ainsworth T, King R, Ball EE, et al. The acute transcriptional response of the coral Acropora millepora to immune challenge: expression of GiMAP/IAN genes links the innate immune responses of corals with those of mammals and plants. BMC Genomics. 2013;14:400.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  20. Libro S, Kaluziak ST, Vollmer SV. RNA-seq profiles of immune related genes in the staghorn coral Acropora cervicornis infected with white band disease. PLoS One. 2013;8:e81821.

    Article  PubMed Central  PubMed  Google Scholar 

  21. Gladfelter EH. Circulation of Fluids in the Gastrovascular System of the Reef Coral Acropora cervicornis. Biol Bull. 1983;165:619–36.

    Article  Google Scholar 

  22. Pearse VB, Muscatine L. Role of Symbiotic Algae (Zooxanthellae) in Coral Calcification. Biol Bull. 1971;141:350–63.

    Article  CAS  Google Scholar 

  23. Taylor DL. Intra-colonial transport of organic compounds and calcium in some Atlantic reef corals. In: Proceedings, 3rd Int Coral Reef Symp. 1977. p. 431–6.

    Google Scholar 

  24. Fang LS, Chen YWJ, Chen CS. Why does the white tip of stony coral grow so fast without zooxanthellae? Mar Biol. 1989;103:359–63.

    Article  Google Scholar 

  25. Oren U, Rinkevich B, Loya Y. Oriented intra-colonial transport of 14C labeled materials during coral regeneration. Mar Ecol Prog Ser. 1997;161:117–22.

    Article  Google Scholar 

  26. Roff G, Hoegh-Guldberg O, Fine M. Intra-colonial response to Acroporid “white syndrome” lesions in tabular Acropora spp. (Scleractinia). Coral Reefs. 2006;25:255–64.

    Article  Google Scholar 

  27. Meyer E, Aglyamova GV, Matz MV. Profiling gene expression responses of coral larvae (Acropora millepora) to elevated temperature and settlement inducers using a novel RNA-Seq procedure. Mol Ecol. 2011;20:3599–616.

    CAS  PubMed  Google Scholar 

  28. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008;9:559.

    Article  PubMed Central  PubMed  Google Scholar 

  29. Roth MS, Deheyn DD. Effects of cold stress and heat stress on coral fluorescence in reef-building corals. Sci Rep. 2013;3:1421.

    Article  PubMed Central  PubMed  Google Scholar 

  30. Smith-Keune C, Dove S. Gene expression of a green fluorescent protein homolog as a host-specific biomarker of heat stress within a reef-building coral. Mar Biotechnol. 2007;10:166–80.

    Article  PubMed  Google Scholar 

  31. Rodriguez-Lanetty M, Harii S, Hoegh-Guldberg O. Early molecular responses of coral larvae to hyperthermal stress. Mol Ecol. 2009;18:5101–14.

    Article  CAS  PubMed  Google Scholar 

  32. DeSalvo MK, Voolstra CR, Sunagawa S, Schwarz J, Stillman JH, Coffroth M, et al. Differential gene expression during thermal stress and bleaching in the Caribbean coral Montastraea faveolata. Mol Ecol. 2008;17:3952–71.

    Article  CAS  PubMed  Google Scholar 

  33. Bay LK, Ulstrup KE, Nielsen HB, Jarmer H, Goffard N, Willis BL, et al. Microarray analysis reveals transcriptional plasticity in the reef building coral Acropora millepora. Mol Ecol. 2009;18:3062–75.

    Article  CAS  PubMed  Google Scholar 

  34. Fujita T. Evolution of the lectin-complement pathway and its role in innate immunity. Nat Rev Immunol. 2002;2:346–53.

    Article  CAS  PubMed  Google Scholar 

  35. Kvennefors ECE, Leggat W, Hoegh-Guldberg O, Degnan BM, Barnes AC. An ancient and variable mannose-binding lectin from the coral Acropora millepora binds both pathogens and symbionts. Dev Comp Immunol. 2008;32:1582–92.

    Article  CAS  PubMed  Google Scholar 

  36. Kvennefors ECE, Leggat W, Kerr CC, Ainsworth TD, Hoegh-Guldberg O, Barnes AC. Analysis of evolutionarily conserved innate immune components in coral links immunity and symbiosis. Dev Comp Immunol. 2010;34:1219–29.

    Article  CAS  PubMed  Google Scholar 

  37. Brown T, Bourne D, Rodriguez-Lanetty M. Transcriptional Activation of c3 and hsp70 as Part of the Immune Response of Acropora millepora to Bacterial Challenges. PLoS One. 2013;8:e67246.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  38. Matsushita M, Fujita T. Activation of the Classical Complement Pathway by Marmose-binding Protein in Association with a Novel Cls-like Serine Protease. J Exp Med. 1992;176:1497–502.

    Article  CAS  PubMed  Google Scholar 

  39. Closek CJ, Sunagawa S, DeSalvo MK, Piceno YM, DeSantis TZ, Brodie EL, Weber MX, Voolstra CR, Andersen GL, Medina M: Coral transcriptome and bacterial community profiles reveal distinct Yellow Band Disease states in Orbicella faveolata. ISME J 2014;8:2411–2422.

  40. Kornberg HL, Krebs HA. Synthesis of Cell Constituents from C2-Units by a Modified Tricarboxylic Acid Cycle. Nature. 1957;179:988–91.

    Article  CAS  PubMed  Google Scholar 

  41. Kondrashov FA, Koonin EV, Morgunov IG, Finogenova TV, Kondrashova MN. Evolution of glyoxylate cycle enzymes in Metazoa: evidence of multiple horizontal transfer events and pseudogene formation. Biol Direct. 2006;1:31.

    Article  PubMed Central  PubMed  Google Scholar 

  42. Polato NR, Altman NS, Baums IB. Variation in the transcriptional response of threatened coral larvae to elevated temperatures. Mol Ecol. 2013;22:1366–82.

    Article  CAS  PubMed  Google Scholar 

  43. Wanner L, Keller F, Matile P. Metabolism of radiolabelled galactolipids in senescent barley leaves. Plant Sci. 1991;78:199–206.

    Article  CAS  Google Scholar 

  44. DeSalvo M, Sunagawa S, Voolstra C, Medina M. Transcriptomic responses to heat stress and bleaching in the elkhorn coral Acropora palmata. Mar Ecol Prog Ser. 2010;402:97–113.

    Article  CAS  Google Scholar 

  45. Lõhelaid H, Teder T, Tõldsepp K, Ekins M, Samel N. Up-Regulated Expression of AOS-LOXa and Increased Eicosanoid Synthesis in Response to Coral Wounding. PLoS One. 2014;9:e89215.

    Article  PubMed Central  PubMed  Google Scholar 

  46. Rotjan R, Lewis S. Impact of coral predators on tropical reefs. Mar Ecol Prog Ser. 2008;367:73–91.

    Article  Google Scholar 

  47. Fang L-S, Shen P. A living mechanical file: the burrowing mechanism of the coral-boring bivalve Lithophaga nigra. Mar Biol. 1988;97:349–54.

    Article  Google Scholar 

  48. Bythell JC, Bythell M, Gladfelter EH. Initial results of a long-term coral reef monitoring program : Impact of Hurricane Hugo at Buck Island Reef National Monument, St. Croix, U.S. Virgin Islands. J Exp Mar Bio Ecol. 1993;172:171–83.

    Article  Google Scholar 

  49. Bythell JC, Hillis-Starr ZM, Rogers CS. Local variability but landscape stability in coral reef communities following repeated hurricane impacts. Mar Ecol Prog Ser. 2000;204:93–100.

    Article  Google Scholar 

  50. Palmer CV, Traylor-Knowles NG, Willis BL, Bythell JC. Corals use similar immune cells and wound-healing processes as those of higher organisms. PLoS One. 2011;6:e23992.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  51. Bosch TCG. Why polyps regenerate and we don’t: towards a cellular and molecular framework for Hydra regeneration. Dev Biol. 2007;303:421–33.

    Article  CAS  PubMed  Google Scholar 

  52. Parks WC, Wilson CL, López-Boado YS. Matrix metalloproteinases as modulators of inflammation and innate immunity. Nat Rev Immunol. 2004;4:617–29.

    Article  CAS  PubMed  Google Scholar 

  53. Burge C, Mouchka ME, Harvell CD, Roberts S. Immune response of the Caribbean sea fan, Gorgonia ventalina, exposed to an Aplanochytrium parasite as revealed by transcriptome sequencing. Front Physiol. 2013;4:1–9.

    Article  Google Scholar 

  54. Armstrong PB, Quigley JP. a 2 -macroglobulin : an evolutionarily conserved arm of the innate immune system. Dev Comp Immunol. 1999;23:375–90.

    Article  CAS  PubMed  Google Scholar 

  55. Tentori E, Allemand D. Light-enhanced calcification and dark decalcification in isolates of the soft coral Cladiella sp. during tissue recovery. Biol Bull. 2006;211:193–202.

    Article  CAS  PubMed  Google Scholar 

  56. Clausen CD, Roth AA. Effect of Temperature and Temperature Adaptation on Calcification Rate in the Hermatypic Coral Pocillopora damicomis. Mar Biol. 1975;33:93–100.

    Article  Google Scholar 

  57. Kleypas JA. Geochemical Consequences of Increased Atmospheric Carbon Dioxide on Coral Reefs. Science. 1999;284:118–20.

    Article  CAS  PubMed  Google Scholar 

  58. Goffredo S, Vergni P, Reggi M, Caroselli E, Sparla F, Levy O, et al. The skeletal organic matrix from Mediterranean coral Balanophyllia europaea influences calcium carbonate precipitation. PLoS One. 2011;6:e22338.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  59. Weiner S, Addadi L. Acidic macromolecules of mineralized tissues: the controllers of crystal formation. Trends Biol Sci. 1991;16:252–6.

    Article  CAS  Google Scholar 

  60. Hwang JS, Takaku Y, Momose T, Adamczyk P, Özbek S, Ikeo K, et al. Nematogalectin, a nematocyst protein with GlyXY and galectin domains, demonstrates nematocyte-specific alternative splicing in Hydra. Proc Natl Acad Sci U S A. 2010;107:18539–44.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  61. Moya A, Tambutté S, Bertucci A, Tambutté E, Lotto S, Vullo D, et al. Carbonic anhydrase in the scleractinian coral Stylophora pistillata: characterization, localization, and role in biomineralization. J Biol Chem. 2008;283:25475–84.

    Article  CAS  PubMed  Google Scholar 

  62. Furla P, Galgani I, Durand I, Allemand D. Sources and mechanisms of inorganic carbon transport for coral calcification and photosynthesis. J Exp Biol. 2000;203(Pt 22):3445–57.

    CAS  PubMed  Google Scholar 

  63. Allemand D, Ferrier-pagès C, Furla P, Houlbrèque F, Darwin C, Fitz C, et al. Biomineralisation in reef-building corals : from molecular mechanisms to environmental control. Gen Paleontol. 2004;3:453–67.

    Google Scholar 

  64. Oren U, Benayahu Y, Lubinevsky H, Loya Y. Colony integration during regeneration in the stony coral Favia favus. Ecology. 2001;82:802–13.

    Article  Google Scholar 

  65. Wang GL, Semenza GL. General involvement of hypoxia-inducible factor 1 in transcriptional response to hypoxia. Proc Natl Acad Sci U S A. 1993;90:4304–8.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  66. D’Angelo G, Duplan E, Boyer N, Vigne P, Frelin C. Hypoxia up-regulates prolyl hydroxylase activity: a feedback mechanism that limits HIF-1 responses during reoxygenation. J Biol Chem. 2003;278:38183–7.

    Article  PubMed  Google Scholar 

  67. Cummins EP, Berra E, Comerford KM, Ginouves A, Fitzgerald KT, Seeballuck F, et al. Prolyl hydroxylase-1 negatively regulates IkappaB kinase-beta, giving insight into hypoxia-induced NFkappaB activity. Proc Natl Acad Sci U S A. 2006;103:18154–9.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  68. Beeden R, Willis BL, Raymundo LJ, Page CA, Weil E. Underwater Cards for Assessing Coral Health on Indo-Pacific Reefs. 2008.

    Google Scholar 

  69. Sweet MJ, Bythell JC, Nugues MM. Algae as reservoirs for coral pathogens. PLoS One. 2013;8:e69717.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  70. Meyer E, Aglyamova G V, Wang S, Buchanan-Carter J, Abrego D, Colbourne JK, Willis BL, Matz M V: Sequencing and de novo analysis of a coral larval transcriptome using 454 GSFlx. BMC Genomics 2009, 10:219.

  71. Putnam NH, Srivastava M, Hellsten U, Dirks B, Chapman J, Salamov A, et al. Sea anemone genome reveals ancestral eumetazoan gene repertoire and genomic organization. Science. 2007;317:86–94.

    Article  CAS  PubMed  Google Scholar 

  72. Dunlap WC, Starcevic A, Baranasic D, Diminic J, Zucko J, Gacesa R, van Oppen MJH, Hranueli D, Cullum J, Long PF: KEGG orthology-based annotation of the predicted proteome of Acropora digitifera: ZoophyteBase - an open access and searchable database of a coral genome,14:509.

  73. Altschul SF, Madden TL, Schaffer AA, Zhang JH, Zhang Z, Miller W, et al. Gapped BLAST and PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res. 1997;25:3389–402.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  74. Tatusov RL, Fedorova ND, Jackson JD, Jacobs AR, Kiryutin B, Koonin EV, et al. The COG database: an updated version includes eukaryotes. BMC Bioinformatics. 2003;4:41.

    Article  PubMed Central  PubMed  Google Scholar 

  75. R Core Team. R: A Language and Environment for Statistical Computing. Vienna: Vienna R Found Stat Comput; 2014. ISBN 3–900.

  76. Love MI, Huber W, Anders S: Moderated estimation of fold change and dispersion for RNA-Seq data with DESeq2 Moderated estimation of fold change and dispersion for RNA-Seq data with DESeq2. bioRxiv 2014, doi:10.11.

  77. Kauffmann A, Gentleman R, Wolfgang H. arrayQualityMetrics–a bioconductor package for quality assessment of microarray data. Bioinformatics. 2009;25:415–6.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  78. Benjamini Y, Hochberg Y. Controlling the False Discovery Rate: A Practical and Powerful Approach to Multiple Testing. J R Stat Soc Ser B. 1995;57:289–300.

    Google Scholar 

  79. Horvath S. Weighted Network Analysis. Applications in Genomics and Systems Biology: Springer; 2011.

    Book  Google Scholar 

  80. Jombart T: adegenet: an R package for the multivariate analysis of genetic markers. Bioinformatics 2008:doi:10.1093/bioinformatics/btr521.

  81. Dixon P. VEGAN, a package of R functions for community ecology. J Veg Sci. 2003;14:927–30.

    Article  Google Scholar 

  82. Voolstra CR, Sunagawa S, Matz MV, Bayer T, Aranda M, Buschiazzo E, et al. Rapid evolution of coral proteins responsible for interaction with the environment. PLoS One. 2011;6:e20392.

    Article  PubMed Central  CAS  PubMed  Google Scholar 

  83. Kolde R. pheatmap: Pretty Heatmaps. 2013. R package version 0.7.7.

    Google Scholar 

Download references

Acknowledgements

We thank C Kenkel, D Stump, and C Palmer for sample collections and photographs. Funding for this study was provided by National Science Foundation grant DEB-1054766 to MVM.

Author information

Authors and Affiliations

  1. Institute for Cell and Molecular Biology, The University of Texas at Austin, Austin, USA

    Rachel M Wright

  2. Department of Integrative Biology, The University of Texas at Austin, Austin, USA

    Galina V Aglyamova & Mikhail V Matz

  3. Department of Zoology, Oregon State University, Corvallis, USA

    Eli Meyer

Authors
  1. Rachel M Wright
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Galina V Aglyamova
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Eli Meyer
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Mikhail V Matz
    View author publications

    You can also search for this author in PubMed Google Scholar

Corresponding author

Correspondence to Rachel M Wright.

Additional information

Competing interests

The authors declare that they have no competing interests.

Authors’ contributions

GVA prepared the libraries for sequencing, EM assembled the transcriptome, RMW conducted the tag-based RNA-seq data analysis. MVM supervised the study, annotated the transcriptome, and performed simulations. MVM and RMW wrote the manuscript. All read and approved the final manuscript.

Additional files

Additional file 1:

Tab-delimited table of unique transcript counts. Rows are genes, columns are samples. The letters in the column names identify disease state (AL, D, and H) and the number identifies individual colony.

Additional file 2:

Venn diagram of DEGs passing the FDR threshold of 10% for each of the contrasts. The largest number of DEGs resulted from the contrast of diseased (D) tissue to healthy (H) tissues, followed by the contrast of diseased to tissues ahead of the lesion (AL). No isogroups passed the FDR cutoff for the AL-healthy contrast.

Additional file 3:

Gene ontology categories enriched by genes up-regulated (red) or down-regulated (blue) in diseased compared to AL samples, summarized by molecular function (MF), biological process (BP), and cellular component (CC). The size of the font indicates the significance of the term as indicated by the inset key. The fraction preceding the GO term indicates the number of genes annotated with the term that pass an unadjusted p-value threshold of 0.05. The trees indicate sharing of genes among GO categories (the categories with no branch length between them are subsets of each other).

Additional file 4:

Gene expression heatmaps for annotated DEGs (adjusted p-value < 0.01) for the disease-healthy contrast. Rows are genes, columns are samples ordered as in the bottom panel: ahead-of-lesion (AL), healthy (H), and diseased (D). The color scale is in log2 (fold change relative to the gene’s mean). The trees are a hierarchical clustering of genes based on Pearson’s correlation of their expression across samples.

Additional file 5:

Gene expression heatmaps for annotated DEGs (unadjusted p-value < 0.01) for the disease-AL contrast. Rows are genes, columns are samples ordered as in the bottom panel: ahead-of-lesion (AL), healthy (H), and diseased (D). The color scale is in log2 (fold change relative to the gene’s mean). The trees are hierarchical clustering of genes based on Pearson’s correlation of their expression across samples.

Additional file 6:

Excel spreadsheet containing genes, annotations, and log fold changes for the disease-healthy contrast.

Additional file 7:

Excel spreadsheet containing genes, annotations, and log fold changes for the disease-AL contrast.

Additional file 8:

Excel spreadsheet containing genes, annotations, and log fold changes for the AL-healthy contrast.

Additional file 9:

Sample dendrogram and outlier heatmap for WGCNA. Sample clustering allows the visualization of how traits (health states and individual genotypes in this case) relate to samples. The dendrogram does not present any obvious outliers, but individual “D4” is called as an outlier based on a standardized connectivity test.

Additional file 10:

Heatmap of module-trait correlations. The strength of the correlations between traits (health states or individual corals) and gene coexpression modules are indicated by the intensity of color. The numbers in the cells give Pearson’s correlation between the module eigengene and the trait and the p-value according to the correlation test. Red boxes mark the three modules that are highly and specifically correlated to each of the health states.

Additional file 11:

Dependency between individual genes’ module membership (correlation with module’s eigengene) and significance for the disease state (correlation of the gene with the disease state). The grey region encompasses genes with both module membership and gene significance scores higher than 0.6. Pearson correlation values and the p-value of the correlation test are indicated in the lower-right hand of each scatterplot.

Additional file 12:

Simulation quality plots compare real and simulated DESeq2 datasets. (A) Size factors are nearly identical between real (dds) and simulated (sims) data sets. (B-C) Dispersion estimates of the real data set strongly agree with simulated dispersion estimates. (D) Dispersions are nearly identical between real (dds) and simulated (sims) data sets. (E-F) MA plots of log-fold change by mean expression value of real and simulated datasets are nearly identical.

Additional file 13:

Calculating empirical false discovery rates using simulated data. The x-axis is the Wald test p-value and the y-axis is the number of differentially expressed genes (DEGs) passing this p-value cutoff. The black line corresponds to DEGs discovered in the real dataset and red line corresponds to DEGs discovered in the simulated dataset (false positives). The vertical dashed line indicates the empirical 10% false discovery rate (FDR) cutoff. (A) Healthy-AL comparison. (B) Healthy-Disease. (C) AL-Disease comparison.

Additional file 14:

Shuffling sample-condition assignments dissolves module-trait relationships in WGCNA. On both panels, rows are module eigengenes (MEs) and the columns are traits (AL – ahead of lesion, H – healthy, and D –diseased). The color scale reflects the correlation of the module’s eigengene with the trait. The numbers in the cells are the Pearson r and the p-value of the correlation test (in parentheses). Strong correlations between modules of gene co-expression and health states were identified when the samples were assigned to correct conditions (left). These correlations disappeared when condition designations were shuffled among samples (right), providing evidence of a true biological relationship between the identified gene coexpression modules and health states.

Rights and permissions

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly credited. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Wright, R.M., Aglyamova, G.V., Meyer, E. et al. Gene expression associated with white syndromes in a reef building coral, Acropora hyacinthus . BMC Genomics 16, 371 (2015). https://doi.org/10.1186/s12864-015-1540-2

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-015-1540-2

Keywords

Comments

BMC Genomics

ISSN: 1471-2164

Contact us