Print
Choice of DNA extraction method affects DNA metabarcoding of unsorted invertebrate bulk samples
expand article infoMarkus Majaneva§, Ola H. Diserud|, Shannon H.C. Eagle§, Mehrdad Hajibabaei§, Torbjørn Ekrem
‡ Norwegian University of Science and Technology, Trondheim, Norway
§ University of Guelph, Guelph, Canada
| Norwegian Institute for Nature Research, Trondheim, Norway
Open Access

Abstract

Characterisation of freshwater benthic biodiversity using DNA metabarcoding may allow more cost-effective environmental assessments than the current morphological-based assessment methods. DNA metabarcoding methods where sorting or pre-sorting of samples are avoided altogether are especially interesting, since the time between sampling and taxonomic identification is reduced. Due to the presence of non-target material like plants and sediments in crude samples, DNA extraction protocols become important for maximising DNA recovery and sample replicability. We sampled freshwater invertebrates from six river and lake sites and extracted DNA from homogenised bulk samples in quadruplicate subsamples, using a published method and two commercially available kits: HotSHOT approach, Qiagen DNeasy Blood & Tissue Kit and Qiagen DNeasy PowerPlant Pro Kit. The performance of the selected extraction methods was evaluated by measuring DNA yield and applying DNA metabarcoding to see if the choice of DNA extraction method affects DNA yield and metazoan diversity results. The PowerPlant Kit extractions resulted in the highest DNA yield and a strong significant correlation between sample weight and DNA yield, while the DNA yields of the Blood & Tissue Kit and HotSHOT method did not correlate with the sample weights. Metazoan diversity measures were more repeatable in samples extracted with the PowerPlant Kit compared to those extracted with the HotSHOT method or the Blood & Tissue Kit. Subsampling using Blood & Tissue Kit and HotSHOT extraction failed to describe the same community in the lake samples. Our study exemplifies that the choice of DNA extraction protocol influences the DNA yield as well as the subsequent community analysis. Based on our results, low specimen abundance samples will likely provide more stable results if specimens are sorted prior to DNA extraction and DNA metabarcoding, but the repeatability of the DNA extraction and DNA metabarcoding results was close to ideal in high specimen abundance samples.

Key Words

Macroinvertebrates, Biomonitoring, Freshwater

Introduction

Modern biomonitoring and water quality assessments of freshwater ecosystems rely on standardised sampling and identification of benthic macroinvertebrates (Bonada et al. 2006). Samples are usually collected using kick-nets or surber-nets in lotic environments, while different ways of grab- or core sampling typically are used in deeper lentic environments (European Community 2000). Sufficient sampling effort is important to obtain an appropriate representation of the fauna (Bongard et al. 2011; references therein). In most cases, the effort needed to sample a large proportion of the community is substantial, resulting in massive volumes of samples. For an undisturbed boreal highland river in Norway, for instance, it is estimated that > 3000 specimens of Ephemeroptera, Plecoptera and Trichoptera (EPT) must be sampled in order to cover 95% of the EPT diversity (Bongard et al. 2011). This corresponds to eight-ten minutes of active kick-sampling. Therefore, the following sorting and morphological identification of specimens is very time-consuming since samples are very large and especially if non-EPT taxa are included. Subsampling is often implemented to expedite the sample handling and identifications rely on taxonomic expertise and knowledge (e.g. availability of taxonomic keys for a specific group and location); therefore, the same samples can produce different taxonomic lists depending on the identifier and the subsample examined (Haase et al. 2010). In this situation, approaches based on DNA metabarcoding promise more efficient and standardised sample processing (Hajibabaei et al. 2011, Yu et al. 2012), especially since large or replicate samples from the same location will not require significant additional effort.

In DNA metabarcoding, short homologous DNA fragments – e.g. part of the cytochrome c oxidase I (COI) gene – from bulk samples are sequenced in parallel and the resulting reads are matched against a reference library to identify the taxa present in the samples (Hajibabaei et al. 2011, Taberlet et al. 2012). The approach has proven to give more detailed taxon lists than morphological identification (e.g. Ji et al. 2013, Gibson et al. 2015, Emilson et al. 2017) and several studies have demonstrated the suitability of the method in monitoring of freshwater macroinvertebrates (Hajibabaei et al. 2011, Carew et al. 2013, Elbrecht et al. 2017b). Although some methodological limitations, such as estimation of biomass (Elbrecht and Leese 2015), detection of rare species (Leray and Knowlton 2017) and incomplete reference libraries (Hajibabaei et al. 2016; Hering et al. 2018) still exist, DNA metabarcoding based monitoring of freshwater macroinvertebrates has attracted great interest and effort has been made towards large-scale implementation (Baird and Hajibabaei 2012, Hajibabaei et al. 2016) i.e. in accordance with the EU Water Framework Directive (Leese et al. 2016; Leese et al. 2018).

Efficient sample homogenisation and DNA extraction is a prerequisite for any successful DNA-based study. Conventionally, DNA has been extracted using desalting and organic solvents like phenol-chloroform (Tan and Yiap 2009). Commercial solid-phase (spin-column based) methods make use of the finding that DNA binds to a silica membrane in the presence of a high concentration of chaotropic salts (Boom et al. 1990) and non-target substances are rinsed off while the target DNA is bound to the silica membrane. Although more expensive, the spin-column methods have advantages compared to the earlier ‘artisanal’ and often toxic DNA extraction methods, since the former are standardised, easy to perform and can save time. Therefore, many freshwater DNA metabarcoding studies are performed using commercial spin column kits for DNA extraction (e.g. Carew et al. 2013; Vivien et al. 2016; Blackman et al. 2017; Emilson et al. 2017), although some research groups have used salt extraction (Elbrecht et al. 2017a; Elbrecht et al. 2017b). However, DNA may be made readily available for downstream analyses by simply using alkaline lysis and neutralisation in a single tube (Truett et al. 2000). Applying such a simple DNA extraction method to DNA metabarcoding samples would further facilitate the reduction of sample processing time.

