The chromosome-level genome assembly of the slug Deroceras laeve facilitates its use as a comparative model of regeneration

Abstract

The genome of the land slug Deroceras laeve was sequenced, assembled up to the chromosome level, and annotated for non-coding RNAs and protein-coding genes. Due to the small size of this pulmonate species, ease of laboratory culture, cosmopolitan distribution, as well as recently released anatomical and histological resources, this genomic resource creates new opportunities for the investigation of the largely unexplored mechanisms of regeneration in mollusks. Moreover, it also makes this slug an attractive model for functional genomics and evolutionary biology.

Keywords: Genome Assembly; Hi-C; PacBio; Stylommatophora; cognate chromosomes; mollusk.

Fig. 1.

Hi-C-guided chromosome level scaffolding. a) Heatmap of the normalized proximity signal from the Hi-C analysis. Blue line squares enclose the scaffolds that are defined by an enriched proximity signal along the diagonal. Scaffolds are arranged in descending order by length. The Hi-C coverage is shown along the y and x axes. b) Hi-C coverage markedly differs between the main assembly and the unplaced scaffolds. c) Whole genome self-alignment of the D. laeve assembly showing long stretches of inter-chromosomal synteny as off-diagonal alignments. The vertical and horizontal blue lines indicate the boundary between D. laeve chromosomes and unplaced scaffolds. d) Conserved one-to-one synteny between the 31 chromosomes of 2 Deroceras gastropods.

Fig. 2.

Map of Deroceras laeve mitochondrial genome. a) Circular map of the mitochondrial genome with features annotated. The sense of arrows indicate the direction of transcription. b) Alignment of NC_072953 to Hi-C scaffold 1563. c) Linear map of the mitochondrial chromosome. Gaps in the segments in the shaded box represent mutations and small notches, or indels. d) Phylogenetic reconstruction with full mitochondrial sequences from the limacids D. lasithionense, D. reticulatum, Ambigolimax valentianus, a recently reported D. laeve sequence (Dl Chin) and the one obtained in this study (Dl Mex).

Fig. 3.

The repeat content landscape occupies 70% of the slug genome. a) General repeat divergence landscape showing the evolutionary distribution of the most common kind of repeats and sub-classifications. b) Repeat divergence landscape of DNA elements. c) Divergence landscape of LINEs. d) Divergence landscape of LTRs, RC/Helitron and SINEs. e–i) Distribution of selected types of repeats in Chromosomes 1, 4, 10, 23, and 30 and j) in Unplaced scaffolds sorted from longest to shortest. For each chromosome, panels show Top: TTAGG repeats found in a given tandem element. Middle: Span of tandem repeats colored by GC enrichment (lighter) or depletion (darker) compared to the genome average GC% (42%). Red triangles show the location of piRNA clusters. Bottom: stacked densities of different classes of interspersed repeats in 1 Mb windows: DNA elements, Rolling Circle/Helitron, LINEs, LTR, SINEs and Unknown classification.

Fig. 4.

Metrics of D. laeve genome annotation quality. a) Benchmarking Universal Single-Copy Orthologs (BUSCO) for the most general (eumetazoan) and the most specific (mollusca) datasets available for the organism. EVM, genes obtained with EvidenceModeller; EVM post PASA, EVM genes enriched with PASA assemblies. b) OMArk results for completeness (percent out of 2,373 lophotrochozoa conserved Hierarchical Orthologous Groups (HOGs) found in D. laeve proteome; left) and consistency (percentage of 24,337 D. laeve proteins that correspond to a gene family known to exist in the selected lineage; right) c). Comparison of OMArk percentages of D. laeve to other gastropod proteomes present in the OMArk database.

Fig. 5.

Functional overview of the D. laeve proteome and evolution within mollusks via gene duplication and gene loss. a) Proportions and number of transcripts assigned to 5 categories according to their KEGG annotation: Organismal Systems, Metabolism, Genetic Information Processing, Cellular Processes, and Environmental Information Processing. b) Central circle area: Syntenic relationships between chromosomes revealed by close paralogy connections of genes annotated by eggNOG Mapper. Stronger interchromosomal connections highlight notable correspondences. Rings from inside-out: Intrachromosomal arches connect paralogous genes that are more than 5 Mb apart. Smaller arches connect paralogous genes located less than 5 Mb apart (most are so close to each other that only a tick can be observed). Solid notches in the following outer ring mark clusters of 5 or more related genes within segments of 5 Mb. Solid bars along the circumference represent chromosomes and tickmarks indicate intervals of 5 Mb. The outer-most ring indicates the abundance of telomeric motifs as in Fig. 3 and the chromosome number. The eggNOG functional descriptions of some genes are provided in the outer region along with the location of the Hox gene clusters in chromosomes 6 and 25. c) Frequency of the top GO terms associated with the expanded gene families found in a CAFE analysis of mollusk proteomes. d) Reconciliation tree depicting inferred protein gains and losses in extant mollusk species and their ancestors. Numbers in parentheses at speciation events indicate de novo generated proteins while those at the leaf level represent proteins with no identified orthologs.

