The Genome of the Cauliflower Coral Pocillopora verrucosa.

Climate change and ocean warming threaten the persistence of corals worldwide. Genomic resources are critical to study the evolutionary trajectory, adaptive potential, and genetic distinctiveness of coral species. Here, we provide a reference genome of the cauliflower coral Pocillopora verrucosa, a broadly prevalent reef-building coral with important ecological roles in the maintenance of reefs across the Red Sea, the Indian Ocean, and the Pacific Ocean. The genome has an assembly size of 380,505,698 bp with a scaffold N50 of 333,696 bp and a contig N50 of 75,704 bp. The annotation of the assembled genome returned 27,439 gene models of which 89.88% have evidence of transcription from RNA-Seq data and 97.87% show homology to known genes. A high proportion of the genome (41.22%) comprised repetitive elements in comparison to other cnidarian genomes, in particular in relation to the small genome size of P. verrucosa.

Significance

The ongoing destruction of coral reefs worldwide heightens the need to better understand the genomic underpinnings of coral resilience. One critical resource to assist such efforts is the generation and provision of reference genomes to enable population genomics approaches and adaptive evolution studies. Here, we generated the genome of the common coral Pocillopora verrucosa, a broadly studied and ecologically important species. Our study demonstrates that reference genomes can be obtained utilizing improved assembly algorithms and the availability of highly assembled genomes from close-by species to fill critical gaps for species that lack genomic resources.

Introduction

Coral reefs are among the most biodiverse ecosystems on earth, providing habitat for about a third of all marine species (Spalding et al. 2001; Plaisance et al. 2011). Besides their ecological importance, coral reefs provide numerous goods and ecosystem services to millions of people (Moberg and Folke 1999; Cesar et al. 2003). The well-being of these ecosystems relies on the health of reef-building corals (order Scleractinia), which comprise the foundation species and main ecosystem architects of coral reefs. Despite their importance, corals are severely threatened by local (overfishing, pollution, industrialization) and global (climate change, ocean warming) anthropogenic impacts that continue to decrease coral cover at an alarming rate (IPCC 2018). Therefore, better knowledge on coral biology is critical to conceive strategies to mitigate coral reef demise.

The cauliflower coral  Pocillopora verrucosa (Ellis and Solander 1786; NCBI Taxonomy ID: 203993; supplementary fig. S1, Supplementary Material online) is one of eight species currently described in the genus Pocillopora (Schmidt-Roach et al. 2014) and belongs to the family Pocilloporidae together with the genera Madracis, Seriatopora, and Stylophora. P. verrucosa is widely distributed and is prevalent in the Red Sea, across the Indian Ocean, and across the Pacific Ocean, up to the Eastern Tropical Pacific, colonizing an extensive range of depths that cover shallow to mesophotic settings (Schmidt-Roach et al. 2014; De Palmas et al. 2018; Soto et al. 2018). Pocillopora verrucosa is an important reef framework builder and has an important ecological role as an early colonizer, aiding the recovery of reefs after natural disturbances (Pearson 1981; Bianchi et al. 2006). Despite the cosmopolitan distribution of P. verrucosa, which makes it a model for experimental and biological studies covering the fields of histology (Kruger and Schleyer 1998; Hirose et al. 2000), physiology (Ziegler et al. 2014; Sawall et al. 2015; Edmunds and Burgess 2016), phylogenetics (Flot et al. 2008, 2011; Pinzón et al. 2013; Schmidt-Roach et al. 2014), and population genetics (Combosch and Vollmer 2015; Robitzch et al. 2015), a reference genome is not available.

Here, we report on the sequencing, assembly, and annotation of the genome of the cauliflower coral P. verrucosa from the Red Sea (supplementary fig. S1, Supplementary Material online), which we obtained through initial assembly of a high-throughput single paired-end short-read library and subsequent reference scaffolding to the available genome of Pocillopora damicornis (Cunning et al. 2018). The generated draft genome is ∼380 Mb in size (scaffold N50 of 333 kb, contig N50 of 75 kb) and features a total of 27,439 predicted protein-coding genes. The genome of P. verrucosa will be of value for comparative studies (Voolstra et al. 2015, 2017; Bhattacharya et al. 2016; Cunning et al. 2018), as a genomic reference to enable population genomic approaches (e.g., RAD-Seq), and as a foundation for molecular biology studies (e.g., RNA-Seq).