Extensive scientific literature on comparisons of DNA extraction methods exists, each piece of literature concentrating on different target taxa and environments and showing differences in DNA yield and PCR success amongst methods (e.g. Miller et al. 1999; Fredricks et al. 2005; Mahmoudi et al. 2011; Psifidi et al. 2015; Mäki et al. 2017; Schiebelhut et al. 2017). High-quality DNA was best achieved using bead-based homogenisation (Miller et al. 1999; Mahmoudi et al. 2011) and best PCR success was achieved with efficient inhibitor removal (Miller et al. 1999; Whitehouse and Hottel 2007; Mahmoudi et al. 2011; Psifidi et al. 2015). However, performance may vary by taxon, for instance in the case of fungi and macroinvertebrates (Fredricks et al. 2005; Schiebelhut et al. 2017). Although no universal best method for DNA extraction exists, optimal methods might be available in literature or be developed based on existing approaches for specific applications.

Most freshwater DNA metabarcoding studies rely on separating biomass of organisms from the non-target organic and inorganic material and on grinding the organisms before adding lysis buffer (e.g. Gibson et al. 2015; Elbrecht et al. 2017b; Emilson et al. 2017). Separating target biomass reduces sample size before grinding and efficiently removes the main contributors to PCR inhibition, the phenolic compounds (Wilson 1997). In addition, there are some benefits of sorting by specimen size before DNA extraction for metabarcoding (Elbrecht et al. 2017a). However, there is still the issue of time used for separating the target biomass of samples, especially for sufficiently large samples (Bongard et al. 2011), including all macroinvertebrate groups.

In various soil studies, DNA has been successfully extracted directly, without first separating organisms (Ogram et al. 1987; Van Elsas et al. 1997; Delmont et al. 2011; Taberlet et al. 2012). The use of such unsorted bulk samples for DNA extraction of macroinvertebrates is appealing due to the dramatic reduction in sample processing time (i.e. no specimen sorting), but only a few freshwater macroinvertebrate studies have used such techniques (Hajibabaei et al. 2011, Gardham et al. 2014, Dowle et al. 2016, Andújar et al. 2018). The sample homogenisation and DNA extraction steps become crucial when using unsorted bulk samples because uneven specimen distribution in the sample and presence of non-target organic matter and derived PCR inhibitors might influence the DNA metabarcoding success.

Here, we compared the success of DNA metabarcoding for assessing a set of unsorted freshwater benthic macroinvertebrate samples, processed using different DNA extraction methods. Rather than being an exhaustive comparison, we chose to compare methods that have low toxicity, are easy to use and require little handling time. Thus, no conventional time-consuming or toxic desalting and organic solvent-based methods were included.

Material and methods

Samples were collected from three sites in Norway: River Atna at Skranglehaugen (10 August 2015, 61.98186°N, 09.80454°E, 1117 m above sea level) and at Vollen (10 August 2015, 61.98471°N, 10.02823°E, 720 m above sea level), as well as Lake Jonsvatn at Jonsborg (28 September 2015, 63.39569°N, 10.55370°E, 150 m above sea level). River Atna originates in the Rondane National Park in Central Norway and is a well-documented Nordic freshwater ecosystem (Aagaard et al. 2004). Lake Jonsvatn, near Trondheim, is a moderately large lake (area 15 km2, maximum depth 97 m) and is the main source for drinking water for the city of Trondheim.

Four minute benthic kick samples were collected from the River Atna (Bongard et al. 2011) and 0.01-m2 Van Veen grab samples from the benthos of the Lake Jonsvatn. The Van Veen samples were collected from 0.25, 2, 7.5 and 15 m depth, from the shoreline to approximately 50 m from the shore. Samples are called Skranglehaugen, Vollen, 15 m, 7.5 m, 2 m and 0.25 m hereafter (Fig. 1). Excess plant material was carefully removed from the samples on site and 96% ethanol was added. The next day, the ethanol was changed to fresh 96% ethanol and the samples were stored at +4 °C. After storing the samples three to four months (until 16 December 2015), excess ethanol was carefully collected from the samples and the remaining material including all plant material, animals, sediment and ethanol was homogenised, using a hand-held blender. Between each sample, the blender was sterilised by submerging it first in deionised water, then in a 5% bleach solution for 15 min. After 15 min, the bleach was washed away from the blender by submerging it in deionised water then in double-deionised water. Finally, the blender was dried under UV light for 30 min. From each of these homogenised samples, excess ethanol was removed and twelve 27.6–130.5-mg (mean 69.4 mg, standard deviation 25.8 mg; Fig. 2) subsamples per sample were collected into microcentrifuge tubes, using forceps sterilised in the same manner as the blender. Samples were shaken vigorously in between subsampling.

Figure 1.

Experimental setup. The 4-minute kick samples were collected from two sites along the River Atna and the Van Veen grab samples were collected from a depth gradient in the Lake Jonsvatn. From each sample, 12 subsamples were taken and DNA was extracted using three extraction methods (four subsamples per extraction method); HotSHOT extraction, Qiagen DNeasy Blood & Tissue Kit and MO BIO PowerPlant Pro Kit. An extraction blank was done for each extraction method.

Figure 2.

Subsample sizes (a) and DNA yields (b). Two-way ANOVA followed by Tukey’s HSD was used to test differences amongst the extraction methods and F-statistic value (F) and significance (p) are given. The small letters denote significantly different extraction methods at each sampling site based on site and extraction interaction factors.

The subsamples were randomly divided into three groups exposed to three different DNA extraction protocols (Fig. 1). Prior to extraction using the DNeasy Blood and Tissue Kit (Qiagen GmbH, Hilden, Germany), the subsamples were homogenised, using Tissue Lyser (Qiagen) with 19.0/5 frequency for 2 min and incubated in 20 µl of proteinase K overnight at 56 °C on a rocking platform (200 rpm). The PowerPlant Pro DNA Isolation Kit (MO BIO Laboratories, Carlsbad, CA, USA) was used according to manufacturer’s protocol (currently named Qiagen DNeasy PowerPlant Pro Kit). The third DNA extraction protocol was based on the HotSHOT approach (Truett et al. 2000), which is a cheap and quick DNA extraction method based on lysis of the tissues in an alkaline sodium hydroxide and neutralisation with an acidic Tris buffer. The extracted DNA was quantified with Qubit 2.0 (Invitrogen, CA, USA), using the dsDNA HS Assay Kit according to the manufacturer’s protocol.