Fig. 6.

Analysis of RNA-silencing pathways in Deroceras laeve. a) Small RNAs from different organs. Plots show frequencies of 18–32 nucleotide RNAs sequenced from whole juvenile slugs, as well as head, foot, and ovotestes from adult animals. Each plot represents an average of 3 biological replicates. b) Dicer-dependent small RNA biogenesis relies on a single Dicer homolog with canonical domain organization. Below the domain organization, the motifs important for ATP-dependent helicase activity implicated in processive long dsRNA cleavage are highlighted in red (Kidwell et al. 2014). Species name initials: Av, Arion vulgaris; Bg, Biomphalaria glabrata; Ce, Caenorhabditis elegans; Dl, D. laeve; Dm, Drosophila melanogaster; Hs, Homo sapiens; Mm, Mus musculus; Sp, Schizosaccharomyces pombe. c) Analysis of AGO-clade Argonaute proteins. Two AGO protein-coding genes were found in the D. laeve genome assembly. Both proteins have the typical AGO protein architecture and their closest homolog in other species is AGO2. Hence, they were denoted AGO2a and AGO2b. Below the domain composition, the conservation of the 4 amino acid residues in the PIWI domain that form the “catalytic tetrad” in the functional nuclease is highlighted by shaded boxes (Nakanishi et al. 2012).

Fig. 7.

Analysis of the piRNA pathway in Deroceras laeve. a) Two PIWI protein-coding genes were found in the D. laeve genome assembly. Both proteins have the typical PIWI protein architecture and are considerably divergent (BLAST alignment shows 41% identity and 61% similarity). Below the domain composition, conservation of the 4 amino acid residues in the PIWI domain that form the “catalytic tetrad” is highlighted in shaded boxes. PIWIa has all residues conserved, PIWIb carries one E>D aminoacid substitution. b) A representative genome track depicting a typical unidirectional piRNA cluster (piC10) with 27–30 nucleotide small RNAs mapped to it. The 27–30 nucleotide reads mapping to the $+$ (sense) and $-$ (antisense) strands are shown separately. CPM, counts per million. c) Ranked annotated piRNA clusters displaying the cumulative fraction (CF) of piRNAs. The top 21 piRNA clusters ranked by reads per million (RPM) represent 90% of piRNAs mapping to all clusters. d) A sequence logo of all 27–30 nucleotide reads mapped to the genome (top) and reads mapping to annotated piCs (bottom). A typical 1U pattern is observed for the entire 27–30 nt RNA population and the signal is increased when only considering reads mapping to piCs. The 10A signature typical of the ping-pong mechanism is not apparent. e) Analysis of 27–30 nucleotide RNA overlaps. The most common overlap is for 10 nucleotides suggesting that a fraction of piRNAs originates from piRNA-mediated cleavage, although this fraction is not large enough to manifest as a 10A signature in the sequence logo.

Fig. 8.

Hox clusters of D. laeve. a) Structure of Hox clusters in Drosophila and Vertebrates. Drosophila has an Antenappedia and a Bithorax complex. Vertebrates possess in general 4 clusters with Paralogous Groups (PG) 1 to 13, equivalent to the Drosophila genes. PG2 and PG3 are orthologous to pb and PG9-PG13 are orthologous to AbdB. The inferred ancestral Hox cluster had 6 genes. Both Drosophila and vertebrate clusters possess a conserved miR-10. Adapted from (Mallo and Alonso 2013). b) D. laeve HoxA cluster in chromosome 6 has putative orthologs of the first 3 ancestral PGs and a miR-10-like locus in the expected place in comparison to Drosophila and vertebrates. The second part of this cluster, which is located after approximately 3 Mb on the same chromosome, contains possible derivatives of PG7-8 and another copy of miR-10 upstream of these genes. c) D. laeve HoxB cluster in chromosome 25 seems to contain only bithorax “posterior” genes from PG6-13. Another copy of miR-10-like miRNA is present within this cluster.