Materials and Methods

Coral Sampling and Isolation of Genomic DNA

On April 26, 2018, a medium-size fragment (10 cm length) of a coral colony of P. verrucosa was collected at 6 m depth in the Red Sea Al Fahal forereef (22°15.100N, 38°57.386E) and transferred to the Coastal and Marine Resources Core Lab (CMOR) aquaria facilities in KAUST, Saudi Arabia (supplementary extended materials and methods, Supplementary Material online). The Saudi Coastguard Authority issued sailing permits to the site that included coral collection, following the Nagoya protocol. Coral genomic DNA (gDNA) was extracted using the Blood & Tissue kit (Qiagen, Hilden, Germany).

Pocillopora verrucosa Lineage Determination

Previous phylogenetic analyses have identified two genetic lineages (type 3 and type 7) within P. verrucosa in the Red Sea (Pinzón et al. 2013; Schmidt-Roach et al. 2014) with type 3 more commonly found than type 7. To determine the genetic lineage of our P. verrucosa specimen, we used the primers FATP6.1 (5′‐TTTGGGSATTCGTTTAGCAG‐3′) and RORF (5′‐SCCAATATGTTAAACASCATGTCA‐3′) (Flot et al. 2008) to amplify the mitochondrial open reading frame (mORF) region. Phylogenetic analysis showed that our P. verrucosa specimen clustered within the mitochondrial open reading frame lineage type 3 sequences (Pinzón et al. 2013), Clade 2 (Schmidt-Roach et al. 2014). Hence, the genome represents the most common Red Sea P. verrucosa lineage (supplementary fig. S2, file S1, and extended materials and methods, Supplementary Material online).

Physical Genome Size Estimation

We determined the physical genome size by measurement of the P. verrucosa nuclei size using chicken erythrocyte nuclei as a reference (DNA QC Particles kit, BD Biosciences, San Jose, CA; supplementary extended materials and methods and fig. S3, Supplementary Material online). The coral genome size estimation was based on the diploid DNA content of chicken erythrocytes of 2.5 pg±0.04 per cell and was calculated as follows: sample genome size [pg]=1.25 x/y (x, fluorescence intensity of unknown sample; y, fluorescence intensity of chicken erythrocyte nuclei). After calculating the mean DNA content per copy of genetic information (1C), the genome size was determined by considering that 1 pg DNA equals 978 Mb (Doležel et al. 2003).

Sequencing Library, Read Filtering, Genome Assembly, and Scaffolding

A single gDNA library was constructed using the NEBNext Ultra II DNA Library Prep Kit for Illumina according to the manufacturer’s instructions (New England Biolabs, Ipswich, MA). The paired-end library was sequenced on the Illumina HiSeq2500 platform (rapid run—500 cycles) at the KAUST Bioscience Core Lab (BCL). The sequencing yielded 141,993,203 read pairs (>174× coverage assuming a genome size of ∼400 Mb), which were trimmed, quality assessed, and filtered to remove unwanted sequences according to the pipeline detailed in the Supplement.

The filtered paired-end reads were assembled with DISCOVAR de novo (Weisenfeld et al. 2014). To assess the level of heterozygosity, the filtered paired-end reads were used to obtain k-mer counts of lengths 25 and 31 using the software jellyfish version 2.2.6 (Marçais and Kingsford 2011) (details available in the Supplement). Based on k-mer distributions of lengths 31 and 25, the obtained heterozygosity rate was 1.21% and 1.32% (supplementary table S1, Supplementary Material online), respectively. This is in line with the high heterozygosity rates (1–2%) reported for other coral genomes (Bellis et al. 2016; Helmkampf et al. 2019; Robbins et al. 2019). Consequently, we applied a hierarchical filtering strategy as follows: 1) Circular scaffolds as well as contigs of mitochondrial origin were identified and removed, 2) contigs containing potential sequences from dinoflagellate, bacterial, or viral origin were identified with BLASTN (Altschul et al. 1990; Camacho et al. 2009) applying a 90% identity cutoff over 50% query length and removed (Voolstra et al. 2017). After that, the assembly was scaffolded using a reference-based approach employing CSAR (Chen et al. 2018) and the available genome of the closely related species P. damicornis (Cunning et al. 2018). Gaps were filled with GapFiller version 1.11 (Boetzer and Pirovano 2012) using the filtered reads of the paired-end sequencing library (see above). As a final filtering step 3), the scaffolded assembly was processed with the Haplomerger2 pipeline v20180514 (Huang et al. 2017) to further improve the assembly.