We amplified three fragments of the mitochondrial COI gene. The full barcode region (approximately 648 bp) was amplified, using the standard Folmer et al. (1994) LCO1490 (GGTCAACAAATCATAAAGATATTGG) and HC02198 (TAAACTTCAGGGTGACCAAAAAATCA) primers. The approximately 239bp long F230 fragment at the 5’ end of the barcode region was amplified, using the standard Folmer LCO1490 primer and the 230_R (Gibson et al. 2015) reverse primer (CTTATRTTRTTTATICGIGGRAAIGC). The approximately 314bp long BE fragment at the 3’ end of the barcode region was amplified, using the forward primer B CCIGAYATRGCITTYCCICG (Hajibabaei et al. 2011) and the reverse primer R5 GTRATIGCICCIGCIARIACIGG (Gibson et al. 2014). The first-round PCR amplifications were undertaken once with attached Illumina forward (TCGTCGGCAGCGTCAGATGTGTATAAGAGACAG) and reverse (GTCTCGTGGGCTCGGAGATGTGTATAAGAGACAG) adapters.

The first-round PCRs (Suppl. material 1: Table 1) had a final volume of 25μl containing 2 μl DNA template, 17.8 μl molecular biology grade water, 2.5 μl 10× reaction buffer (200 mM Tris HCl, 500 mM KCl, pH 8.4), 1 μl MgCl2 (50 mM), 0.5 μl dNTPs mix (10 mM), 0.5 μl forward primer (10 mM), 0.5 μl reverse primer (10 mM) and 0.2 μl Invitrogen’s Platinum Taq polymerase (5 U/μl; Thermo Fisher Scientific, Waltham, MA, USA). Slightly modified PCRs were necessary for 15 m (all HotSHOT subsamples, F230 and BE Blood & Tissue Kit subsamples and one F230 PowerPlant Kit subsample), 7.5 m (all HotSHOT) and 2 m (all HotSHOT) samples. These reactions contained increased concentrations of forward primer (1 μl), reverse primer (1 μl) and Taq (0.5 μl). All PCRs included negative control reactions (no DNA template). The PCR conditions for the F230 and BE fragments were, with a heated lid, 94 °C for 5 min, followed by a total of 35 cycles of 94 °C for 40 s, 46 °C for 1 min and 72 °C for 20 s and a final extension at 72 °C for 2 min and hold at 10 °C. The modified PCRs were amplified, using 40 cycles. The PCR conditions for the full barcode region were, with heated lid, at 94 °C for 1 min, followed by 5 cycles of 94 °C for 40 s, 45 °C for 40 s and 72 °C for 1 min, 35 cycles of 94 °C for 40 s, 51 °C for 40 s and 72 °C for 1 min and a final extension at 72 °C for 5 min and hold at 10 °C. The DNA template was 1:50 diluted due to PCR inhibition (except some subsamples; see Suppl. material 1: Table 1). Samples that were amplified with increased conditions were purified, using Qiagen’s MinElute PCR purification Kit according to manufacturer’s instructions. PCR products were visualised on a 1.5% agarose gel to check the amplification success and a subset of samples was quantified using PicoGreen (Thermo Fisher Scientific) according to the manufacturer’s protocol. In the second-round PCRs, the Illumina tailed amplicons were dual indexed, using Nextera XT Indices (FC-131-1002, Illumina, San Diego, CA, USA) in a reduced-cycle PCR according to the manufacturer’s protocol. The amplification reactions contained the same reagent concentrations as above with three modifications: 3 µl of amplicon, 1 µl of each primer and 0.25 µl of Platinum Taq (Invitrogen). Indexed amplicons were pooled into a library and sequenced on a standard flow cell with three extraction controls, using the 600-cycle V3 Illumina MiSeq sequencing kit (MS-102-3003).

The resulting raw amplicon reads (available in the ENA SRA repository with study name PRJEB26589) were processed, using mothur v.1.36.1 (Schloss et al. 2009, see Suppl. material 5: Methods for details). The forward and reverse reads of the F230 and BE fragments were merged keeping only the reads with no mismatch in primer sequence, using the command make.contigs. The forward and reverse fragments overlapped on average with 282 bases (F230 reads) and 238 bases (BE reads). All merged F230 reads shorter than 250 bases and longer than 340 bases and all merged BE reads shorter than 330 and longer than 375 bases, as well as merged reads with ambiguous bases were removed, using the command screen.seqs. Primer sequences were removed, using the command trim.seqs. The full barcode region reads were processed separately since the forward and reverse reads could not be merged reliably (only few bases overlapping). All reads with mismatches in primer sequences, with ambiguous bases, with homopolymers longer than 12 bases and shorter than 200 bases were removed; primer sequences were removed from the retained reads and the reads were truncated when the average quality score decreased to under 30 in 30 bases sliding window, using the command trim.seqs. After the quality control, the resulting good quality reads were processed further, using usearch v8.1.1831_win32 (Edgar 2013). Exact duplicates were removed, using the command -derep_fulllength and the reads were de-noised, using a sequence identity threshold of 98% in the command -cluster_otus. The above steps were undertaken for each subsample and the resulting operational taxonomic unit (OTU) fasta files were pooled to F230, BE, Folmer forward and Folmer reverse files, using merge.files in mothur. Then, the resulting pooled OTU files were de-replicated, using -derep_fulllength and the final OTUs (each fragment separately) were clustered, using -cluster_otus at 97% similarity level. The abundance of each OTU in each sample was searched, using -usearch_global against the pooled OTU-files. All OTUs that had less than 6 reads/sample were removed. Also, OTUs present only or predominantly in three negative samples were removed (42 F230 and BE OTUs, 38 full barcode forward OTUs and 12 full barcode reverse OTUs).

