Annotated Draft Genomes of Two Caddisfly Species Plectrocnemia conspersa CURTIS and Hydropsyche tenuis NAVAS (Insecta: Trichoptera).

Members of the speciose insect order Trichoptera (caddisflies) provide important ecosystem services, for example, nutrient cycling through breaking down of organic matter. They are also of industrial interest due to their larval silk secretions. These form the basis for their diverse case-making behavior that allows them to exploit a wide range of ecological niches. Only five genomes of this order have been published thus far, with variable qualities regarding contiguity and completeness. A low-cost sequencing strategy, that is, using a single Oxford Nanopore flow cell per individual along with Illumina sequence reads was successfully used to generate high-quality genomes of two Trichoptera species, Plectrocnemia conspersa and Hydropsyche tenuis. Of the de novo assembly methods compared, assembly of low coverage Nanopore reads (∼18×) and subsequent polishing with long reads followed by Illumina short reads (∼80-170× coverage) yielded the highest genome quality both in terms of contiguity and BUSCO completeness. The presented genomes are the shortest to date and extend our knowledge of genome size across caddisfly families. The genomic region that encodes for light (L)-chain fibroin, a protein component of larval caddisfly silk was identified and compared with existing L-fibroin gene clusters. The new genomic resources presented in this paper are among the highest quality Trichoptera genomes and will increase the knowledge of this important insect order by serving as the basis for phylogenomic and comparative genomic studies.

Introduction

