Methodology article | Open | Published:
Improving the performance of true single molecule sequencing for ancient DNA
BMC Genomicsvolume 13, Article number: 177 (2012)
Second-generation sequencing technologies have revolutionized our ability to recover genetic information from the past, allowing the characterization of the first complete genomes from past individuals and extinct species. Recently, third generation Helicos sequencing platforms, which perform true Single-Molecule DNA Sequencing (tSMS), have shown great potential for sequencing DNA molecules from Pleistocene fossils. Here, we aim at improving even further the performance of tSMS for ancient DNA by testing two novel tSMS template preparation methods for Pleistocene bone fossils, namely oligonucleotide spiking and treatment with DNA phosphatase.
We found that a significantly larger fraction of the horse genome could be covered following oligonucleotide spiking however not reproducibly and at the cost of extra post-sequencing filtering procedures and skewed %GC content. In contrast, we showed that treating ancient DNA extracts with DNA phosphatase improved the amount of endogenous sequence information recovered per sequencing channel by up to 3.3-fold, while still providing molecular signatures of endogenous ancient DNA damage, including cytosine deamination and fragmentation by depurination. Additionally, we confirmed the existence of molecular preservation niches in large bone crystals from which DNA could be preferentially extracted.
We propose DNA phosphatase treatment as a mechanism to increase sequence coverage of ancient genomes when using Helicos tSMS as a sequencing platform. Together with mild denaturation temperatures that favor access to endogenous ancient templates over modern DNA contaminants, this simple preparation procedure can improve overall Helicos tSMS performance when damaged DNA templates are targeted.
Ancient DNA has been shown to survive in fossil material for over thousands to even hundreds of thousands of years [1–3]. However, post-mortem DNA damage reactions, which fragment the DNA backbone into short pieces and generate hydrolytic and oxidative base derivatives, often limit the amount of DNA templates preserved [4–6]. Together with DNA-DNA crosslinks [7, 8] and by-products of the Maillard reaction , these modifications could limit the success of PCR amplification and thus the ability to incorporate ancient DNA fragments into DNA libraries, such as those required for sequencing on second-generation sequencing platforms . As a result of such DNA degradation, the overwhelming majority of DNA within ancient DNA extracts often consists of exogenous contaminant DNA templates, mainly originating from environmental microbes , with endogenous DNA levels rarely greater than a few percent of the sequenceable templates [12, 13].
Notable exceptions have been reported. For example, ancient keratinous materials such as hair, horn and nail, are relatively exempt from exogenous contamination [14, 15], as is the occasional bone that has been exceptionally well preserved in both permafrost  and temperate environments .
Several different targeted re-sequencing techniques have been demonstrated to be able to preferentially access ancient DNA templates over exogenous contaminants, but these have been limited to the characterization of complete mitochondrial genomes [18, 19], complete plasmids  and complete exomes . The exact efficiency of such methods is also unclear with reported enrichment values ranging over several orders of magnitude and subject to whether or not clonal collapse is considered in the calculation . So far, no DNA capture method has been successfully used for recovering gigabases of sequence information and complete genomes. As a result, the characterization of ancient genomes is still principally dependent on shotgun sequencing [17, 23–25]. The introduction of methods that would improve the balance of endogenous over exogenous sequences obtained would therefore be especially beneficial in the emerging field of paleogenomics, since such methods would serve to improve sequencing depth and consequently the quality of the genome drafts generated [17, 23, 25]. Additionally, since most DNA extraction approaches are destructive, improving current methods would help to minimize damage to valuable samples, a point of critical importance in cases where the fossil record is limited (e.g. for Denisovans, a newly discovered group of archaic hominins, currently known by no more than three fossils).
Recently, Green and collaborators have shown that enzymatic restriction digestion can be used to successfully deplete the environmental microbe component within Illumina sequencing libraries, thus enriching the fraction of sequenced material of endogenous origin (in their case, Neandertal) . Similarly, we have recently demonstrated that true Single Molecule DNA Sequencing (tSMS), as performed on the Helicos HeliScope platform, advantageously complemented Illumina GAIIx platforms when sequencing DNA extracted from a Pleistocene horse bone preserved in permafrost . Indeed, a simple modification to the conventional tSMS DNA template preparation procedure, involving a modified denaturation temperature prior to template 3’ poly-A tailing, increased the efficiency of endogenous sequence recovery, suggesting that the efficiency of tSMS could be further improved for sequencing ancient DNA extracts . Inspired by this observation, in this study, we further explore the potential of the tSMS procedure when sequencing ancient DNA extracts through two new modifications, namely oligonucleotide spiking and DNA phosphatase treatment. We show that both approaches increase tSMS performance, with up to 3.3-fold enrichment in the amount of endogenous sequence information recovered per sequencing channel. In addition, we confirm the existence of molecular preservation niches in bones [10, 26], from which DNA can be preferentially extracted, leading to enriched fractions of endogenous DNA. Such niches should be preferentially extracted for characterizing the complete genomes of individuals and species from the past.
Results and discussion
Two horse bone fossils preserved in permafrost were tSMS-sequenced on the Helicos platform using from 110 to 1100 FOV (fields of view) per channel. Overall sequence yield was expected to be directly proportional to the number of FOV. The first bone, CA, originated from Thistle Creek, Yukon Territory, Canada, and was associated with infinite radiocarbon dates (>50,300BP). The second bone, TP, was excavated in the Taymir Peninsula and was radiocarbon dated at 13,389 ± 52 BP. Each bone was extracted in duplicate, henceforth referred to as CA1 and CA2 for sample CA, and TP1 and TP2 for sample TP. In addition, for sample TP, the residual bone powder that remained undigested after the initial DNA extraction was subjected to redigestion until completely digested. DNA was subsequently purified from this second fraction, generating two additional extracts, referred to below as TP1RE and TP2RE.
DNA denaturation temperatures
Ancient DNA molecules are short, damaged, and showing overhangs at both termini. In contrast to the denaturation temperature at 95°C used in standard template preparation for the Helicos platform, mild denaturation temperatures limit further fragmentation of ancient DNA molecules . For all extracts considered, when compared to a high denaturation temperature of 95°C, a more mild denaturation temperature of 80°C was found to yield: (1) a higher fraction of endogenous sequences (Tables 1 and 2); (2) longer reads (Figure 1, top); (3) lower overall %GC contents (Figure 1, middle), and; (4) higher levels of post-mortem cytosine deamination (Figure 1, bottom). Such enrichment could not be explained by the difference in volumes of CA extracts sequenced following denaturation at the respective temperatures since this would be expected to yield differences in the overall total number of sequences recovered but not in the relative proportions of the different fractions of reads present in the extracts. These observations are in agreement with previous observations derived from a more limited number of bone extracts . This suggests that the presence of deaminated cytosine residues reduces the number of hydrogen bonds between DNA strands and thereby likely favor a preferential denaturation of endogenous ancient templates, especially compared to longer, hence more thermostable, exogenous templates. Interestingly, the younger sample (TP) exhibited lower levels of cytosine deamination than the older sample (CA; Additional file 1, Figure S1).
The cytosine deamination levels observed in sequencing reads recovered from the extracts generated by redigesting residual pellet that was left-over from the first extraction (TP1RE and TP2RE), were found to be significantly lower than those observed in the sequencing reads generated from the initial extraction (TP1 and TP2) (Figure 2, bottom). In addition, re-extracted samples showed higher endogenous sequence yield regardless of the DNA denaturation temperature used during the template preparation procedure (Table 2). This confirmed our previous observations with other bone material . We hypothesize that the initially undigested bone powder likely consisted of larger bone crystals, hence required additional time before full digestion could be completed. Crystal growth, and the related diagenetic increase in bone crystallinity observed post-mortem, could contribute to the entrapment and protection of endogenous DNA molecules in larger crystals, as observed with the TP1RE and TP2RE extracts. Conversely, an increase in bone porosity, and related fluid circulation, could result in microbial attack and the dissolution of crystals, which could potentiate the deposition of exogenous DNA templates on the fraction of bone powder that appeared easily decalcified (here, TP1 and TP2 extracts) [27, 28]. We anticipate that the preferential extraction of larger bone crystals would offer a simple strategy for acquiring a higher proportion of endogenous DNA, thus improving the performance of shotgun sequencing.
We additionally explored the potential of two other sequencing template preparation approaches, namely, oligonucleotide spiking and DNA phosphatase treatment, to improve the performance of tSMS on ancient DNA extracts. Such improvement could be either direct by increasing the relative ratio of endogenous content, or indirect, through higher total sequence recovery rates per channel.
In comparison to fresh material, degraded bone samples yield relatively low DNA concentrations that are up to several orders of magnitude below those observed in fresh bone extracts [29, 30]. While samples with extremely low DNA yields are expected to perform poorly on any sequencing platform, such samples are still sequenceable on the Helicos platform when they are spiked with oligonucleotides. In this system, a camera scans the sequencing flow cell and records images of each field of view. After each nucleotide addition, the images must be aligned with each other in order to properly assign sequence information. If the DNA is at a very low concentration the images cannot be aligned properly, but the addition of spikes increases the total yield of DNA sequences, allowing more accurate image alignment so that it becomes possible to obtain sequences from extremely small numbers of endogenous DNA molecules . This method has been shown to work effectively with modern DNA samples present at very low concentrations.
In order to test this idea on ancient DNA extracts, we focused on extracts CA1 and CA2. These extracts showed relatively low levels of endogenous DNA content (ca. 2.4% – 2.9% with 95°C as the denaturation temperature) and delivered low numbers of sequences per channel on a 1100 FOV platform (Table 1). For each extract, DNA templates prepared using 80°C and 95°C as denaturation temperatures were pooled, and spiked just prior to sequencing with a 40 pM solution consisting of 30 oligonucleotides just prior to sequencing (see Methods). The reads that derived from oligonucleotides were identified and filtered by mapping against those 30 known sequences. The results were then compared to the sequencing performance observed when no oligonucleotide spiking was performed. The fractions of extracts that were spiked with a mixture of exogenous oligonucleotides exhibited lower endogenous sequence ratios than unspiked extracts (1.11% – 1.27% versus 3.79% – 7.61%; Table 1). If the spiking procedure was equivalent to the standard procedure then a larger number of endogenous sequences would have been expected (see Methods). In particular, 142,076 – 143,299 and 179,788 – 181,159 horse reads would have been expected following spiking for extracts CA1 and CA2, respectively, in contrast to the 57,946 and 54,266 actually observed (Table 1). This difference shows that oligonucleotide spiking reduced our capability to access endogenous ancient reads. BLAST analyses revealed that the subset of sequence reads that was filtered for the sequence of oligonucleotide spikes presented higher proportions of unidentified sequences (97.6% and 97.8% for samples CA1 and CA2, respectively) compared to unspiked controls (88.1% and 93.3%). Likewise, we detected a 6.3-fold to 14.2-fold reduction in the proportion of BLAST hits assigned to bacteria following oligonucleotide spiking. This change in sequence composition suggests that a significant fraction of the reads consisted in oligonucleotide sequences that could not be filtered, possibly due to the presence of substantial levels of sequencing errors and/or the formation of chimeric sequencing reads.
Further analyses of the sequences revealed that the population of horse DNA sequences recovered following oligonucleotide spiking exhibited a skewed %GC contents (Figure 1, middle). This lower %GC-content (average = 38.5%) might have resulted from the presence of unfiltered oligonucleotide sequences. Following oligonucleotide spiking, we found levels of cytosine deamination similar to the ones observed in the absence of spiking at a denaturation temperature of 80°C (Figure 1, bottom). The observed misincorporation pattern was in agreement with that expected for ancient DNA molecules, ie. an excess of G to A mismatches at the 5’-ends of reads (Figure 3, bottom).
Despite its ability to access genuine ancient sequence information, the oligonucleotide spiking approach tested here could not be recommended for sequencing ancient DNA extracts unless further developments are performed in order to optimize the filtering of oligonucleotide sequences and to prevent the shifts observed both in base composition and sequence diversity.
DNA phosphatase treatment
A previous comparison of tSMS shotgun sequencing of modern versus ancient (Pleistocene horse bone) DNA extracts revealed no major difference in the base composition of the genomic region preceding sequencing reads . In particular, blocking sites were found to be enriched in guanine residues due to the preferential incorporation of dGTP Virtual Terminator during the locking reaction (Additional file 1, Figure S2). Additionally, genomic coordinates located upstream of blocking sites were enriched for thymine residues most likely resulting from the preferential capture of poly-A rich genomic regions by oligo-dT probes linked to the surface of the flow cell  (Additional file 1, Figure S2). In contrast, Illumina sequencing of the same extracts has confirmed the model of DNA fragmentation through depurination presented by Briggs and colleagues, with excessive proportions of guanine and adenine residues at genomic positions preceding the sequencing start site .
The absence of a clear signal for depurination in tSMS reads has been proposed to result from the terminal transferase enzyme used for poly-A tailing in the Helicos template preparation procedure, which might show non-optimal activity on modified 3’-ends of depurinated sites . After depurination the DNA backbone undergoes further cleavage at the site of base loss through β-elimination reactions, leaving DNA fragments with modified 3’-ends, which could be incompatible with poly-A tailing and further tSMS sequencing [4, 33]. Phosphorylation is one of the many modifications possibly present at 3’-ends of ancient templates. Hence, we decided to test whether dephosphorylation of ancient DNA extracts could restore access to a significant fraction of ancient templates and therefore could improve the overall performance of tSMS by increasing either the overall number of sequences generated per channel or the fraction of endogenous sequences recovered. As mild denaturation temperatures had been shown to deliver higher endogenous sequence yields (see above), we selected an enzyme with 3’ DNA phosphatase activity that could be heat inactivated at 65°C. In addition, the enzymatic buffer was selected to be compatible with a further polyA-tailing reaction in order to avoid extra template purification steps. For each of the four extracts generated from sample TP a fraction was treated with DNA phosphatase and the template preparation procedure was performed using 80°C as the denaturation temperature. Controls were run using an identical template preparation procedure but with no phosphatase treatment.
In all but one duplicate (extract TP2RE) DNA phosphatase treatment generated a 7.5 – 9.7-fold enrichment in the overall number of sequences produced per sequencing channel (Table 2). This enrichment was paralleled by a 2.0 – 2.9-fold decrease in the ratio of endogenous sequences (Table 2), suggesting that exogenous templates benefited more from the DNA phosphatase treatment than endogenous DNA molecules. This was observed for the four independent experiments performed and could result from a saturation effect whereby DNA phosphatase restored all ancient molecules with the phosphorylated 3’-ends that were present in the relatively small volume of extracts analyzed (0.5 μl).
Except for extract TP2RE, DNA phosphatase treatment resulted in a higher coverage of the horse genome per channel, as the increase in the overall number of sequences outweighed the loss resulting from the decrease in relative endogenous content. In this way the DNA phosphatase treatment allowed us to retrieve 15.7 Mb of sequence information from four 110 FOV channels, while only 6.7 Mb could be recovered in controls with a similar DNA denaturation temperature but no DNA phosphatase treatment. One sample (TP2RE) showed, however, lower absolute sequence yield after phosphatase treatment, most likely due to a channel-specific drop in sequencing quality. Even though all 110 FOV channels were sequenced on the same run, a series of known factors can lead to channel-to-channel variation, including 1) inhomogeneous deposition of primers (causing DNA reads to be located too close to each other at one extreme or not enough reads per image at the other extreme); 2) air bubbles in the channel that cause inefficient distribution of reagents, and; 3) laser focusing issues caused by inhomogeneity of the flow cell surface. With the channel considered, the overall yield of aligned sequences was lower than all the other channels, phosphatased or not; additionally, the rates of endogenous sequences appears lower than non-phosphatased samples of the same type, in agreement with all the other phosphatased samples (Table 2). This channel-specific drop in sequencing quality has limited the number of sequences recovered, but has not changed misincorporation patterns, cumulative damage distributions and %GC content, that all appear similar to what observed for sample TP1RE.
Following DNA phosphatase treatment the base composition at the genomic position preceding the sequencing start was found to be enriched in pyrimidines (position −1; Figure 4, middle). This position corresponds to the blocking site, which was previously shown to preferentially consist of guanine residues . In addition to affecting the base composition at the blocking site, the DNA phosphatase treatment was found to change the base composition at position −2 from the sequencing start, which corresponds to the position of the first nucleotide following the 3’-end of the ancient DNA strand (Additional file 1, Figure S2). After DNA phosphatase treatment this position was found to be enriched in cytosine residues, and depleted in equivalent proportions of adenine residues (Figure 4, top). As Helicos sequencing reads provide the sequence of the strand complementary to the ancient DNA template, this suggests that the first nucleotide following the 3’-end of the ancient DNA strand consisted of an excess of guanine residues, in agreement with the model of DNA fragmentation through depurination [34, 35], and with previous analyses using Illumina sequencing that have reported a preferential enrichment in guanine residues at the genomic coordinate preceding the beginning of sequencing reads . Hence, the DNA phosphatase treatment restored the signal of DNA depurination which was nearly absent from previous tSMS analyses of ancient DNA templates, confirming that modification of 3’-ends, in particular phosphorylation, precluded efficient poly-A tailing of ancient DNA templates with terminal transferase.
Furthermore, we noted that, following DNA phosphatase treatment, the base composition of the first nucleotide sequenced was enriched in cytosine residues (Figure 4, bottom). This enrichment occurred at the expense of adenine residues. Overall, horse reads recovered from the DNA phosphatase treatment were found to be slightly shorter than reads recovered without DNA phosphatase treatment (Figure 2, top) and the base composition of these reads was similar to those recovered following the use of a 95°C denaturation temperature (Figure 2, middle). The fraction of horse reads recovered after treatment of TP1 and TP2 extracts with DNA phosphatase exhibited G→A cumulative misincorporation rates that were lower than those observed without DNA phosphatase treatment, suggesting that DNA phosphatase allowed preferential access to ancient DNA templates with lower levels of cytosine deamination (Figure 2, bottom). More specifically, G→A cumulative misincorporation rates ranged from 38.1% – 45.1% in the first 25 nucleotides sequenced following DNA phosphatase treatment in contrast to 45.9% – 61.0% in the absence of this treatment. No significant difference in overall G→A cumulative misincorporation rates were found in horse sequences recovered from extracts TP1RE and TP2RE, despite different DNA denaturation temperatures and/or DNA phosphatase treatment.
We should caution that DNA phosphatase treatments of ancient DNA extracts could be deleterious if second-generation sequencing methods and not tSMS were used as sequencing approaches, since the former require building ancient templates into DNA libraries through adapter ligation at phosphorylated 5’-ends. This has been demonstrated by Briggs et al.  who have shown that the fraction of mammoth reads was drastically reduced following incubation with Calf Intestinal Alkaline Phosphatase in the absence of further rephosphorylation of 5’ ends with PNK. Conversely, treatment with DNA phosphatase increased the amount of templates accessible to DNA polymerases used for filling-in adapter sequences by restoring 3’ hydroxyl termini. Therefore the pre-treatment of ancient DNA extracts with DNA phosphatase before library building, followed by a proper end repair step rephosphorylating 5’ ends, could be used for improving our ability to insert ancient DNA fragments into DNA libraries and consequently for improving the performance of 454 and Illumina sequencing of ancient DNA templates. The exact benefit of this approach for second-generation sequencing, however, remains to be determined. We anticipate that it will prove highly dependent on the relative efficiency of Antarctic phosphatase over T4 Polynucleotide Kinase, as the latter is part of the standard end-repair enzymatic cocktail used for library building and could also catalyse the removal of 3’-phosphoryl groups from polynucleotides.
Overall, we have demonstrated that sequencing of DNA purified from redigested bone powder, and the use of DNA template preparation procedures combining a mild DNA denaturation temperature and DNA phosphatase treatment, significantly increased the amount of endogenous sequences retrieved from ancient bones using tSMS. Following these findings, we decided to perform further tSMS sequencing on the TP1RE extract using our modified DNA template preparation procedure (Table 2). A total of eight channels were used on a 550 FOV platform, generating 46.4 Mb of horse DNA sequences. When added to the information recovered from the twelve 110 FOV channels this yielded 75.2 Mb of horse DNA sequence for the 13ky-old sample (sample TP). A total of 42,253 nucleotides could be aligned with no ambiguity on the horse mitochondrial genome, corresponding to a global coverage of 2.53X (Accession Nb. NC_001640).
We subsequently used the mitochondrial information to assess further the quality of tSMS sequences. We first focused on nucleotide positions (sites) covered with a depth of 2 and calculated the proportion of sites showing identical bases among reads as a conservative estimate of the quality of the data generated by tSMS on ancient DNA templates. In addition, we further considered sites showing a minimal sequencing depth of 3 and calculated the proportion of these with >50% base agreement among reads. Ignoring indels, we found identical sequences in 97.8% of the sites covered at a depth of 2. Furthermore, levels of accuracy ranging from 99.97% – 100.0% were found at sequence depths equal to or greater than 3 (a total of 7,510 sites), suggesting that a strict consensus approach considering sites with minimal sequencing depth of 2 should deliver high quality mitochondrial data. We consequently built a Maximum Likelihood phylogenetic tree using such sites (11,316 in total) together with donkey and horse mitogenomic data available from the literature [36–40]. As expected, the late Pleistocene mitogenome sequence data recovered from sample TP was found to cluster among horse sequences with maximal bootstrap support (Figure 5). Interestingly, one haplotype characterized from the wild Przewalski’s horse exhibited maximum sequence similarity to the late Pleistocene sequence. However, the latter grouping was not supported by high bootstrap values, but is in agreement with the notion of this haplotype originating from a Przewaski’s horse specimen that is thought to have no recent domestic horse genetic contribution to its known pedigree .
In this study we have confirmed that large bone crystals that resist a first digestion step in EDTA-rich buffers provide access to molecular preservation niches in ancient bones, resulting in high endogenous sequence contents as well as lower levels of cytosine deamination. Pending further confirmation for a wider range of preservation conditions, such niches could be targeted as a means of optimizing the characterization of complete ancient genome sequences through shotgun sequencing. This further shows that non-destructive extraction methods consisting of smooth digestion in decalcifying buffers do not provide access to the optimal fraction of DNA templates preserved in bone specimens. In addition, we show that even simple temperature shifts and/or enzymatic treatments, have a major impact on the population of molecules retrieved from tSMS. Overall, this suggests that the amount of ancient DNA sequences that are accessible using current methods could potentially be increased by improving extraction procedures and restoring accessibility to 3’ termini. We therefore suggest that DNA from large bone crystals should be preferentially extracted, and when tSMS is used as a sequencing method, DNA templates should be prepared following DNA phosphatase treatment and using mild denaturation temperatures. Spiking of ancient DNA extracts with oligonucleotides of known sequences also appears to be a promising means to increase tSMS yields further.
The DNA extract used for testing the extent to which oligonucleotide spiking could improve the performance of tSMS was prepared from a horse fossil bone excavated at Thistle Creek in Pleistocene permafrost deposits. The sample consists of two pieces of bone (CA1 and CA2) collected from a larger bone (TC21 YG148.2) that has already been successfully analyzed by tSMS and that has been associated with an infinite radiocarbon date (OxA-23933 >50,300BP; and UBA-16493 and UBA-17013 >50,505BP) . The DNA phosphatase experiments were performed on DNA extracts from a humerus bone fragment from another horse fossil (CGG10023 XA08-53) that has been preserved in the permafrost (Taymir peninsula) and that has been associated with a finite radiocarbon date (UBA-16479 13,389 ± 52BP). All fossil bones were extracted using ancient DNA facilities at the Centre for GeoGenetics, Copenhagen, Denmark. Samples were drilled to powder at low speed and this bone powder was then extracted at 37°C under agitation using a 48 hours digestion in 5 ml of EDTA buffer (EDTA 0.5 M, pH = 8.0; N-lauryl Sarcosyl 0.5%; Proteinase K 1 mg/ml). A total of 0.45 and 0.46 grams of bone powder were digested for samples CA1 and CA2, whereas the digestion was replicated twice on 4.0 and 2.2 grams of bone powder for sample CGG10023. The CGG10023 sample derived from the initial digestion of 4.0 grams of bone powder will be referred to in the main text as TP1 and that from the 2.2 grams sample as TP2. The supernatants of these first digestions, recovered after spinning the digestion solution at 2,000 rpm for 5 minutes, were transferred to 30 kDa Amicon centrifugal filter units (Millipore) and spun at 3,000 rpm for a minimum of 30 minutes. A volume of ca. 200 μl was recovered from each Centricon column and further concentrated to 150 μl in an EB buffer using Minelute columns (QIAgen) following the manufacturer’s instructions except that the elution was performed 15 minutes after loading the EB buffer on the column and incubating the column at 37°C. The bone pellets that were left undigested after this first extraction of TP1 and TP2 were further extracted by adding 5 ml of EDTA buffer and following the same procedure as above. These re-extractions will be referred to TP1RE and TP2RE in the main text.
Helicos tSMS reactions were performed at Helicos BioSciences Corporation facilities in Cambridge, Massachusetts, USA. For samples CA1 and CA2, we followed the same methodology as described in  using 8 μl of DNA extract and two different DNA denaturation temperatures while preparing DNA templates for tSMS (80°C and 95°C, respectively). Two microliters of the 95°C extracts and 4 μl of the extracts denatured at 80°C were loaded per channel on a 1100 FOV sequencing platform using a total number of four channels (Table 1). In addition, 6 μl of the CA1 95°C extract and 4 μl of the CA1 80°C extract were pooled, spiked with polyA-tailed oligonucleotides (final concentration 40 pM) and loaded on one channel of the 1100 FOV platform. A similar procedure was followed for the CA2 extracts. This procedure aimed at ensuring that there were enough DNA molecules being sequenced to provide good image registration.
Since different volumes of the 80°C and 95°C denatured extracts were pooled, we estimated the expected minimal endogenous ratio Rmin for spiking experiments using a weighted average of the endogenous ratios observed without spiking. For each spiking experiment, we further assumed a binomial distribution with Rmin as the probability of success for estimating the number of expected endogenous sequences given the total number of sequences generated. Five percent and 95% quantiles were estimated using the function qbinom in R .
The effects of different denaturation temperatures and of DNA phosphatase treatment were tested on the four extracts recovered from sample CGG10023 (TP1, TP1RE, TP2 and TP2RE) using 110 FOV and 550 FOV. For the phosphatase reactions, 0.5 μl of the DNA extracts were mixed with 8 μl of nuclease-free water, 1 μl of NEB Antarctic Phosphatase 10X Reaction buffer and 2.5 units of NEB Antarctic Phosphatase. Reactions were incubated for 30 minutes at 37°C. Reactions were stopped by inactivating the enzyme at 65°C for 15 minutes. For the poly-A tailing reaction, the heat-inactivated phosphatase reactions were mixed with 0.8 μl of nuclease-free water, 2 μl of NEB Terminal Transferase 10X buffer and 2 μl of a 2.5 mM CoCl2 solution. This mixture was heated at 80°C for 5 minutes in a thermocycler for denaturation. For DNA extracts that were not phosphatased, 0.5 μl of the DNA extracts were mixed with 10.3 μl of nuclease-free water, 2 μl of NEB Terminal Transferase 10X buffer and 2 μl of a 2.5 mM CoCl2 solution. This was then heated at 80°C or 95°C for 5 minutes for denaturation. All the extracts were rapidly cooled on ice to minimize re-annealing of the denatured DNA strands. For all poly-A tailing reactions the volume of the previous mixes was increased to 20 μl through the addition of 5 units of NEB Terminal Transferase, NEB BSA (to 1X final concentration) and dATP at a final concentration of 10 μM in the 20 μl reaction. The single stranded DNA molecules were poly-A tailed for 1 hr at 37°C. Reactions were stopped by inactivating the enzyme at 70°C for 10 minutes. As the DNA is prone to re-annealing during the tailing step, heating at 80°C or 95°C followed by rapid cooling was repeated before 10 μl of 3’-end blocking mastermix (NEB Terminal Transferase buffer, 250 μM CoCl2, 5U NEB Terminal Transferase, 10 μM ddATP) was added to the tailing reaction volume. The 3’-end blocking reactions were performed for 1 hr at 37°C and stopped by denaturing the enzyme at 70°C for 20 minutes. DNA that may have reannealed during the blocking reaction was converted back to single strands by repeating the previous heating-rapid cooling conditions prior to loading the samples on the flowcell. 20 μl of sample consisting of 1 μl of the tailed and blocked DNA reaction, 9 μl of nuclease-free water and 10 μl of 2X hybridization buffer was added to each channel and allowed to hybridize for 1 hr at 37°C. The buffer was then rinsed away and the extra bases of the poly-A tailed reaction were filled in with dTTP and then locked in place with the first non-TTP base . Sequencing was carried out using Virtual Terminator© Nucleotides as described in . The resulting sequence reads were then filtered for minimal sizes of 25 nucleotides and for artifactual sequences that closely matched the order of nucleotide addition. Sequence reads starting with a minimum number of 2 Ts were further trimmed as the latter could arise from incomplete fill and lock reactions. The remaining set of filtered reads was then analyzed.
The DNA sequence data analyzed in this study are available on the NCBI Sequence Read Archive (SRA Accession Number SRA045862). Helicos reads were mapped using BWA against two equid mitochondrial genomes (Accession numbers: NC_001640 for the horse, and NC_001788 for the donkey), the human reference genome (hg19) and the horse reference genome (equCab2, filtered for the mitochondrial genome and chromosome Un), each available at the UCSC Genome Bioinformatic website: http://genome.ucsc.edu/). Seed was disallowed for mapping using the –l option and providing a value larger than the size of the longest sequence reads. In addition, as tSMS shows high error rates, especially indels [43, 44], gaps were allowed all along the reads, including at termini, by setting the -i option to zero . The BWA output was further converted to sam format with the samse command and filtered for mapping quality scores higher than 25 using the view command in samtools. Reads mapping uniquely to the horse reference genome but showing no hit to the human reference were retained. This conservative approach was selected in order to filter for possible contamination with human DNA and remove spurious hits that could result from short sequencing reads. BLAST analyses of sequences obtained for samples CA1 and CA2, after filtering of spiking oligonucleotides for spiking experiments were performed using blastn  against the NCBI non-redundant database with the following options: -F F –e 0.01 –q −2, i.e., no masking for repeat regions, expect value threshold 0.01 and mismatch penalty adjusted to −2. Base compositions, fragmentation and nucleotide misincorporation patterns were generated using SAM files as input and the map Damage package . %GC contents and cumulative damage rates were plotted using mapDamage outputs and in-house-made R scripts .
Maximum Likelihood phylogenetic inferences were performed using PhyML 3.0  and the robustness of nodes was estimated from 1,000 bootstrap pseudo-replicates. The best-fit model of substitution (HKY + G8) was selected in jModeltest  using the Akaike Information Criterion (AICc). The ML phylogeny was plotted using FigTree v1.3.1 available at http://tree.bio.ed.ac.uk/software/figtree.
True Single Molecule DNA Sequencing
Field of View.
Willerslev E, Hansen AJ, Ronn R, Brand TB, Barnes I, Wiuf C, Gilichinsky D, Mitchell D, Cooper A: Long-term persistence of bacterial DNA. Curr Biol. 2004, 14: R9-R10. 10.1016/j.cub.2003.12.012.
Willerslev E, Cappellini E, Boosma W, Nielsen R, Hebsgaard MB, Brand TB, Hofreiter M, Bunce M, Poinar HN, Dahl-Jensen D, Johnsen S, Steffensen JP, Bennike O, Schwenninger JL, Nathan R, Armitage S, de Hoog CJ, Alfimov V, Christl M, Beer J, Muscheler R, Barker J, Sharp M, Penkman KE, Haile J, Taberlet P, Gilbert MT, Casoli A, Campani E, Collins MJ: Ancient biomolecules from deep ice cores reveal a forested southern Greenland. Science. 2007, 317: 111-114. 10.1126/science.1141758.
Valdiosera C, Garcia N, Dalen L, Smith C, Kahke RD, Liden K, Anderbjorn A, Arsuaga JL, Gotherstrom A: Typing single polymorphic nucleotides in mitochondrial DNA as a way to access Middle Pleistocene DNA. Biol Lett. 2006, 22: 601-603.
Lindahl T: Instability and decay of the primary structure of DNA. Nature. 1993, 362: 709-715. 10.1038/362709a0.
Gilbert MT, Hansen AJ, Willerslev E, Rudbeck L, Barnes I, Lynerup N, Cooper A: Characterization of genetic miscoding lesions caused by postmortem damage. Am J Hum Genet. 2003, 72: 48-61. 10.1086/345379.
Paabo S, Poinar H, Serre D, Jaenicke-Despres V, Hebler J, Rohland N, Kuch M, Krause J, Vigilant L, Hofreiter M: Genetic analyses from ancient DNA. Annu Rev Genet. 2004, 38: 645-679. 10.1146/annurev.genet.37.110801.143214.
Hansen AJ, Mitchell DL, Wiuf C, Paniker L, Brand TB, Binladen J, Gilichinsky DA, Ronn R, Willerslev E: Crosslinks rather than strand breaks determine access to ancient DNA sequences from frozen sediments. Genetics. 2006, 173: 1175-1179. 10.1534/genetics.106.057349.
Heyn P, Stenzel U, Briggs AW, Kircher M, Hofreiter M, Meyer M: Road-blocks on paleogenomes—polymerase extension profiling reveals the frequency of blocking lesions in ancient DNA. Nucleic Acids Res. 2010, 38: e161-10.1093/nar/gkq572.
Poinar HM, Hofreiter M, Spaulding WG, Martin PS, Stankiewicz BA, Bland H, Evershed RP, Possnert G, Paabo S: Molecular coproscopy: dung and diet of the extinct ground sloth Nothrotheriops shastensis. Science. 1998, 281: 402-406.
Orlando L, Ginolhac A, Raghavan M, Vilstrup J, Rasmussen M, Magnussen K, Steinmann KE, Kapranov P, Thompson JF, Zazula G, Froese D, Moltke I, Shapiro B, Hofreiter M, Al-Rasheid KA, Gilbert MT, Willerslev E: True Single-molecule DNA sequencing of a Pleistocene horse bone. Genome Res. 2011, 21: 1705-1719. 10.1101/gr.122747.111.
Noonan JP, Hofreiter M, Smith D, Priest JR, Rohland N, Rabeder G, Krause J, Detter JC, Paabo S, Rubin EM: Genomic sequencing of Pleistocene cave bears. Science. 2005, 309: 597-599. 10.1126/science.1113485.
Green RE, Krause J, Ptak SE, Briggs AW, Ronan MT, Simons JF, Du L, Egholm M, Rothberg JM, Paunovic M, Paabo S: Analysis of one million base pairs of Neanderthal DNA. Nature. 2006, 444: 330-336. 10.1038/nature05336.
Noonan JP, Coop G, Kudaravalli S, Smith D, Krause J, Alessi J, Chen F, Platt D, Paabo S, Pritchard JK, Rubin EM: Sequencing and analysis of Neanderthal genomic DNA. Science. 2006, 314: 1113-1118. 10.1126/science.1131412.
Gilbert MT, Tosho LP, Rendulic S, Packard M, Drautz DI, Sher A, Tikhonov A, Dalen L, Kuznetsova T, Kosintsev P, Campos PF, Higham T, Collins MJ, Wilson AS, Shidlovskiy F, Buigues B, Ericson PG, Germonpre M, Gotherstrom A, Iacumin P, Nicolaev V, Nowak-Kemp M, Willerslev E, Knight JR, Irzyk GP, Perbost CS, Fredrikson KM, Harkins TT, Sheridan S, Miller W, Schuster SC: Whole-genome shotgun sequencing of mitochondria from ancient hair shafts. Science. 2007, 317: 1927-1930. 10.1126/science.1146971.
Willerslev E, Gilbert MT, Binladen J, Ho SY, Campos PF, Ratan A, Tomsho LP, da Fonseca RR, Sher A, Kuznetsova TV, Nowak-Kemp M, Roth TL, Miller W, Schuster SC: Analysis of complete mitochondrial genomes from extinct and extant rhinoceroses reveals lack of phylogenetic resolution. BMC Evol Biol. 2009, 9: 95-10.1186/1471-2148-9-95.
Poinar HN, Schwarz C, Qi J, Shapiro B, MacPhee RDE, Buigues B, Tikhonov A, Huson DH, Tomsho LP, Auch A, Rampp M, Miller W, Schuster SC: Metagenomics to Paleogenomics: Large-scale sequencing of mammoth DNA. Science. 2006, 311: 392-394. 10.1126/science.1123360.
Reich D, Green RE, Kircher M, Krause J, Patterson N, Durand EY, Viola B, Briggs AW, Stenzel U, Johnson PL, Maricic T, Good JM, Marques-Bonet T, Alkan C, Fu Q, Mallick S, Li H, Meyer M, Eichler EE, Stoneking M, Richards M, Talamo S, Shunkov MV, Derevianko AP, Hublin JJ, Kelso J, Slatkin M, Paabo S: Genetic history of an archaic hominin group from Denisova Cave in Siberia. Nature. 2010, 468: 1053-1060. 10.1038/nature09710.
Briggs AW, Good JM, Green RE, Krause J, Maricic T, Stenzel U, Lalueza-Fox C, Rudan P, Brajkovic D, Kucan Z, Gusic I, Schmidz R, Doronichev VB, Golovanova LV, de la Rasilla M, Frotea J, Rozas A, Paabo S: Targeted retrieval and analysis of five Neandertal mtDNA genomes. Science. 2009, 325: 318-321. 10.1126/science.1174462.
Krause J, Briggs AW, Kircher M, Maricic T, Zwyns N, Derevianko A, Paabo S: A complete mtDNA genome of an early modern human from Kostenki, Russia. Curr Biol. 2010, 20: 231-236. 10.1016/j.cub.2009.11.068.
Schueneman VJ, Bos K, Dewitte S, Schmedes S, Jamieson J, Mittnik A, Forrest S, Coombes BK, Wood JW, Earn DJ, White W, Krause J, Poinar HN: Targeted enrichment of ancient pathogens yielding the pPCP1 plasmid of Yersinia pestis from victims of the Black Death. Proc Natl Acad Sci USA. 2011, 108: E746-10.1073/pnas.1105107108.
Burbano HA, Hodges E, Green RE, Briggs AW, Krause J, Meyer M, Good JM, Maricic T, Johnson PLF, Xuan Z, Rooks M, Bhattacharjee A, Brizuela L, Albert FW, de la Rasilla M, Fortea J, Rosas A, Lachmann M, Hannon GJ, Paabo S: Targeted investigation of the Neanderthal genome by array-based sequence capture. Science. 2010, 328: 723-725. 10.1126/science.1188046.
Fordyce SL, Ávila-Arcos MC, Rockenbauer E, Børsting C, Frank-Hansen R, Torp Petersen F, Willerslev E, Hansen AJ, Morling N, Gilbert MTP: Highthroughput sequencing of core STR loci for forensic genetic investigations using the Roche Genome Sequencer FLX platform. Biotechniques. 2011, 51: 127-133.
Miller W, Drautz DI, Ratan A, Pusey B, Qi J, Lesk AM, Tomsho LP, Packard MD, Zhao F, Sher A, Tikhonov A, Raney B, Patterson N, Lindblad-Toh K, Lander ES, Knight JR, Irzyk GP, Frederikson KM, Harkins TT, Sheridan S, Pringle T, Schuster SC: Sequencing the nuclear genome of the extinct woolly mammoth. Nature. 2008, 456: 387-390. 10.1038/nature07446.
Rasmussen M, Li Y, Lindgreen S, Pedersen JS, Albrechtsen A, Moltke I, Metspalu M, Metspalu E, Kivisild T, Gupta R, Bertalan M, Nielsen K, Gilbert MT, Wang Y, Raghavan M, Campos PF, Kamp HM, Wilson AS, Gledhill A, Tridico S, Bunce M, Lorenzen ED, Binladen J, Guo X, Zhao J, Zhang X, Zhang H, Li Z, Chen M, Orlando L, Kristiansen K, Bak M, Tommerup N, Bendixen C, Pierre TL, Gronnow B, Meldgaard M, Andreasen C, Fedorova SA, Osipova LP, Higham TF, Ramsey CB, Hansen TV, Nielsen FC, Crawford MH, Brunak S, Sicheritz-Ponten T, Villems R, Nielsen R, Krogh A, Wang J, Willerslev E: Ancient human genome sequence of an extinct paleo-eskimo. Nature. 2010, 463: 757-762. 10.1038/nature08835.
Green RE, Krause J, Briggs AW, Maricic T, Stenzel U, Kircher M, Patterson N, Li H, Zhai W, Fritz MH, Hansen NF, Durand EY, Malaspinas AS, Jensen JD, Marques-Bonet T, Alkan C, Prüfer K, Meyer M, Burbano HA, Good JM, Schultz R, Aximu-Petri A, Butthof A, Höber B, Höffner B, Siegemund M, Weihmann A, Nusbaum C, Lander ES, Russ C, Novod N, Affourtit J, Egholm M, Verna C, Rudan P, Brajkovic D, Kucan Z, Gusic I, Doronichev VB, Golovanova LV, Lalueza-Fox C, de la Rasilla M, Fortea J, Rosas A, Schmitz RW, Johnson PL, Eichler EE, Falush D, Birney E, Mullikin JC, Slatkin M, Nielsen R, Kelso J, Lachmann M, Reich D, Pääbo S: A draft sequence of the neandertal genome. Science. 2010, 328: 710-722. 10.1126/science.1188021.
Salamon M, Tuross N, Arensburg B, Weiner S: Relatively well preserved DNA is present in the crystal aggregates of fossil bones. Proc Natl Acad Sci Unit States Am. 2005, 102: 13783-13788. 10.1073/pnas.0503718102.
Hedges REM: Bone diagenesis: an overview of processes. Archaeometry. 2002, 44: 319-328. 10.1111/1475-4754.00064.
Collins MJ, Nielsen-Marsh CM, Hiller Jm Smith CI, Roberts JP, Prigodich RV, Weiss TJ, Csapo J, Millard AR, Turner-Wlaker G: The survival of organic matter in bone: a review. Archaeometry. 2002, 44: 383-394. 10.1111/1475-4754.t01-1-00071.
Campos PF, Craig OE, Turner-Walker G, Peacock E, Gilbert MT: DNA in ancient bones - where is it located and how should we extract it?. Ann Anat. 2011, 194: 7-16.
Handt O, Krings M, Ward RH, Paabo S: The retrieval of ancient human DNA sequences. Am J Hum Genet. 1996, 59: 368-376.
Hart C, Lipson D, Ozsolak F, Raz T, Steinmann K, Thompson J, Milos PM: Single-molecule sequencing: sequence methods to enable accurate quantitation. Methods Enzymol. 2010, 472: 407-430.
Briggs AW, Stenzel U, Johnson PLF, Green RE, Kelso J, Prufer K, Meyer M, Krause J, Ronan MT, Lachmann M, Paabo S: Patterns of damage in genomic DNA sequences from a Neandertal. Proc Natl Acad Sci Unit States Am. 2007, 104: 14616-14621. 10.1073/pnas.0704665104.
Mitchell D, Willerslev E, Hansen A: Damage and repair of ancient DNA. Mutat Res. 2005, 571: 265-276. 10.1016/j.mrfmmm.2004.06.060.
Brotherton P, Endicott P, Sanchez JJ, Beaumont M, Barnette R, Austin J, Cooper A: Novel high-resolution characterization of ancient DNA reveals C > U-type base modification events as the sole cause of post-mortem miscoding lesions. Nucleic Acids Res. 2007, 35: 5717-5728. 10.1093/nar/gkm588.
Briggs AW, Stenzel U, Meyer M, Krause J, Kircher M, Paabo S: Removal of deaminated cytosines and detection of in vivo methylation in ancient DNA. Nucleic Acids Res. 2010, 38: e87-10.1093/nar/gkp1163.
Xu X, Arnason U: The complete mitochondrial DNA sequence of the horse Equus caballus: extensive heteroplasmy of the control region. Gene. 1994, 148: 357-362. 10.1016/0378-1119(94)90713-7.
Xu X, Gullberg A, Arnason U: The complete mitochondrial DNA (mtDNA) of the donkey and mtDNA comparisons among four closely related mammalian species-pairs. J Mol Evol. 1996, 43: 438-446. 10.1007/BF02337515.
Xu S, Luosang J, Hua S, He J, Ciren A, Wang W, Tong X, Liang Y, Wang J, Zheng X: High altitude adaptation and phylogenetic analysis of Tibetan horse based on the mitochondrial genome. J Genet Genom. 2007, 34: 720-729. 10.1016/S1673-8527(07)60081-2.
Jiang Q, Wei Y, Huang Y, Jiang H, Guo Y, Lan G, Liao DJ: The complete mitochondrial genome and phylogenetic analysis of the Debao pony (Equus caballus). Mol Biol Rep. 2010, 38: 593-599.
Goto H, Ryder OA, Fisher AR, Schultz B, Pond SK, Nekrutenko A, Makova KD: A massively parallel sequencing approach uncovers ancient origins and high genetic variability of endangered Przewalski’s horses. Genome Biol Evol. 2011, 3: 1096-1106. 10.1093/gbe/evr067.
R development Core Team: R: a language and environment for statistical computing. Genome Biol Evol. R Foundation for Statistical Computing, Vienna, Austria. 2010. ISBN 3-900051-07-0, URL http://www.R-project.org
Lipson D, Raz T, Kieu A, Jones DR, Giladi E, Thayer E, Thompson JF, Letovsky S, Milos P, Causey M: Quantification of the yeast transcriptome by single-molecule sequencing. Nat Biotechnol. 2009, 27: 652-658. 10.1038/nbt.1551.
Bowers J, Mitchell J, Beer E, Buzby PR, Causey M, Efcavitch JW, Harosz M, Krzymanska-Olejnik E, Kung L, Lipson D, Lowman GM, Marappan S, MaInerney P, Platt A, Roy A, Siddiqi SM, Steinmann K, Thompson JF: Virtual terminator nucleotides for next-generation DAN sequencing. Nat Methods. 2009, 6: 593-595. 10.1038/nmeth.1354.
Thompson JF, Steinmann KE: Single molecular sequencing with a HeliScope genetic analysis system. Curr Protoc Mol Biol. 2010, 7: unit 7.10-
Li H, Durbin R: Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics. 2009, 25: 1754-1760. 10.1093/bioinformatics/btp324.
Altschul SF, Madden TL, Schäffer AA, Zhang J, Zhang Z, Miller W, Lipman DJ: Gapped BLAST and PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res. 1997, 25: 3389-3402. 10.1093/nar/25.17.3389.
Ginolhac A, Rasmussen M, Gilbert MTP, Willerslev E, Orlando L: mapDamage: testing for damage patterns in ancient DNA sequences. Bioinformatics. 2011, 27: 2153-2155. 10.1093/bioinformatics/btr347.
Guindon S, Dufayard JF, Lefort V, Anisimova M, Hordijk W, Gascuel O: New algorithms and methods to estimate maximum-likelihood phylogenies: assessing the performance of PhyML 3.0. Syst Biol. 2010, 59: 307-321. 10.1093/sysbio/syq010.
Posada D: jModeltest: phylogenetic model averaging. Mol Biol Evol. 2008, 25: 1253-1256. 10.1093/molbev/msn083.
We thank Tina Brand and the laboratory technicians at the Danish High-throughput DNA Sequencing Centre for technical assistance; Anders Krogh and Thomas Sicheritz-Ponten for access to computation facilities; Anders Krogh, Stinus Lindgreen and Mikkel Schubert for fruitful discussions related to read mapping; Maanasa Raghavan for constructive suggestions regarding manuscript editing; two anonymous reviewers for valuable input on a previous version of the manuscript. This work was supported by the Danish Council for Independent Research, Natural Sciences (FNU); the Danish National Research Foundation; a Marie-Curie Career Integration Grant CIG- 293845; the National Science Foundation ARC-0909456; the Searle Scholars Program; King Saud University Distinguished Scientist Fellowship Program (DSFP). This project was supported in part by American Recovery and Reinvestment Act (ARRA) funds through grant number RC2 HG005598 to JFT from the National Human Genome Research Institute, National Institutes of Health. Its contents are solely the responsibility of the authors and do not necessarily represent the official views of NHGRI.
LO and JT designed research. LO, JV and JS extracted ancient DNA. LO and AG analyzed tSMS data. LO and JV performed phylogenetic analyses. DF, GZ, BS, MS contributed fossils. KS and JT performed tSMS. KAH, KS, JT, TG, EW contributed methods, reagents and tools, and valuable input for data interpretation. LO wrote the manuscript. All authors read and approved the final manuscript.