The OTUs were assigned taxonomically in two steps. First, they were searched against the NCBI non-redundant nucleotide database, using the BLAST 2.3.0+ (Zhang et al. 2000) (F230 and BE fragments searched 17 November 2016 and forward and reverse fragments of the full barcode region 7 April 2017). Taxonomic assignment of an OTU was undertaken, using the lowest common ancestor algorithm in MEGAN 6.5.10 (minimum bit score 100, top percentage 8.0 and minimum support 1; Huson et al. 2016). The Metazoa OTUs were translated to amino acids and aligned, using Mafft online tool (Katoh and Standley 2013) to remove all non-COI OTUs and OTUs with stop codons in the mitochondrial invertebrate code. Secondly, the Metazoa OTUs were assigned taxonomically using the BOLD v.4 identification engine and Species Level Barcode Records (Ratnasingham and Hebert 2007; F230 and BE fragments 24 November 2016 and forward and reverse fragments of the full barcode region 21 April 2017; > 97% matches to a species in BOLD were considered as an assignment). The F230, BE and the forward and reverse full barcode fragment OTU tables were combined, but only OTUs that matched with an invertebrate species in the BOLD v.4 Species Level Barcode Records were kept (889 OTUs). All OTUs assigned to the same taxonomic name were merged leaving 305 invertebrate species. The different DNA extraction methods were compared, using both non-normalised and normalised species level reads (36991 reads/sample).

To find statistically significant differences (p < 0.05) in DNA yield and in metazoan diversity amongst sites, one-way Kruskal-Wallis followed by pair-wise comparisons was used. To find significant differences in subsample size, in DNA yield and in metazoan diversity amongst different DNA extraction methods, two-way ANOVA followed by Tukey’s pair-wise comparisons was used. Partitioned pair-wise betadiversity values (βsim and βsne, Baselga 2010) were calculated using the R-package betapart (Baselga et al. 2018, R Core Team 2016) to see how much of the variation in species richness was due to nestedness (βsne, subsamples with less species include a subset of the species in subsamples with more species) and turnover (βsim, change in species composition). To study repeatability of the DNA extraction methods, similarity in invertebrate community composition was analysed by fitting the bivariate Poisson-lognormal species abundance distribution to pairs of subsamples (Engen et al. 2002, Engen et al. 2011, Diserud et al. 2013), using the R-package poilog (Grøtan and Engen 2015, R Core Team 2016). Species abundances from a community were assumed to follow a lognormal distribution and the sampling process to be described by Poisson sampling. For a species with a given abundance, the expected number of reads observed will then equal the abundance scaled by the sampling intensity. When considering two community subsamples simultaneously, the log abundances in the pair of communities were assumed to follow the binormal distribution, with the estimated correlation from this bivariate Poisson-lognormal species abundance distribution serving as a robust measure of similarity between subsamples. This approach utilises all the species abundance information, accounts for the sampling process and can deal with over-dispersion relative to the Poisson, that is, the variance in the observed log abundances is greater than the mean (Engen et al. 2002). When community subsamples are grouped according to the DNA extraction method, the dissimilarity between pairs of subsamples, after Poisson sampling has been accounted for, can be attributed to over-dispersion and the extraction method. Assuming the variance component caused by different sampling sites to be the same for all extraction methods, the different methods can be evaluated by how consistently they represent the community structure, as quantified by the bivariate correlation. If the dissimilarity between subsamples is large for a given extraction method, this method’s representation of the community structure will vary more than for an extraction method with lower dissimilarity.

Results

Total DNA was extracted from 28–131 mg subsamples, four subsamples/extraction method from each site (Figs 12). The weights of subsamples were similar overall (Fig. 2a; two-way ANOVA, extraction factor F2,5 = 0.55, site factor F2,5 = 1.04, p > 0.05) in spite of a crossover interaction (two-way ANOVA, interaction factor F2,5 = 3.60, p < 0.01), which was due to HotSHOT extraction subsamples being significantly heavier in the 2 m samples and significantly lighter in the 0.25 m samples compared to the PowerPlant Kit subsamples (Tukey’s HSD on extraction at each site, p < 0.05). The DNA yield correlated weakly with the subsample weights overall (Spearman’s ρ = 0.31, p < 0.01). This correlation could be attributed to a significant and strong correlation of the PowerPlant Kit subsample weight and DNA yield (Spearman’s ρ = 0.77, p < 0.001). Weights of Blood & Tissue Kit and HotSHOT extraction subsamples did not correlate with the DNA yield (Spearman’s ρ, p > 0.05). All HotSHOT extraction, 75% of Blood & Tissue Kit and 25% of PowerPlant Kit subsamples suffered from PCR inhibition as only diluted templates worked in PCR (Suppl. material 1: Table 1). There was no correlation between the amount of initial material extracted and inhibition.

The DNA yield was lower in 15 m, 7.5 m and 2 m samples (42–495 ng) than in 0.25 m, Skranglehaugen and Vollen samples (169–4280 ng; Kruskal-Wallis χ2=50.99, p < 0.001, pair-wise Mann-Whitney with sequentially Bonferroni corrected p < 0.05; Fig. 2b). The DNA yield of PowerPlant Kit was significantly higher than that of HotSHOT extraction in the 15 m and 0.25 m samples and significantly higher than that of Blood & Tissue Kit in the 15 m, 7.5 m and 0.25 m samples; also the DNA yield of HotSHOT extraction was significantly higher than that of Blood & Tissue Kit in the 15 m and Skranglehaugen samples (two-way ANOVA, interaction factor F2,5 = 3.99, Tukey’s HSD on interaction, p < 0.001). Overall, the highest DNA yield was achieved with the PowerPlant Kit, also when the subsample weight was taken into consideration (Suppl. material 4: Figure 1a).