Basic summary statistics of the initial DISCOVAR de novo assembly as well as the putative haploid final genome assembly were estimated by QUAST version 5.0.2 (Gurevich et al. 2013) (supplementary table S2, Supplementary Material online). QUAST version 5.0.2 was also used to estimate genome statistics for the genomes of P. damicornis (Cunning et al. 2018) and Stylophora pistillata (Voolstra et al. 2017). Both genomes were downloaded from the reefgenomics.org database (Liew et al. 2016) at http://pdam.reefgenomics.org/ and http://spis.reefgenomics.org, respectively.

Genome Assembly Completeness Analysis

Completeness of the P. verrucosa genome was assessed by searching for 978 universal metazoan Single-Copy Orthologs (metazoa_odb9) using BUSCO version 3 in mode 'genome'. For comparative purposes, we performed the same genome completeness assessment for the genomes of P. damicornis and S. pistillata (supplementary table S3, Supplementary Material online).

Annotation of Repetitive Elements

Repetitive elements in the genome assembly of P. verrucosa were identified de novo using RepeatScout version 1.0.5 (Price et al. 2005) with an l-mer size of 16 bp and default settings (supplementary extended materials and methods and table S4, Supplementary Material online).

Reference Transcriptome Sequencing and Assembly

An RNA-Seq library was generated from 500 ng of total RNA using the TruSeq Stranded mRNA Library prep Kit (Illumina, San Diego, CA) as per the manufacturer’s instructions (further details in the Supplement). The constructed library (median size 340 bp) was indexed and paired-end sequenced on the Illumina HiSeq4000 platform (300 cycles) at the KAUST BCL. RNA sequencing yielded 121,056,948 read pairs (i.e., 242 mio. paired-end reads), which were trimmed and further processed as detailed in the Supplement.

Gene Model Prediction

Coral transcripts from the reference transcriptome (n = 49,384, supplementary extended materials and methods, Supplementary Material online) were mapped to the genome assembly and filtered by PASA version 2.3.3 (Haas et al. 2008) to create a training set for AUGUSTUS version 2.5.5 (Stanke and Morgenstern 2005).

Protein Set Annotation

The final set of predicted proteins was generated by selecting the longest isoform per gene. The protein sequences were annotated following a hierarchical approach using the UniProt (SwissProt and TrEMBL) and the NCBI “nr” databases sensu Baumgarten et al. (2015). This approach led to the annotation of 19,427 proteins in SwissProt, 7,145 proteins in TrEMBL, and 283 proteins using the NCBI “nr” database; 584 proteins had no hits in either of the reference databases (supplementary extended materials and methods and file S2, Supplementary Material online).

Results and Discussion

Genome Size Estimation

We estimated the genome size of P. verrucosa using FACS of propidium-iodide-stained nuclei in addition to a k-mer based approach. The physical estimate yielded a genome size of ∼407 Mb (supplementary fig. S3, Supplementary Material online) showing a good agreement with the k-mer based estimate, which obtained a size of 406 and 418 Mb for k-mers 25 and 31, respectively (supplementary table S1, Supplementary Material online). This genome size is larger than the estimated genome size of P. damicornis (349 Mb assessed using k-mer counts, Cunning et al. 2018) and smaller than that of S. pistillata (434 Mb assessed using FACS, Voolstra et al. 2017), which are the two other coral species in the family Pocilloporidae with published genomes.

Genome Assembly and Statistics