With >16,000 extant species (Morse: Trichoptera World Checklist; available online: http://entweb.sites.clemson.edu/database/trichopt/; last accessed December 4, 2019) the holometabolous insect order Trichoptera (caddisflies) is the seventh most speciose order of all insects (Adler and Foottit 2017). Members of Trichoptera provide important ecosystem services (Macadam and Stockan 2015; Morse et al. 2019). They stabilize gravel bed sediments, play an important role in food webs, and are important biological indicator organisms used for assessing and monitoring water quality (Morse 2017). Similar to lepidopteran caterpillars, caddisfly larvae produce silk in specially modified labial glands. Their silk secretion forms the basis for their diverse case-making behavior which, in turn, may allow caddisfly larvae to exploit a wide range of ecological niches. Caddisfly silk is an important target for interdisciplinary research because it is a tough, energy dissipating fiber that adheres to diverse substrates completely submerged in water (Stewart and Wang 2010; Ashton et al. 2013, 2016), making it an attractive model for tough and adhesive synthetics, including materials with potential applications in medicine (e.g., Brooks 2015; Lane et al. 2015). We report the first draft genomes of two caddisfly species, each generated with a single Oxford Nanopore flow cell along with Illumina sequencing data. This low-cost sequencing strategy was tested on two Trichoptera species from families for which genomic resources were previously not available: Plectrocnemia conspersa CURTIS (Polycentropodidae) and Hydropsyche tenuis NAVAS (Hydropsychidae) of the suborder Annulipalpia. The most closely related species for which a high-quality assembled genome is available is Stenopsyche tienmushanensis Hwang. It shared a common ancestor with Hydropsyche ∼177 Ma and with Plectrocnemia ∼166 Ma, representing large gaps in evolutionary time among groups (Malm et al. 2013; Luo et al. 2018). Despite their ecological importance and potential technological applications, caddisflies are underrepresented among publicly available genomes. Only five genomes have been published thus far, most with low-quality continuity and completeness (table 1). The new genomic resources presented in this paper will serve as a foundation for future genomic and applied research on caddisflies.

Table 1

Comparison of Genome Assemblies among the Seven Published Caddisfly Genomes

Species	Accession	Suborder	Sequencing Platform (Coverage a )	Assembly Length (bp)	Scaffold N50 (kb)	BUSCOs Present b (%)	Estimated Haploid Genome Size (1C)
Hydropsyche tenuis	VTFK00000000	Annulipalpia	Nanopore+Illumina	229,663,394	2190	98.3	222 Mbc
(this study)			(16.5× + 167.6×)a				(FCM: 258 Mbd)
Plectrocnemia conspersa	VTON00000000	Annulipalpia	Nanopore+Illumina	396,695,105	869	98.6	315 Mbc
(this study)			(17.1× + 82.9×)b				(FCM: 455 Mb)
Stenopsyche tienmushanensis	v1	Annulipalpia	PacBio+Illumina	451,494,475	1297	98.7	407c to 453 Mbe
(Luo et al. 2018)			(153× + 150×)
Limnephilus lunatus	Llun_2.0	Integripalpia	Illumina	1369,180,260	69.1	94.3	n.a.
(i5k Consortium)			(80.1×)
Glossosoma conforme	ASM334726v1	Annulipalpia	Illumina	604,293,666	16.7	96.8	1.52 Gbe
(Weigand et al. 2018)			(53×)
Sericostoma sp.	ASM300347v1	Integripalpia	Illumina	1015,727,762	3.1	76.2	616 Mbe
(Weigand et al. 2017)			(43×)
Glyphotaelius pellucidus	n.a.	Integripalpia	Illumina	757,289,448	1.47	62.2	n.a.
(Ferguson et al. 2014)			(8.12×)

Coverage of data used for genome assembly.

N Insecta = 1,658; present = complete + fragmented.

Based on GenomeScope (Vurture et al. 2017).

Based on closely related species.

Based on 17-mer analysis.

FCM, flow cytometry; n.a., not available.

Materials and Methods

Genomic DNA was extracted from larvae (H: tenuis: 49°3′29.31″N, 13°17′44.02″E, P. conspera: 50°15′28.58″N, 9°25′1.23″E). Illumina paired-end sequencing was performed on a HiSeq 2000 sequencer. The Oxford Nanopore library was sequenced on a single flow cell using the MinION portable DNA sequencer. Further details and sequence read processing are described in supplementary note 1, Supplementary Material online.

Genome size estimates were conducted using flow cytometry (FCM) according to Otto (1990) using Lycopersicon esculentum cv. Stupické polnítyčkové rané (2C = 1.96 pg; Doležel et al. 1992) as internal standard and from the short-read sequence data with GenomeScope 2.0 (Ranallo-Benavidez et al. 2019) and backmap.pl v0.1 (Schell et al. 2017; supplementary note 2, Supplementary Material online). GC content was estimated based on comparison of nuclei fluorescence sample/standard ratios stained with the DNA intercalating propidium iodide and AT-specific DAPI using the protocols and the GC content calculation tool by Šmarda et al. (2008). We confirm with our genomic data that this method, which is usually applied in plant research can be successfully applied on animals.

Different de novo genome assembly methods (short-, long-read, and hybrid assembly approaches) were tested with the P. conspersa data set (supplementary note 3, Supplementary Material online). The pipeline, which yielded the best genome for P. conspersa, was used to assemble the H. tenuis genome. Quality of assemblies was evaluated based on continuity [QUAST v5.0.2 (Gurevich et al. 2013)] and BUSCO completeness (Simão et al. 2015; Waterhouse et al. 2018). Genome assemblies were screened for potential contaminations with BlobTools v1.0 (supplementary note 4, Supplementary Material online).

Genome annotation is described in supplementary note 5, Supplementary Material online. In short, repeatModeler v1.0.11 (Smit and Hubley 2008–2015) was used to identify de novo repeat families. Masking and annotating repetitive elements were conducted using RepeatMasker v4.0.7 (Smit et al. 2013–2015). BLAST-like Alignment Tool v3.6 (BLAT, Kent 2002) was used to align the transcriptome of Plectrocnemia sp. (INSqiqTBXRAAPEI-18) to the P. conspersa and the transcriptome of Arctopsyche palpata (INSeqtTAARAAPEI-18) to the H. tenuis draft genome sequence. Transcriptomes were provided by 1KITE (http://www.1kite.org/; last accessed November 11, 2019). Protein-coding genes were annotated with AUGUSTUS v3.3 (Stanke et al. 2008). All protein sequences were assigned putative names by BlastP Protein–Protein BLAST 2.2.30+ searches (Camacho et al. 2009) and were functionally annotated using InterProScan v5 (Zdobnov and Apweiler 2001; Quevillon et al. 2005) and command line Blast2Go v1.3.3 (Götz et al. 2008).

Sequences of L-fibroin were extracted from the functional annotations and aligned with publicly available L-fibroin genes in Geneious v8.0.5. Distance and maximum likelihood trees were reconstructed as described in supplementary note 7, Supplementary Material online.

Results and Discussion

Genome Assembly

These two Trichoptera genomes were sequenced using a low-cost sequencing approach. Even so, the genome assembly statistics are comparable to the previously published S.tienmushanensis genome, which was sequenced at a much higher depth of coverage, and they exhibited improved continuity (based on N50) when compared with all other Trichoptera genome assemblies (table 1). Of the different de novo assembly methods compared for P. conspersa (supplementary table 1, Supplementary Material online), sequence assembly of Nanopore reads (∼17.5× and 16.5× coverage for P. conspersa and H. tenuis, respectively) using wtdbg2 (Ruan and Li 2019), followed by polishing with long Nanopore reads [Racon (Vaser et al. 2017)], and Illumina short reads [Pilon (Walker 2014); 82.9× and 167.6× coverage for P.conspersa and H. tenuis, respectively] performed best with respect to both genome quality and computing time. The draft P. conspersa assembly comprised 1,614 scaffolds (1,613 ≥ 1,000 bp) with a cumulative length of 396,695,105 bp. It exhibited a contig N50 of 868,980, a contig L50 of 141 and a GC content of 30.1%. The genome contained 97.8% (complete: 93.5%, fragmented: 4.3%) of an Endopterygota core gene collection indicating an almost complete coverage of known single copy orthologs in the coding fraction. Remapping the preprocessed Illumina reads showed that 99.5% could be placed within the assembly (supplementary fig. 1, Supplementary Material online). The wtdbg2-Racon-Pilon pipeline also yielded excellent results for the H. tenuis data, which had even lower coverage Nanopore data (∼16.5× coverage). After contaminated contigs were filtered out using BlobTools (supplementary fig. 2B, Supplementary Material online), the assembly consisted of 403 scaffolds of ≥1,000 bp with a cumulative length of 229,663,394 bp. The assembly had a contig N50 of 2,190,134 bp, the largest known contig N50 for a caddisfly, a contig L50 of 33 and a GC content of 33%. In total, 97.7% (complete: 94.2%, fragmented: 3.5%) of the Endopterygota core gene collection was detected by BUSCO. Mapping the preprocessed Illumina reads showed a mapping rate of 99.3% (supplementary fig. 1, Supplementary Material online). FCM-based estimation revealed a GC content of 30.39% for Hydropsyche siltalai and 29.77% for P. conspersa, which shows that GC rich regions are not underrepresented in the genome assemblies and thus adds toward quality of the assemblies.

Genome Size Estimations

Regarding P. conspersa (supplementary note 2, Supplementary Material online; table 1), the mean haploid genome size estimate resulting from the mapping-based genome size estimation with backmap.pl was ∼365 Mb, FCM yielded an estimate of 455.20 ± 8.39 Mb, and GenomeScope predicted a size of ∼316 Mb. A high level of repeats, that fail to assemble correctly may account for this discrepancy. Differences between k-mer-based methods and FCM have been observed before and can be explained by the presence of repetitive elements, which can affect k-mer estimates (Austin et al. 2017; Sánchez-Herrero et al. 2019). Generally, our results indicate that researchers should not rely on a single method for estimating genome size. For H. tenuis, the genome size estimates based on sequencing data ranged from ∼222 to ∼229 Mb (supplementary note 2, Supplementary Material online). This was in line with FCM-based estimates of closely related species of Hydropsyche (H. siltalai: 257.72 ± 7.24 Mb, H. pellucidula 283.94 ± 4.16 Mb, H. saxonica 241.65 ± 4.41 Mb), and with the assembly length. To our knowledge, we present the first genome size estimates for Trichoptera based on FCM. Varying genome sizes (based on sequence data) have been observed between suborders, that is, between retreat making Annulipalpia [Hydropsyche: 241.65–283.94 Mb; Stenopsyche: 407–453 Mb (Luo et al. 2018), Plectrocnemia: ∼ 455 Mb, Glossosoma conforme: 1.5 Gb (Weigand et al. 2018)] and case-building Integripalpia [Sericostoma sp.: 616 Mb (Weigand et al. 2017), assembly length of Limnephillus lunatus: ∼1.4 Gb (i5k Consortium), table 1], but also among families within suborders (table 1). Further studies are needed to investigate genome size evolution in Trichoptera.

Genome Annotation and Functional Annotation of Protein-Coding Genes

A total of 34.36% of the P. conspersa genome assembly was masked as repeats. A 32.39% of the annotated repeats were interspersed repeats. Almost half of the interspersed repeats, 17.58%, could not be classified by comparison with known repeat databases, and therefore may be specific for Trichoptera. Only 19.12% of the H. tenuis genome assembly was masked as repeats. A 17.16% of the repeats were classified as interspersed repeats. A large portion, 10.89%, of the interspersed repeats remained unclassified. In comparison, repeat content in the previously published S.tienmushanensis is 36.76 (10.56% unclassified; Luo et al. 2018). Details on the repeat classes are given in supplementary note 6, Supplementary Material online. The annotation of the genomes resulted in the prediction of 19,104 proteins and 18,586 genes and 15,844 proteins and 15,873 genes for P. conspersa and H. tenuis, respectively. In total, BUSCO detected 95.6% of the Endopterygota core gene collection in the predicted proteins of P. conspersa (complete: 88.1%, fragmented: 7.5%) and H. tenuis (complete: 90.5%, fragmented: 5.1%), respectively. In comparison, 15,658 proteins and 14,672 genes were annotated in S. tienmushanensis (Luo et al. 2018). More details can be found in supplementary table 2, Supplementary Material online. Of the annotated proteins 70% and 62% showed significant sequence similarity to entries in the NCBI nr database for P. conspersa and H. tenuis, respectively. When searched against InterPro databases, 14,476 (P. conspersa: 75.7%) and 12,949 (H. tenuis: 81.7%) annotated proteins include motifs/domains were identified by InterProScan of which 7,192 (P. conspersa) and 6,889 (H. tenuis) genes were assigned to Gene Ontology IDs with a corresponding InterPro entry (supplementary figs. 3 and 4, Supplementary Material online). The major biological process found in all three caddisfly genomes was a cellular process (P. conspersa: 1,639 genes, H. tenuis: 1,460 genes, S. tienmushanensis: 1,523 genes). Although binding was the largest subcategory in molecular function in H. tenuis (1,338 genes) and S. tienmushanensis 1,506 genes), catalytic activity was the largest subcategory in P. conspersa (1,415 genes). Regarding the cellular component category, most genes were assigned to the cell subcategory in H. tenuis (1,017 genes) and S. tienmushanensis (1,008 genes) and to the membrane subcategory in P. conspersa (1,123 genes). InterProScan repeat distribution revealed that leucine-rich repeats ranked highest among all repeats for P. conspersa and H. tenuis, whereas ankyrin-rich repeats ranked highest in S. tienmushanensis. A detailed distribution of InterProScan repeats can be found in supplementary figure 5, Supplementary Material online.

Identification of the L-fibroin Gene and Evolutionary Implications

L-fibroin, a protein component of larval caddisfly silk was successfully assembled, annotated, and compared with previously published L-fibroin sequences. Translation of the L-fibroin sequence yielded a protein of 280 and 276 amino acid residues for Hydropsyche tenuis and Plectrocnemia conspersa, respectively. Amino acid sequences of the L-fibroin proteins exhibit 56.5% pairwise identity and 26.7% identical sites among the eight caddisfly species (fig. 1). The L-fibroins of Trichoptera exhibit a similar organization to the L-fibroin of Lepidoptera (19.2% pairwise identity; 18.7% identical sites; fig. 1: Bombyx mori) suggesting that the ancestral L-fibroin gene evolved prior to the separation of the two orders and was retained in both (Yonemura et al. 2009). Moreover, distribution of residues with characteristic properties is similar in both orders, for example, that of Cystein (polar, hydrophobic, and neutral) which has proved to play an important role for L-fibroin function in B. mori (Tanaka et al. 1999). Within Trichoptera, the conservation of regions with alternating groups of hydrophobic and hydrophilic residues can be observed. The two reconstructed trees of the L-fibroin exon sequences agree closely with previously reported phylogenetic trees (Kjer et al. 2016; Frandsen et al., in preparation). The distance tree is separated into two large clades that correspond to the suborders Annulipalpia (fig. 1I) and Integripalpia (fig. 1II). Within Annulipalpia, Hydropsyche (Hydropsychidae) and Stenopsyche (Stenopsychidae) form a clade with Plectrocnemia (Polycentropodidaea) as sister. The maximum likelihood tree revealed a similar topology (fig. 1) and high bootstrap support (>90) for most of the nodes. However, the position of Plectrocnemia remains unclear. Although Hydropsychidae are typically characterized by a relative flat filter net built adjacent to their actual retreat (Clarke 1883; Schuhmacher 1970), larvae of Plectrocnemia are “sit-and-wait” predators with silken retreats attached to aquatic plants with many irregularly radiating threads forming capture nets adapted for predation of organisms in the benthos which become entangled in its nets (Wesenberg-Lund 1911; Townsend and Hildrew 1976, Lillpopp et al. 1998). Retreats/capture nets of Stenopsychidae are loosely woven (Morse et al. 2019). The capture-net portion consists of three different parts which serve as a cover, feeding mesh, and collecting seston (Tanida 2000). Current phylogenetic evidence suggests that retreat-making behavior followed evolutionary transitions from 1) organisms with a silken retreat that housed the organism with a separated filter net outside of the larval retreat to 2) larvae with retreats that also served as filter nets to 3) more advanced retreats (such as trumpet nets) that enable behaviors such as collecting fine particulate matter, gardening, scraping perphyton, and predation (Morse et al. 2019).