A total of 10.3, 9.8 and 1.3 million reads were generated from the F230, BE and full barcode fragment, respectively (Suppl. material 3: Table 3). After merging and quality filtering, 3.9 and 3.0 million good-quality paired-end F230 and BE reads as well as 0.5 million good-quality full barcode fragments remained, which were clustered into 2741, 4832 and 2289 OTUs at 97% similarity level, respectively. From those OTUs, 1426, 1261 and 1430 OTUs were affiliated with Metazoa in an NCBI BLAST search and from those, 616, 489 and 588 OTUs were taxonomically assigned to a metazoan species in BOLD (> 97% matches to a species in BOLD). After merging the taxonomic information from all three fragments, 305 OTUs were taxonomically assigned to an invertebrate taxon (named DNA species hereafter), which were used for comparing the extraction methods. The majority of DNA species were Insecta (78%). Clitellata (8%), Maxillopoda (4%) and Branchiopoda (2%) were the next richest groups.

The number of DNA species was the highest using PowerPlant Kit in the 15 m and 2 m samples, using Blood & Tissue Kit in the 7.5 m and Vollen samples and using HotSHOT extraction in the Skranglehaugen samples (Fig. 3; two-way ANOVA, interaction factor F2,5 = 5.09, p < 0.001, Tukey’s HSD on extraction at each site p < 0.05) and, thus, the number of DNA species did not give a clear indication of which extraction method performed best overall.

Figure 3.

Number of DNA species in each subsample. Two-way ANOVA followed by Tukey’s HSD was used to test differences amongst the extraction methods and F-statistic value (F) and significance (p) are given. The small letters denote significantly different extraction methods at each sampling site based on site and extraction interaction factors.

The portion of shared DNA species amongst the extraction methods was lower for the 15 m, 7.5 m and 2 m samples (Fig. 4a–c) than for the 0.25 m Skranglehaugen and Vollen samples (Fig. 4g–i). The beta diversity (measured as Sørensen dissimilarity) of the Blood & Tissue Kit was significantly higher than that of the HotShot extraction in the 15 m, 7.5 m and 2 m samples and significantly higher than that of the PowerPlant Kit in the 15 m, 7.5 m, 2 m and 0.25 m samples (Fig. 4d–f, j–l; mean beta diversity value ± 95% confidence interval higher). The beta diversity of the HotShot extraction samples was also significantly higher than that of the PowerPlant Kit samples in the 7.5 m, 2 m and 0.25 m samples. The significantly lower beta diversity values of the PowerPlant Kit samples than those of Blood & Tissue Kit and HotShot samples indicate less change in community composition using PowerPlant Kit in the lake site. The beta diversity values were similar amongst the extraction methods in the river site (Fig. 4k–l). However, partitioning of Sørensen dissimilarity showed that the change in community composition due to nestedness of PowerPlant Kit samples was significantly higher than that of the Blood & Tissue Kit and HotShot extraction in Vollen sample (Fig. 4l; mean βsne value ± 95% confidence interval higher). This indicates that the PowerPlant Kit subsamples were more subsets of each other and less affected by species replacement than the Blood & Tissue Kit and HotShot extraction subsamples. Overall, the PowerPlant Kit performed best in the lake site and showed indications of better repeatability in the river site.

Figure 4.

The number of shared DNA species amongst extraction methods illustrated using Venn diagrams (a–c; g–i) and pairwise beta diversity values within extraction methods at each site measured as Sørensen dissimilarity (βsor, which is partitioned to nestedness, βsne and species turnover, βsim) (d–f; j–l), lines denote group means and 95% confidence intervals for βsor and βsne. The small letters denote significantly different extraction methods at each sampling site based on non-overlapping 95% confidence intervals for βsor and βsne.

When evaluating repeatability of extraction methods in more detail, based on the bivariate Poisson-lognormal correlations, the PowerPlant Kit outperformed both Blood & Tissue Kit and HotSHOT extraction in most sites (Fig. 5; mean bivariate correlation value ± 95% confidence interval higher than that with other methods at each site except 15 m where Blood & Tissue Kit and PowerPlant Kit confidence intervals overlap). The mean bivariate correlation value was higher than 0.80 in all except the 15 m samples using PowerPlant Kit, while using Blood & Tissue Kit and HotSHOT extraction, it was higher than 0.80 only in river samples. This indicates that the Blood & Tissue Kit and HotSHOT extraction subsamples failed to give a consistent community composition in the lake samples.

Figure 5.

Bivariate correlation values, using the full dataset. Circles denote values for 6 pairs of samples within each extraction method, lines denote group means and 95% confidence intervals and all 66 bivariate correlation values within the sampling site are summarised in the box plots. The small letters denote significantly different extraction methods at each sampling site based on non-overlapping 95% confidence intervals.

Discussion

In this study, we compared three DNA extraction methods on unsorted bulk samples of lentic and lotic macroinvertebrates. DNA metabarcoding, based on unsorted samples, holds great potential for large scale monitoring of freshwaters because of reduced sample processing time, especially for samples with a high invertebrate to debris ratio. Our results show that the DNA extraction method employed influences the success of DNA metabarcoding for assessing macroinvertebrate biodiversity and is dependent on environmental factors such as the benthic substrate, type of vegetation and specimen abundance.

Our sampling sites differed markedly in their physical and biological composition. The Vollen and Skranglehaugen sites are characterised by abundant boreal and alpine species living in variable lotic microhabitats (Aagaard et al. 2004), which were thoroughly sampled with our four-minute kick-samples. The kick-samples resulted in two to three litres of sampled material, including both macroinvertebrates and some living plant material. Our sampling of the lake Jonsvatn depth gradient was more cursory as four 0.01-m2 Van Veen grab samples in ca. 50 m line represented each depth. Thus, different microhabitats were not as extensively sampled as in the river. In addition, the grab samples included dead plant material and excess sediments as well as some living plant material. We took parallel samples at each site for sorting and morphological identification (Majaneva et al. in prep.) and, as an example, a lake grab sample at 15 m included only 17 specimens compared to over 10000 specimens in the river kick-samples. The lentic samples included fewer macroinvertebrates (Figs 3, 4) and less DNA (Fig. 2) and they had lower consistency of the subsamples (Fig. 5) than the lotic samples, leading to uncertainty of the metabarcoding-based community composition results for the lentic site. Consequently, more stable results may be achieved if such low-specimen lentic site samples were sorted prior to DNA extraction and DNA metabarcoding, as is the usual practice currently undertaken in DNA metabarcoding studies.

