Identifying genes and regulatory pathways associated with the scleractinian coral calcification process

Reef building corals precipitate calcium carbonate as an exo-skeleton and provide substratum for prosperous marine life. Biomineralization of the coral’s skeleton is a developmental process that occurs concurrently with other proliferation processes that control the animal extension and growth. The development of the animal body is regulated by large gene regulatory networks, which control the expression of gene sets that progressively generate developmental patterns in the animal body. In this study we have explored the gene expression profile and signaling pathways followed by the calcification process of a basal metazoan, the Red Sea scleractinian (stony) coral, Stylophora pistillata. When treated by seawater with high calcium concentrations (addition of 100 gm/L, added as CaCl2.2H2O), the coral increases its calcification rates and associated genes were up-regulated as a result, which were then identified. Gene expression was compared between corals treated with elevated and normal calcium concentrations. Calcification rate measurements and gene expression analysis by microarray RNA transcriptional profiling at two time-points (midday and night-time) revealed several genes common within mammalian gene regulatory networks. This study indicates that core genes of the Wnt and TGF-β/BMP signaling pathways may also play roles in development, growth, and biomineralization in early-diverging organisms such as corals.


INTRODUCTION
Coral calcification is a process of skeletal elements extension and is accompanied by the development of the coral's tissues. In general, the biomineralization and development processes are generated by the differential gene expressions with specialized cellular properties, as demonstrated in various organisms (Belcher et al., 1996;Gardner et al., 2011;Kuballa & Elizur, 2008;Reyes-Bermudez et al., 2009). These gene expression characteristics of many cellular processes, among them are biomineralization and tissue growth, which are coordinated by networks of regulatory genes.
The process of skeleton growth in corals has been suggested to be controlled by a specialized tissue called calicoblastic epithelium, which is comprised of ectodermal cells within the aboral coral layer, opposite the oral side (Gattuso, Allemand & Frankignoulle, 1999). The transport of calcium ions to the site of the growing skeleton is one feature for characterizing the controlled process of coral calcification. The mechanism behind the transport of calcium is under debate, and three pathways have been proposed: (i) active transcellular transport of calcium through calicoblastic cells, (ii) passive paracellular diffusion of calcium or seawater between calicoblastic cells, and (iii) a combination of transcellular and paracellular pathways (Reviewed in Allemand et al., 2011 andAllison et al., 2014).
Another aspect of the controlled process is the molecular mechanism which is estimated to include genes in several molecular pathways which aid the skeletogenesis. Various genes have been identified to be associated with the calcification mechanism and are involved in the transport of CaCO 3 reactants. The plasma membrane calcium-activated ATPase gene was isolated and cloned from a coral cDNA and ESTs (Zoccola et al., 2004;Zoccola et al., 2015), and its activity was measured and characterized (Furla, Allemand & Orsenigo, 2000;Ip, Lim & Lim, 1991;Isa, Ikehara & Yamazato, 1980). A voltage-dependent calcium channel was cloned and sequenced from cDNA, and its protein product was labeled within the calicoblastic ectoderm (Zoccola et al., 1999). Another gene that is involved in the transport of CaCO 3 reactants is the carbonic anhydrase, which has been identified and characterized by molecular biology, biochemistry and pharmacology techniques (reviewed in Bertucci et al., 2013). Several genes involved in the calcification mechanism have also been discovered as proteins found in the coral's organic matrix with protein extraction, separation and characterization techniques (Drake et al., 2013;Ramos-Silva et al., 2013), or with RNA extraction and transcriptome assembly methods (Mass et al., 2013;Moya et al., 2012). The coral's organic matrix proteins are characterized as highly acidic proteins that are abundant with aspartic and glutamic acids. Still, the molecular aspects of the skeletogenesis process are not fully understood.