We obtained 141,993,203 paired-end read pairs (2×250 bp) from a single gDNA library. About 88% of read pairs were kept after filtering and used for a de novo genome assembly. Read coverage was estimated to be 153× based on the physical genome size estimate (407 Mb), which was reported to be in the optimal range for DISCOVAR de novo genome assemblies (Love et al. 2016). After assembly, we applied a filtering strategy to remove nonnuclear and noncoral contigs before scaffolding using the reference genome of the closely related coral P. damicornis. To address the issue of high levels of heterozygosity, we applied the Haplomerger2 pipeline to consolidate incompletely merged loci/alleles and produce a haploid reference genome assembly with a total length of ∼380 Mb assembled in 25,605 contigs and 18,268 scaffolds with a corresponding contig N50 of 75,704 bp and scaffold N50 of 333,696 bp, respectively (table 1). The assembly statistics of the scaffolded haploid genome compare well to the genome assemblies of the available genomes in the family Pocilloporiadae, namely P. damicornis (contig N50 = 25,987 bp, scaffold N50 = 326,133 bp) and S. pistillata (contig N50 = 20,518 bp, scaffold N50 = 457,453 bp).

Table 1

Summary Statistics of the Pocillopora verrucosa Genome Assembly

Genome Assembly Statistics	Length		Number
	Contig	Scaffold	Contig	Scaffold
Total (bp)	379,998,154	380,505,698	25,605	18,268
Max (bp)	571,091	2,095,917	—	—
≥5,000 (bp)	353,178,732	358,580,702	8,537	3,932
N50 (bp)	75,704	333,696	1,401	326
N75 (bp)	30,836	125,371	3,336	766
N’s per 100 kb	7.89	134.40

To compare the assembly before and after reference genome scaffolding, we computed basic statistics of the DISCOVAR de novo assembly and of the final reference scaffolded and filtered assembly (after DISCOVAR assembly) (supplementary table S2, Supplementary Material online). As expected, the scaffolded haploid genome was a substantial improvement over the de novo assembly, highlighting the utility of having genomes of closely related species available for the generation of reference genomes from so far unconsidered species.

Genome assembly completeness was assessed by searching for 978 universal metazoan single-copy orthologs using BUSCO. We could identify 902 (92.22%) complete metazoan single-copy orthologs (of which 3.17% were duplicated), 18 orthologs (1.84%) were present but fragmented, and 58 orthologs (5.93%) were missing. Comparison to P. damicornis and S. pistillata showed that P. damicornis had 862 (88.14%) and S. pistillata had 861 (88.03%) complete metazoan single-copy orthologs (supplementary table S3, Supplementary Material online). The number of fragmented universal metazoan single-copy orthologs is comparable among the three pocilloporid species, whereas the number of missing single-copy orthologs is smaller in P. verrucosa (58, 5.93%) than in P. damicornis (87, 8.90%) and S. pistillata (89, 9.10%).

Genomic Repeats

The identification and annotation of repetitive elements showed that 41.22% of the genome comprised repetitive elements. Out of the 15,905 repeat motifs identified, we could annotate 4,915 motifs, which were further classified in six groups encompassing 63 superfamilies/clades (supplementary table S4, Supplementary Material online). Transposable elements (TEs), such as DNA transposons, LTR retrotransposons, and non-LTR retrotransposons (among others) comprised the largest portion of the P. verrucosa genome (17.22%), followed by unclassified interspersed repeats (13.47%) and simple repeats (10.28%). Given that unclassified interspersed repeats comprise putative species-specific TEs, we considered them TEs for the purpose of comparison between P. verrucosa and S. pistillata, as previously done for other coral genomes (Ying et al. 2018). Overall, we found a slightly higher proportion of TEs in P. verrucosa (30.69% of genome size) in comparison to S. pistillata (28.43% of genome size), despite the smaller assembly size of P. verrucosa. This was unexpected given the known positive correlation between TEs and genome size (Kidwell 2002; Hua-Van et al. 2005; Biscotti et al. 2015), although some margin of error exists.

At the moment, we can only speculate on the significance of the increased TE content in the genome of P. verrucosa in comparison to S. pistillata. TE expansion can denote signatures of genus radiation (Wong et al. 2019), and TE activity has been associated with response to environmental stress and phenotypic plasticity in plants (Negi et al. 2016), arguably of relevance for corals in the context of climate change and ocean warming. From this perspective, an increased TE content in P. verrucosa fits with the species’ broad geographical and depth distribution, and arguably, broad physiological plasticity (Ziegler et al. 2014; Sawall et al. 2015; Voolstra and Ziegler 2020; Ziegler et al. 2019).