. 1.

—(A) Jukes–Cantor distance tree. I: Annulipalpia; II: Integripalpia. (B) Maximum likelihood tree. Bootstrap values ≥90% are given for each node. L-fibroins of the two sequenced genomes are indicated in red. Genbank accessions: Bombyx mori: NM_001044023.1 (Suetsugu et al. 2013), Rhyacophila obliterata: AB354690.1 (Yonemura et al. 2009), Limnephilus decipienns: AB214510.1 (Yonemura et al. 2006), Hesperophylax occidentalis: KM384738.1 (Wang et al. 2014), Hydropsyche tenuis: this study, Hydropsyche angustipennis: AB214508.1 (Yonemura et al. 2006), Plectrocnemia conspersa: this study, Stenopsyche marmorata: LC057252.1 (Bai et al. 2015), and Stenopsyche tienmushanensis: Luo et al. (2018). (C) Aligned amino acid residues of L-fibroin. Each color in the alignment represents a different amino acid. Mean pairwise identity over all pairs in the column: green, 100% identity; greeny-brown, at least 30% and under 100% identity; red, below 30% identity.

The new genomic resources presented in this paper, which includes the full assembly of the L-fibroin are a further step to disentangling strands of silk, that is, to better understand the evolution of silk biosynthesis, silk production, and associated behavior in caddisflies. Further studies will include this data to assess the phylogeny of caddisflies and evolution of silk usage, and to investigating structural differences in genes involved in biosynthesis of silk and silk usage regulation.

Supplementary Material

Supplementary data are available at Genome Biology and Evolution online.

Click here for additional data file.