All DNA extraction methods tested might have been affected by PCR inhibition because of the presence of non-target organic matter in our bulk samples. However, the lotic samples suffered less than the lentic samples from inhibition, probably because the plant material present in the lotic samples was mainly bryophytes that have a relatively low amount of phenolic compounds (Maksimova et al. 2013), while the lentic samples included remains of vascular plants. There are several ways to deal with inhibition after DNA extraction, but previous studies have shown that species-specific detectability is reduced by 25% using column purification and 52% using dilution (McKee et al. 2015, Buxton et al. 2017). Thus, it is preferable to maximise inhibitor removal already in the DNA extraction step. In our test, the PowerPlant Kit, which includes a specific inhibitor removal step, performed best as it suffered least from PCR inhibition and yielded most DNA in strong correlation with the subsample weight (Fig. 2). The high yield observed in our results is in contrast to previous studies where kits that included inhibitor removal steps had lower yields than kits without this step (Mahmoudi et al. 2011; Eichmiller et al. 2016). Obviously, it is possible that the PowerPlant Kit performs better for extracting plant DNA present in the samples and this DNA can contribute to the higher yield.

Some samples were particularly difficult in terms of PCR inhibitors. In the 15 m lake sample, where also the PowerPlant Kit failed to produce DNA extracts unaffected by PCR inhibition, each DNA extraction method resulted in some close to zero and negative bivariate correlation values (Fig. 3d). Correlation close to zero indicates that the communities are independent in structure, while a negative correlation show that dominant species in one community tend to be rare in the other and vice versa (Diserud et al. 2013). An obvious reason for the low consistency in the lentic subsamples is the lower relative and total number of target specimens compared to the lotic kick samples and subsequent inferior target specimen homogenisation of the lentic samples. In such situations, the time spent on sorting small samples is beneficial compared to using unsorted samples for DNA metabarcoding. Another factor that could have improved the consistency in our lentic subsamples is using DNA extraction developed for a greater amount of starting material, e.g. > 1 g instead of 100 mg (Ranjard et al. 2003; Zinger et al. 2016), although some of such extractions methods have been found unreliable for PCR (Mahmoudi et al. 2011). In any case, using homogenisation of bulk samples when abundance and diversity is low is questionable as it is easy to envisage a situation where the ecological status is incorrectly inferred due to chance. This may also apply for other low abundance and diversity samples, for example if using environmental DNA samples to monitor invertebrate biodiversity.

Often only presence/absence data are used for community studies based on DNA metabarcoding due to the low reliability of abundance data (Elbrecht and Leese 2015). Based on our results, no DNA extraction method was best in terms of number of species (Fig. 3). PowerPlant Kit subsamples showed significantly better repeatability both using incidence based sample nestedness (Fig. 4) and abundance based bivariate correlation values (Fig. 5), but as the overall dissimilarity amongst subsamples was high and proportion of shared species amongst the extraction methods was low in the low-abundance lentic samples (Fig. 4), we do not recommend homogenising low-abundance bulk samples prior to sorting.

Based on our results, all tested DNA extraction methods appear suitable for large, high-abundance unsorted samples, particularly from lotic systems. Smaller samples with fewer specimens are less time-consuming to sort and likely provide more stable results if PCR-inhibiting substances are physically removed prior to DNA extraction, although DNA extraction with the inhibitor removal step performed significantly better in the case of low-abundance samples. In conclusion, we urge for continued testing of homogenised unsorted bulk samples and comparison with results from presorted samples since time spent on sorting samples might be more efficiently used for monitoring more localities using DNA metabarcoding.

Acknowledgements

Thanks to Elisabeth Stur and Karstein Hårsaker for field assistance and Rachel Smith for laboratory assistance. The paper is part of the project ‘Environmental Barcoding of Aquatic Invertebrates (EBAI)’ funded by the Research Council of Norway (project no. 243791/E50) and the Norwegian Environment Agency (project no. 15040013). This publication contributes to the EU COST Action DNAqua-Net (EU COST Action CA15219).