At large, however, the proportion of TEs in the genomes of both pocilloporid species were in line with expectation following Kidwell (2002), which predicts 28.29% and 29.71% for P. verrucosa and S. pistillata, respectively. Of note, a direct comparison with P. damicornis was not possible, due to the lack of high-identity repetitive content, which evaded assembly (Cunning et al. 2018).

Transcriptome Assembly, Gene Models, and Protein Annotation

RNA sequencing yielded 121,056,948 read pairs (i.e., 242 mio. paired-end reads) that produced 49,384 coral transcripts available for gene model prediction. Based on BUSCO , we found a high number of complete single-copy orthologs (893 of 978, 91.31%), whereas some orthologs were complete but duplicated (33, 3.37%). Consequently, the number of fragmented (28, 2.86%) and missing (24, 2.45%) single-copy orthologs was rather small, indicating that the reference transcriptome assembly is largely complete.

We identified a final set of 27,439 protein-coding genes with a mean CDS of 1,849.67 bp (table 2). Of these, 96.60% (n = 26,506) have complete ORFs and 89.88% (n = 24,662) are corroborated by RNA-Seq evidence. Despite largely similar gene statistics (table 2), the proportion of intronless gene models in P. verrucosa (10.62%) was comparable to that of P. damicornis (12.26%), but three times higher than in S. pistillata (3.41%). Intronless genes are commonly thought to arise from horizontal gene transfer events of sequences of bacterial origin or through reverse transcription (Baumgarten et al. 2015). As such, the higher number of intronless genes in both Pocillopora genomes may indicate extensive horizontal gene transfer or increased TE element activity, as genes encoding for reverse transcriptases are most commonly found as components of retrotransposons in eukaryotes (Lescot et al. 2016).

Table 2

Genomic Gene Sets of Pocillopora verrucosa (Present Study) in Comparison to the Other Pocilloporid Genomes of Pocillopora damicornis and Stylophora pistillata

	P. verrucosa	P. damicornis	S. pistillata
Genes
Number of protein-encoding genes	27,439	26,077	25,769
Mean gene length (bp)	7,566.63	5,859.80	8,378.20
Genes with annotation (%)	97.87	N/A	88.91
Genes with annotation to SwissProt (%)	70.80	59.70	67.93
Exons
Mean number of exons per gene	7.54	7.17	7.85
Mean exon length (bp)	316.92	244.86	265.66
Mean coding region (CDS) length (bp)	1,849.67	1,365.30	1,833.33
Introns
Percentage of genes with introns (%)	89.38	87.74	96.59
Mean number of introns per gene	7.37	7.03	7.09
Mean intron length (bp)	784.93	665.17	918.25

The majority of the 27,439 protein-coding genes were annotated (supplementary file S2, Supplementary Material online): in total, 97.87% (n = 26,855) had an annotation, the large majority of which (70.80%, n = 19,427) retrieved hits from the SwissProt database, followed by 7,145 proteins (26.04%) with hits in TrEMBL, and 283 proteins (1.03%) that showed similarity to proteins from the NCBI “nr” database. Only a small percentage of proteins ( 2.13%, n = 584) had no hits to either database.

Toward Generation of Model Species Reference Genomes

In this study, we present the first genome assembly for the coral species P. verrucosa. This is an important resource to advance the understanding of coral reef ecology and evolution, particularly because P. verrucosa is a common and globally distributed reef builder. As demonstrated above, we have assembled >90% of the P. verrucosa genome with a contig and scaffold N50 and gene completeness on par with the available genomes of P. damicornis and S. pistillata. We have accomplished this through utilization of an available genome of a closely related species used for scaffolding after initial assembly based on the sequencing of only a single genomic library, with overlapping paired-end reads. We recognize the challenge of assembling a coral genome particularly with regard to their high heterozygosity. We hope that the P. verrucosa genome assembly will be useful for a number of studies ranging from comparative genomics to genetic variants identification (e.g., SNPs and microsatellites).

Supplementary Material

Supplementary data are available at Genome Biology and Evolution online.

Click here for additional data file.