Figure 1 Experimental layout design.
Eight fragments of Stylophora pistillata were incubated in calcium seawater (SW) concentrations of seawater plus 100 mg/L calcium (SW + 100) and plus 200 mg/L calcium (SW + 200) during daytime (11:00) and during night-time (23:00). In each sampling time, four replicates were analyzed for calcification rates and four replicates from SW and SW+100 treatments, for microarray analysis.
The aim of this study is to identify potential key genes and pathways that are related to the calcification apparatus. This will be achieved by exposing the coral to seawater with high calcium concentrations to induce an alteration in gene expression profiles, identified using a DNA microarray analysis. In addition, we explored the differential genes expressions between day and night in the presence of high calcium concentration seawater in order to identify the presence of networks of regulatory genes in a basal metazoan.

Coral collection
Scleractinian coral colonies of Stylophora pistillata were collected by SCUBA diving from 3-6 m depth at the coral reef near the Inter-University Institute (IUI) in Eilat, Aqaba Gulf, Red Sea. The branching coral colonies were fragmented to 5-10 cm long nubbins and transferred to an outdoor shaded running seawater table at the IUI for a month acclimation and recovery before the experiment started. The Israeli Nature and National Parks Protection Authority approved the collection of corals in this study (permit No. 2009-35556).

Calcification rates measurements under different seawater calcium concentrations
High calcium seawater was prepared by adding either 100 mg/L or 200 mg/L calcium (added as CaCl 2 ·2H 2 O) to filtered seawater (0.45 µm).
Four coral fragments of S. pistillata, for each treatment (Fig. 1), were incubated in running seawater tables for 1.5 h for two time points: midday (11:00) and night-time (23:00). The fragments were incubated separately in three seawater calcium concentrations: natural calcium concentration (SW), natural calcium concentration added with 100 mg/L (SW + 100) and natural calcium concentration added with 200 mg/L (SW + 200).
The calcification rate values (µmol CaCO 3 h −1 cm −2 ) were calculated according to the equation: where A T is the difference in Total Alkalinity (A T ) measured between the beginning and the end the incubation period, V is the volume of the chamber or the coral fragment, T is the duration of the incubation and A is the coral surface area. The coral's surface area (A) was calculated by considering the branched fragments as perfect cylinders, hence the height and diameter were measured and used for the cylinder-coral fragment surface area calculation (Mass et al., 2007). Such short period incubation was already been successful for distinguishing differences of calcification rates between treatments (Cohen, Dubinsky & Erez, 2016;Gutner-Hoch et al., 2016).

Microarray experiment set-up and coral fragments sampling
Four coral fragments of S. pistillata were incubated for 1.5 h in each of the two different calcium concentrations (SW, SW + 100), and for two time points: in the mid-day (11:00) and at night-time (23:00). After incubation, coral fragments were sampled, snap froze in liquid nitrogen and stored in −80 C until further handled. The coral fragments that have been subjected to SW + 200, were removed from gene expression analysis because no statistically significant was found between the calcification rates of SW + 100 and SW + 200 treatments (Fig. 2).