References

  • Andújar C, Arribas P, Gray C, Bruce K, Woodward G, Yu DW, Vogler AP (2018) Metabarcoding of freshwater invertebrates to detect the effects of a pesticide spill. Molecular Ecology 27: 146–166. https://doi.org/10.1111/mec.14410.
  • Blackman RC, Constable D, Hahn C, Sheard AM, Durkota J, Hänfling B, Handley LL (2017) Detection of a new non-native freshwater species by DNA metabarcoding of environmental samples – first record of Gammarus fossarum in the UK. Aquatic Invasions 12: 177–189. https://doi.org/10.3391/ai.2017.12.2.06
  • Bongard T, Diserud OH, Sandlund OT, Aagaard K (2011) Detecting invertebrate species change in running waters: An approach based on the sufficient sample size principle. The Open Environmental & Biological Monitoring Journal 4: 72–82. https://doi.org/10.2174/1875040001104010072
  • Boom R, Sol CJ, Salimans MM, Jansen CL, Wertheim-van Dillen PM, van der Noordaa J (1990) Rapid and simple method for purification of nucleic acids. Journal of Clinical Microbiology 28: 495–503.
  • Carew ME, Pettigrove VJ, Metzeling L, Hoffmann AA (2013) Environmental monitoring using next generation sequencing: rapid identification of macroinvertebrate bioindicator species. Frontiers in Zoology 10: 45. https://doi.org/10.1186/1742-9994-10-45
  • Delmont TO, Robe P, Clark I, Simonet P, Vogel TM (2011) Metagenomic comparison of direct and indirect soil DNA extraction approaches. Journal of Microbiological Methods 86: 397–400. https://doi.org/10.1016/j.mimet.2011.06.013
  • Diserud OH, Stur E, Aagaard K (2013) How reliable are Malaise traps for biomonitoring? – A bivariate species abundance model evaluation using alpine Chironomidae (Diptera). Insect Conservation and Diversity 6: 561–571. https://doi.org/10.1111/icad.12012
  • Dowle EJ, Pochon X, Banks JC, Shearer K, Wood SA (2016) Targeted gene enrichment and high-throughput sequencing for environmental biomonitoring: a case study using freshwater macroinvertebrates Molecular Ecology Resources 16: 1240–1254. https://doi.org/10.1111/1755-0998.12488.
  • Eichmiller JJ, Miller LM, Sorensen PW (2016) Optimizing techniques to capture and extract environmental DNA for detection and quantification of fish. Molecular Ecology Resources 16: 56–68. https://doi.org/10.1111/1755-0998.12421
  • Elbrecht V, Leese F (2015) Can DNA-based ecosystem assessments quantify species abundance? testing primer bias and biomass – sequence relationships with an innovative metabarcoding protocol. PloS one10: e0130324. https://doi.org/10.1371/journal.pone.0130324
  • Elbrecht V, Peinert B, Leese F (2017a) Sorting things out: Assessing effects of unequal specimen biomass on DNA metabarcoding. Ecology and Evolution 7: 6918–6926. https://doi.org/10.1002/ece3.3192
  • Elbrecht V, Vamos E, Meissner K, Aroviita J, Leese F (2017b) Assessing strengths and weaknesses of DNA metabarcoding based macroinvertebrate identification for routine stream monitoring. Methods in Ecology and Evolution 8: 1265–1275. https://doi.org/10.1111/2041-210X.12789
  • Emilson CE, Thompson DG, Venier LA, Porter TM, Swystun T, Chartrand D, Capell S, Hajibabaei M (2017) DNA metabarcoding and morphological macroinvertebrate metrics reveal the same changes in boreal watersheds across an environmental gradient. Scientific Reports 7: 12777. https://doi.org/10.1038/s41598-017-13157-x
  • Engen S, Lande R, Walla T, DeVries PJ (2002) Analyzing spatial structure of communities using the two-dimensional Poisson lognormal species abundance model. The American Naturalist 160: 60–73. https://doi.org/10.1086/340612
  • European Community (2000) Directive 2000/60/EC of the European parliament and of the council of 23 October 2000 establishing a framework for community action in the field of water policy. Official Journal of the European Communities L 327: 1–72.
  • Folmer O, Black M, Hoeh W, Lutz R, Vrijenhoek R (1994) DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Molecular Marine Biology and Biotechnology 3: 294–299.
  • Gardham S, Hose GC, Stephenson S, Chariton AA (2014) DNA metabarcoding meets experimental ecotoxicology: Advancing knowledge on the ecological effects of copper in freshwater ecosystems. Advances in Ecological Research 51: 79–104. https://doi.org/10.1016/B978-0-08-099970-8.00007-5
  • Gibson J, Shokralla S, Porter TM, King I, van Konynenburg S, Janzen DH, Hallwachs W, Hajibabaei M (2014) Simultaneous assessment of the macrobiome and microbiome in a bulk sample of tropical arthropods through DNA metasystematics. Proceedings of the National Academy of Sciences of the United States of America 111: 8007–8012. https://doi.org/10.1073/pnas.1406468111
  • Gibson JF, Shokralla S, Curry C, Baird DJ, Monk WA, King I, Hajibabaei M (2015) Large-scale biomonitoring of remote and threatened ecosystems via high-throughput sequencing. PLoS ONE 10: e0138432. https://doi.org/10.1371/journal.pone.0138432
  • Haase P, Pauls SU, Schindehütte K, Sundermann A (2010) First audit of macroinvertebrate samples from an EU Water Framework Directive monitoring program: human error greatly lowers precision of assessment results. Journal of the North American Benthological Society 29: 1279–1291. https://doi.org/10.1899/09-183.1
  • Hajibabaei M, Shokralla S, Zhou X, Singer GAC, Baird DJ (2011) Environmental barcoding: a next-generation sequencing approach for biomonitoring applications using river benthos. PLoS ONE 6: e17497. https://doi.org/10.1371/journal.pone.0017497
  • Hajibabaei M, Baird DJ, Fahner NA, Beiko R, Golding GB (2016) A new way to contemplate Darwin’s tangled bank: how DNA barcodes are reconnecting biodiversity science and biomonitoring. Philosophical Transactions of the Royal Society B: Biological Sciences 371: 20150330. https://doi.org/10.1098/rstb.2015.0330
  • Hering D, Borja A, Jones I, Pont D, Boets P, Bouchez A, Bruce K, Drakare S, Hänfling B, Kahlert M, Leese F, Meissner K, Mergen P, Reyjol Y, Segurado P, Vogler A, Kelly M (2018) Implementation options for DNA-based identification into ecological status assessment under the European Water Framework Directive. Water Research 138: 192–205. https://doi.org/10.1016/j.watres.2018.03.003
  • Huson DH, Beier S, Flade I, Górska A, El-Hadidi M, Mitra S, Ruscheweyh H-J, Tappu R (2016) MEGAN Community Edition - Interactive exploration and analysis of large-scale microbiome sequencing data. PLoS Computational Biology 12: e1004957. https://doi.org/10.1371/journal.pcbi.1004957
  • Ji Y, Ashton L, Pedley SM, Edwards DP, Tang Y, Nakamura A, Kitching R, Dolman PM, Woodcock P, Edwards FA, Larsen TH, Hsu WW, Benedick S, Hamer KC, Wilcove DS, Bruce C, Wang X, Levi T, Lott M, Emerson BC, Yu DW (2013) Reliable, verifiable and efficient monitoring of biodiversity via metabarcoding. Ecology Letters 16: 1245–1257. https://doi.org/10.1111/ele.12162
  • Katoh K, Standley DM (2013) MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Molecular Biology and Evolution 30: 772–780. https://doi.org/10.1093/molbev/mst010
  • Leese F, Altermatt F, Bouchez A, Ekrem T, Hering D et al. (2016) DNAqua-Net: Developing new genetic tools for bioassessment and monitoring of aquatic ecosystems in Europe. Research Ideas and Outcomes 2: e11321. https://doi.org/10.3897/rio.2.e11321
  • Leese F, Bouchez A, Abarenkov K, Altermatt F, Borja A et al. (2018) Why we need sustainable networks bridging countries, disciplines, cultures and generations for aquatic biomonitoring 2.0: A perspective derived from the DNAqua-Net COST Action. Advances in Ecological Research 58: 63–99. https://doi.org/10.1016/bs.aecr.2018.01.001
  • Leray M, Knowlton N (2017) Random sampling causes the low reproducibility of rare eukaryotic OTUs in Illumina COI metabarcoding. PeerJ 5: e3006. https://doi.org/10.7717/peerj.3006
  • Mahmoudi N, Slater GF, Fulthorpe RR (2011) Comparison of commercial DNA extraction kits for isolation and purification of bacterial and eukaryotic DNA from PAH-contaminated soils. Canadian Journal of Microbiology 57: 623–628. https://doi.org/10.1139/w11-049
  • Mäki A, Salmi P, Mikkonen A, Kremp A, Tiirola M (2017) Sample preservation, DNA or RNA extraction and data analysis for high-throughput phytoplankton community sequencing. Frontiers in Microbiology 8: 1848. https://doi.org/10.3389/fmicb.2017.01848
  • Maksimova V, Klavina L, Bikovens O, Zicmanis A, Purmalis O (2013) Structural characterization and chemical classification of some bryophytes found in Latvia. Chemistry & Biodiversity 10: 1284–1294. https://doi.org/10.1002/cbdv.201300014
  • McKee AM, Spear SF, Pierson TW (2015) The effect of dilution and the use of a post-extraction nucleic acid purification column on the accuracy, precision, and inhibition of environmental DNA samples. Biological Conservation 183: 70–76. https://doi.org/10.1016/j.biocon.2014.11.031
  • Miller DN, Bryant JE, Madsen EL, Ghiorse WC (1999) Evaluation and optimization of DNA extraction and purification procedures for soil and sediment samples. Applied and Environmental Microbiology 65: 4715–4724.
  • Psifidi A, Dovas CI, Bramis G, Lazou T, Russel CL, Arsenos G, Banos G (2015) Comparison of eleven methods for genomic DNA extraction suitable for large-scale whole-genome genotyping and long-term DNA banking using blood samples. PLoS ONE 10: e0115960. https://doi.org/10.1371/journal.pone.0115960
  • R Core Team (2016) R: A language and environment for statistical computing. R foundation for Statistical Computing, Vienna. https://www.R-project.org/
  • Ranjard L, Lejon DPH, Mougel C, Schehrer L, Merdinoglu D, Chaussod R (2003) Sampling strategy in molecular microbial ecology: influence of soil sample size on DNA fingerprinting analysis of fungal and bacterial communities. Environmental Microbiology 5: 1111–1120.
  • Schiebelhut LM, Abboud SS, Gómez Daglio LE, Swift HF, Dawson MN (2017) A comparison of DNA extraction methods for high-throughput DNA analyses. Molecular Ecology Resources 17: 721–729. https://doi.org/10.1111/1755-0998.12620
  • Schloss PD, Westcott SL, Ryabin T, Hall JR, Hartmann M, Hollister EB, Lesniewski RA, Oakley BB, Parks DH, Robinson CJ, Sahl JW, Stres B, Thallinger GG, Van Horn DJ, Weber CF (2009) Introducing mothur: Open-source, platform-independent, community-supported software for describing and comparing microbial communities. Applied and Environmental Microbiology 75: 7537–7541. https://doi.org/10.1128/AEM.01541-09
  • Taberlet P, Prud’Homme SM, Campione E, Roy J, Miquel C, Shehzad W, Gielly L, Rioux D, Choler P, Clement J-C, Melodelima C, Pompanon F, Coissac E (2012) Soil sampling and isolation of extracellular DNA from large amount of starting material suitable for metabarcoding studies. Molecular Ecology 21: 1816–1820. https://doi.org/10.1111/j.1365-294X.2011.05317.x
  • Tan SC, Yiap BC (2009) DNA, RNA, and protein extraction: the past and the present. Journal of Biomedicine and Biotechnology 2009: 574398. https://doi.org/10.1155/2009/574398
  • Truett GE, Heeger P, Mynatt RL, Truett AA, Walker JA, Warman ML (2000) Preparation of PCR-quality mouse genomic DNA with hot sodium hydroxide and Tris (HotSHOT). BioTechniques 29: 52–54.
  • Van Elsas JD, Mäntynen V, Wolters AC (1997) Soil DNA extraction and assessment of the fate of Mycobacterium cholorophenolicum strain PC-1 in different soils by 16S ribosomal gene sequence based most probable number PCR and immunofluorescence. Biology and Fertility of Soils 24: 188–195. https://doi.org/10.1007/s003740050230
  • Whitehouse CA, Hottel HE (2007) Comparison of five commercial DNA extraction kits for the recovery of Francisella tularensis DNA from spiked soil samples. Molecular and Cellular Probes 21: 92–96. https://doi.org/10.1016/j.mcp.2006.08.003
  • Wilson IG (1997) Inhibition and facilitation of nucleic acid amplification. Applied and Environmental Microbiology 63: 3741–3751.
  • Yu DW, Ji Y, Emerson BC, Wang X, Ye C, Yang C, Ding Z (2012) Biodiversity soup: metabarcoding of arthropods for rapid biodiversity assessment and biomonitoring. Methods in Ecology and Evolution 3: 613–623. https://doi.org/10.1111/j.2041-210X.2012.00198.x
  • Zinger L, Chave J, Coissac E, Iribar A, Louisanna E, Manzi S, Schilling V, Schimann H, Sommeria-Klein G, Taberlet P (2016) Extracellular DNA extraction is a fast, cheap and reliable alternative for multi-taxa surveys based on soil DNA. Soil Biology and Biochemistry 96: 16–19. https://doi.org/10.1016/j.soilbio.2016.01.008