RNA extraction and microarray design
Total RNA was isolated from the S. pistillata fragments by using TRIzol R reagent (Invitrogen, Carlsbad, CA, USA) according to manufacturer's instructions. RNA quantity and integrity were assessed with a NanoDrop ND-1000 spectrophotometer and an Agilent 2100 Bioanalyzer, respectively. Total RNA was labeled and hybridized against custom S. pistillata microarray; an Agilent two-color gene expression microarray platform with 8 × 15 K probe per slide; microarray samples included four SW+100 and four SW (control) samples for each time point (11:00 and 23:00). Each treatment was combined with a differently labeled control sample from the same time point for hybridization to each compartment of the microarray slide. Oligonucleotide probes (60-mers) were designed based on approximately 12,000 genes predicted to encode proteins retrieved from a recent de novo assembly of 454-sequenced EST libraries of S. pistillata (Karako-Lampert et al., 2014). Labeling and hybridization were conducted using the Agilent Low Input Quick Amp Labeling Kit according to the manufacturer's instructions. The intensity of the emitted fluorescence from a target spot on the array was detected using an Agilent G2565BA microarray scanner. The raw data as well as the processed data of the microarray were deposited under Gene Expression Omnibus (GEO) (http://www.ncbi.nlm.nih.gov/geo) accession number GSE87159. Stylophora pistillata EST data are also stored at the Cnidarian Database of Centre Scientifique de Monaco: http://data.centrescientifique.mc/CSMdata-home.html.

Microarray results analysis
The data from all arrays were first subjected to background correction and LOESS withinarray normalization using Agilent Feature Extraction software (version 9.5.1.1; Agilent Technologies, Santa Clara, CA, USA). The remaining analyses were performed in Partek R Genomics Suite software (version 6.6; Partek Inc., St. Louis, MO, USA). Data from three biological replicates were used to perform a one-way ANOVA. The quintile normalized data were analyzed to identify genes with significantly up-or down-regulated expression (FDR p-value < 0.05) with an arbitrary cutoff of at least a two-fold change. As the microarray contains a diversity of genes originated from different organisms, all genes were transformed to Homo sapiens ones. Networks of highly interconnected proteins were generated using the STRING (Search Tool for the Retrieval of Interacting Genes, Heidelberg, Germany) 9.0 database (Szklarczyk et al., 2011) and GeneMania (Warde-Farley et al., 2010). The gene function clustering in GeneMania was predicted by using the Gaussian field label propagation algorithm. For this study the Homo sapiens dataset of functional association network was used in a network-weighting method of 'Gene-Ontology (GO) based weighting, biological process based. The names of samples within the microarray analysis were labeled C11 and T11 as for control (reference) and treatment at daytime (11:00), respectively. The samples C23 and T23 were labeled as for control (reference) and treatment for night-time (23:00), respectively.

Semi-quantitive real time PCR
The relative expression of genes was analyzed through Qiagen Corbett machine, using a two step real-time PCR. The gene expression from every sample was amplified in triplicates, with a single volume of 10 µl containing: 3 µl of diluted cNDA with 0.5 µl of primer set (see Table 1 for primer sequences), 1.5 µl double distilled water (DDW) and 5 µl of GoTaq qpcr Master Mix Kit (Promega, USA). Cycle threshold values for each time point were compared to the internal reference genes, according the 2 − CT method (Livak & Schmittgen, 2001). The gene β-actin was used for normalization in all RT-PCR analyses (Levy et al., 2007).

Statistical analyses
For statistical analyses, one-way ANOVA was used and followed by Tukey HSD (honestly significant differences) to assess the differences in experimental treatments of calcification rates and genes' real-time PCR results. All statistical analyses were conducted using SPSS 20.0 (IBM, USA), and the results were considered statistically significant at p < 0.05.

Microarray expression analysis
The Stylophora pistillata gene expression profile was analyzed with Partek Genomics Suite TM software (Fig. 3A) in a principal component analysis (PCA). The analysis showed that experiment replicates were profiled in the 3D PCA plot into two different groups-treated coral fragments (T11 and T23, which are fragments that were incubated in SW+100 mg/L at 11:00 and at 23:00) and second group-reference coral fragments (C11 and C23, which are fragments that were incubated in natural SW at 11:00 and at 23:00). The microarray analysis revealed that 2,768 transcripts were differentially expressed (fold > 2) across all of the treatment groups. Among the differentially expressed gene targets, 1924 were found comparing treatment corals versus reference corals from night-time (T23 vs. C23). Comparing treatment corals with reference coral from daytime (T11 vs. C11) revealed 1,957 differentially expressed genes. In addition, comparing both lists of treatments corals with reference corals from night-time (T23 vs. C23) and daytime (T11 vs. C11) 1,113 common genes were identified (Fig. 3B). Within the 811 unique differentially expressed genes during night-time (which did not overlap with daytime) 267 genes were up-regulated and 544 were down-regulated. Within the 844 unique differentially expressed genes during daytime (which did not overlap with night-time) 437 genes were up-regulated and 407 were down-regulated.
In order to focus on genes that share similar biological processes and correlate with the calcification process, the analysis continued for unique up-regulated gene transcripts and total up-regulated gene transcripts that were differentially expressed during daytime (437 genes for unique and 1,333 for total), and during night-time (267 genes for unique and 1,171 for total). The total (found also to be up-regulated during the ''opposite'' time point) and unique genes were clustered into biological processes descriptions with GeneMania website analysis portal.
Among the total gene transcripts that were differentially expressed during daytime six biological processes were found: Wnt signaling pathway, cell junction organization, SMAD binding, regulation of cell morphogenesis, TGF-β signaling pathway, and amino acid synthesis (Fig. 4A). From the unique genes during daytime four biological processes were found: amino acid synthesis, metabolism, Smad signaling, and cartilage and bone biogenesis ( Fig. 4A and Tables 2 and 3). All biological processes that were analyzed from total and unique genes were found with high significance of FDR p-value <0.05, except of ''connective tissue development'' cluster within the biological process cartilage and bone biogenesis, with FDR p-value of 0.068.
Among the total gene transcripts that were differentially expressed during night-time four biological processes were found: protein folding, mitochondrial membrane, innate immune system, and protein targeting (Fig. 4B). From the unique genes during night-time three biological processes were found to be associated with ribosomal metabolism, peptide synthesis/modification/processing, and ATP synthesis and oxidative phosphorylation. (Fig. 4B and Tables 2 and 3). All nighttime biological processes were found with high significance of FDR p-value <0.05. Venn diagram representing the overlap of the total of 2,768 genes of the total 2,768 differentially expressed genes at treatment samples from 11:00 and from 23:00. Among the 844 unique genes that were significantly up-or down-regulated expression during the day, 437 were showed up-regulation of their expression and 407 were down-regulated. Among the 811 unique genes that were significantly up-or downregulated expression during the night, 267 were showed up-regulation of their expression and 544 were down-regulated.
(List of total genes for C11 vs. C23, T11 vs. C11, T23 vs. C23, and T11 vs. T23 can be found in File S2; Biological processes of these total gene transcripts analyses can be found in File S3).
In order to expand the analysis of calcification associated genes, the gene expressions of S. pistillata's orthologs of galaxin and galaxin-like were analyzed through semi-quantitative

Figure 4 Biological processes according to differential expression of total and unique gene transcripts (genes that were up-regulated only during treatment) from treated vs. control during daytime (A), and from treated vs. control during night-time (B).
real-time PCR (Figs. 5A-5B and File S1). The galaxin gene has displayed a high relative expression during midday (relative expression of 59.76 ± 0.89 at 11:00, C11 in Fig. 5A) with statistically significant of p < 0.0005 (one-way ANOVA). While in treatment with addition of 100 mg/L calcium, relative expression of galaxin was dropped to expression level similar to expression during night-time (relative expression of 17.72 ± 2.61 at 23:00, C23 in Fig. 5A). The expression of gene galaxin-like was enhanced when treated with addition of 100 mg/L of calcium with relative expression of 14.14 ± 4.83 compared to relative expression of 3.48 ± 0.96 of the expression in the control corals at 11:00 (statistically Table 2 List of biological processes and number of gene transcripts within each biological process, divided to daytime and night-time. All biological processes were found with FDR p < 0.05.

Ribosomal Metabolism
Protein localization to endoplasmic reticulum 10 SRP-dependent cotranslational protein targeting to membrane 9 Protein targeting to ER 9 Cotranslational protein targeting to membrane 9 Establishment of protein localization to endoplasmic reticulum 9 Peptide Synthesis Modification Processing

Mitochondrial inner membrane 11
Respiratory chain 6 Respiratory electron transport chain 7 Electron transport chain 7 Mitochondrial ATP synthesis coupled electron transport 5 ATP synthesis coupled electron transport 5
Following the microarray analysis, the gene smad1, which was identified to be upregulated during daytime and associated with the Smad signaling pathway was used to validate the microarray results (Fig. 5C). Semi-quantitative real-time PCR results of smad1 displayed that its expression was elevated during the day to 3.00 ± 0.48 (T11 in Fig. 5C) when seawater had added calcium, compared to relative expression of 1.56 ± 0.2 in control corals (C11 in Fig. 5C; statistical significance p < 0.05, one-way ANOVA, n =3). However, there was no difference in the relative expression of the gene during the night-time (C23 and T23 in Fig. 5C).

DISCUSSION
Genetic pathways and regulatory networks are important for numerous cellular functions which enable the organism to survive in the environment through different physiological processes (Karlebach & Shamir, 2008). Exposing a scleractinian coral to a seawater medium with high concentration of calcium ions caused an excessive expression of ortholog genes ( Fig. 4A and 4B) which are assumed to be associated with processes of calcification, growth, and development (Tables 2 and 3). Moreover, using this approach we identified genes which are likely involved in cartilage and bone biogenesis in vertebrates.
Among the genes that were identified to be up-regulated are orthologs that encode for human collagen proteins: ADAMTS2, ADAMTS7 and COL2A1. These proteins are  galaxin-like (B), and smad1 (C). Relative gene expression was determined by the 2 − method and normalized by β-actin's gene expression. Means followed by the same letter are not significantly different, one-way ANOVA followed by Tukey's post-hoc test. Error bars represent SD. n = 3. responsible for providing the infrastructure of a cartilaginous matrix for bone and cartilage growth in mammals. The cartilaginous matrix is composed of collagen and proteoglycans, and is generated by cartilage cells that are called chondrocytes. In addition, one of the genes identified-an ortholog of human COL2A1 is commonly used as a chondrocytic maker in the development of cartilage (Zaucke et al., 2001).
Another identified ortholog gene that emerged from the analysis is bmp (Bone Morphogenetic Protein) which is putatively associated with cartilage and bone promotion. The BMP protein acts as an activator for cell-surface receptors, and initiates cellular processes and signaling cascades. This protein also participates in another gene signaling pathway that was up-regulated in this study: the Smad signaling pathway. Smad genes are part of a transforming growth factor beta superfamily (TGF-β). The Smad genes are translated into multifunctional peptides, which play a role in the regulation of cell-type specification, proliferation, and differentiation of many cell types within a variety of species (Raftery & Sutherl , 1999;Sekelsky et al., 1995). Smad proteins are recruited after being phosphorylated by BMP receptors and assembled into a complex. The complex then migrates into the cell nucleus, and activates the transcription of specific target genes which are involved in physiological and pathological bone formation (Chen, Deng & Li, 2012;Katagiri & Tsukamoto, 2013;Miyazawa et al., 2002).
The ortholog genes that were found to be associated with cartilage development, along with bmp which acts as a key member gene, are known to be involved in the Wnt signaling pathway (Kim et al., 2013;Lin & Hankenson, 2011). This signaling pathway was found to be important in embryonic bone development (Wan et al., 2013), and in osteoblasts (bone cells) regulation (Glass et al., 2005;Harada & Rodan, 2003). The Wnt and TGF-β signaling pathways were identified in the gene expression analysis as being up-regulated during the day (Fig. 4). Lipid and nitrogen metabolism genes also exhibited up-regulation, which likely corresponds to the high metabolic requirements that support the coral growth and development. During the night, coral fragments treated with high calcium concentrations showed active genes that are involved in the mitochondrial synthesis of ATP and ribosomal metabolism. Ribosomal biogenesis and the translation of proteins included within their processing, folding and degradation have shown to be activated during the night and controlled by the circadian clock in past studies (Jouffe et al., 2013;Panda et al., 2002).
The high calcium concentration treatment during the day revealed genes associated with cell junction organization process exhibiting up-regulation (Fig. 4). Indication of reorganization of the cell junctions, which mediate adhesion and communication between tissues and cells may suggest a possible role in growth and development of the organism. Goldberg (2001) and Tambutté et al. (2007) have previously reported the occurrence of anchoring cells, also named desmocythe cells, which were suggested to secure the coral tissue to the extracellular matrix and/or to the skeleton.
Along with the cellular genes emphasized from the microarray analysis, two genes that are suspected to be involved in the scleractinian coral calcification process were tested for expression (Figs. 5A and 5B). These two genes in S. pistillata, galaxin and galaxin-like, were found to be similar to the Acropora millepora's galaxin and galaxin-like1 genes (Karako-Lampert et al., 2014), and have been characterized to be associated with the calcifying organic matrix (Fukuda et al., 2003). Both genes were affected by the calcium concentration treatment, during the day, but not at night. While S. pistillata's galaxin-like gene was up-regulated during the daytime treatment, the galaxin gene was down-regulated, and was similar to night-time expression levels. These results are interesting because galaxin was first found due to its high levels in the organic matrix of scleractinian coral, and it was characterized as a soluble protein (Watanabe et al., 2003) and localized along the calcifying septa (Reyes-Bermudez et al., 2009). According to the findings of Reyes-Bermudez et al. (2009), the galaxin-like gene is considered to be a gene expressed only during planula stages (free swimming planula, pre-settled and post-settled planula in the case study of A. millepora), and not in the adult stages. Our results suggest that the role of galaxin genes in S. pistillata is possibly different from that in A. millepora.
The evidence of ortholog genes expression which likely associated with Smad proteins, Wnt, TGF-β signaling pathways, and cartilage formation in S. pistillata supports previous reports suggesting that gene regulatory pathways of development and biomineralization processes across the metazoan evolution are considerably conserved (Jackson et al., 2007;Marin et al., 2008). This study suggests that the Wnt signaling pathway cooperates with the TGF-β signaling pathway, similar to the mammalian cellular mechanism (Chen, Deng & Li, 2012). Both pathways appear to be responsible for the differentiation and proliferation of cell type similar to bone mass and osteoblasts (Fig. 6). Similarities between invertebrates and vertebrates regarding the elements associated with calcification process come in agreement with previous reports of Atlan et al. (1997) andLiao et al. (1997). Those studies showed that nacre implants from a bivalve mollusk were successfully incorporated and transplanted into a mammalian bone tissue due to similarity between the growth proteins present in nacre and those in bone. Similarity between the osteoblasts and the calicoblastic epithelium cells of scleractinian coral was also shown in the localization of the protein BMP in the tissue of S. pistillata (Zoccola et al., 2009), and for its role in the coral biomineralization (Mass et al., 2016;Zoccola et al., 2009).
Exposing S. pistillata to seawater with an elevated calcium concentration resulted in up-regulation of gene expression during the day relative to the night. In addition, the effect of the high calcium concentration was stronger in the calcification rates during day compared to the night (Fig. 2). Previous studies showed that the coral calcification process has a diel rhythmic cycle of increasing rates towards midday, and then decreasing towards dusk (Gutner-Hoch et al., 2016;Schneider et al., 2009). This could be an indication of a controlled expression of genes responsible for calcification during daytime, as evidenced by the up-regulation of genes associated with the cartilage and bone biogenesis process. Another possible explanation is the activity of the coral endo-symbionts during daytime which acts as important mean of energy to the coral's cellular functions and physiology (Colombo-Pallotta, Rodríguez-Román & Iglesias-Prieto, 2010). Hence, it could be suggested that the small differences in the calcification rates observed between 100 mg/L and 200 mg/L calcium additions during the day could be a result of relatively mediocre photosynthetic activities of the endo-symbionts, but this should be further studied.
Uncovering regulatory genes involved in the scleractinian biomineralization process and their coded proteins is truly a challenge because genetic manipulation methods of mutagenesis are not yet accessible in scleractinian corals. However, the use of antagonists for the TGF-β/BMP and Wnt pathways associated with the coral calcification process has not been tested. Studying coral calcification with antagonists can provide further information and knowledge to better understand the role and function of genes that participate in the cellular control mechanism, and in the coral organic matrix.