PaperlayerView exact revision on arXiv

Source-linked AI summary

Chromosome-scale shotgun assembly using an in vitro method for long-range linkage

Nicholas H. Putnam, Brendan O'Connell, Jonathan C. Stites, Brandon J. Rice, Andrew Fields, Paul D. Hartley, Charles W. Sugnet, David Haussler, Daniel S. Rokhsar, Richard E. Green

arXiv:1502.05331v1q-bio.GNq-bio.BM

TL;DR

Long-range, accurate genome assembly from short reads remains challenging. This paper introduces Chicago, an in vitro method for generating long-range mate-pair data, and shows that it can substantially improve de novo genome scaffolding.

  • Problem

    Accurate, long-range scaffolding of de novo assemblies from high-throughput short-read sequencing remains a central genomics challenge.

  • Method

    The authors generate Chicago long-range mate-pair data in vitro and use it with short-fragment Illumina sequencing to scaffold vertebrate genomes.

  • Results

    Chicago data dramatically improved scaffolding of de novo assembled vertebrate genomes using short-fragment Illumina sequence plus Chicago libraries.

  • Takeaways & Limitations

    The approach provides a practical route to highly contiguous vertebrate genome scaffolding from short-read sequencing and in vitro libraries.

  • Takeaways & Limitations

    The authors report that patents have been filed on the technology described in the manuscript.

Abstract

from arXiv · show

Long-range and highly accurate de novo assembly from short-read data is one of the most pressing challenges in genomics. Recently, it has been shown that read pairs generated by proximity ligation of DNA in chromatin of living tissue can address this problem. These data dramatically increase the scaffold contiguity of assemblies and provide haplotype phasing information. Here, we describe a simpler approach ("Chicago") based on in vitro reconstituted chromatin. We generated two Chicago datasets with human DNA and used a new software pipeline ("HiRise") to construct a highly accurate de novo assembly and scaffolding of a human genome with scaffold N50 of 30 Mb. We also demonstrated the utility of Chicago for improving existing assemblies by re-assembling and scaffolding the genome of the American alligator. With a single library and one lane of Illumina HiSeq sequencing, we increased the scaffold N50 of the American alligator from 508 kb to 10 Mb. Our method uses established molecular biology procedures and can be used to analyze any genome, as it requires only about 5 micrograms of DNA as the starting material.

Results

Chicago libraries provided long-range genomic links, enabling HiRise to improve human assembly scaffolding and alligator assembly contiguity while supporting haplotype phasing and structural-variant detection.

  • Human genome assembly: HiRise produced human genomic scaffolds that were longer and had fewer global mis-assemblies than published MERACULOUS and APLG assemblies.HiRise used a likelihood model and Chicago links to break and re-scaffold contigs.
  • American alligator assembly: 10.3 Mbp: The HiRise-scaffolded American alligator assembly achieved this scaffold N50 from a de novo assembly with N50 81 Kbp.A single Chicago library was sequenced with 210.7 million reads.
  • American alligator assembly: 98.1%: BAC end pairs in the alligator HiRise assembly were aligned on the same scaffold in correct relative orientation.The assembly had a global density of misjoins of less than 1 per 8.36 Mbp.
  • Haplotype phasing: 99.83%: Haplotype-informative read pairs spanning 10 Kbp to 150 Kbp agreed with the known GM12878 haplotype phase.Chicago data therefore provided haplotype-phasing information across these separations.
  • Structural variation: 0.81 sensitivity and 0.76 specificity: Chicago data detected 5 Kbp inversions.The analysis also revealed two structural differences when mapped to GRCh38.

Discussion

Chicago is a simple in vitro library and bioinformatic approach that improves de novo genome scaffolding while avoiding living material and confounding biological proximity signals. The method supports long-range, highly contiguous scaffolding, genome-variation discovery, and progress toward haplotype-resolved chromosome reconstruction.

  • Method advantages: Chicago library construction requires no living biological material and uses 5.0 micrograms of input DNA per library.The libraries described were each generated from 5.0 micrograms of input DNA.
  • Method advantages: In vitro chromatin proximity ligation avoids persistent biological proximity signals that can confound genome assembly.The authors report low background noise and a virtual absence of persistent, spurious read pairs.
  • Scaffolding performance: Chicago libraries support highly contiguous vertebrate-genome scaffolding using short-fragment Illumina sequence, with pair separation limited by input-DNA molecular weight.Chicago successfully scaffolded assemblies with 30 Kbp N50 input contigs.
  • Method advantages: The approach eliminates the need to create and sequence combined long-range mate-pair and fosmid libraries or use specialized high-molecular-weight DNA equipment.This simplifies library preparation compared with methods requiring shearing or size selection.
  • Implications: Chicago generates longer-range genome-assembly scaffolds than existing methods and also supports genome-variation discovery.The authors frame these results as progress toward accurate reconstruction of full-length, haplotype-resolved chromosome sequences.

Methods … Biotinylation and restriction digestion

The method extracted genomic DNA, assembled and validated chromatin in vitro, then biotinylated, fixed, dialyzed, restriction-digested, and captured the chromatin on streptavidin beads. These steps prepared chromatin-linked DNA fragments for downstream processing.

  • DNA Preparation: Genomic DNA was extracted by lysing cells, isolating nuclei, digesting with Proteinase K and RNAse A, purifying on a Qiagen column, and resuspending the pellet in 200 µL TE.DNA was washed, eluted, precipitated in isopropanol, pelleted by centrifugation, dried, and resuspended.
  • Chromatin assembly: Chromatin was assembled overnight at 27°C from genomic DNA using the Active Motif in vitro Chromatin Assembly kit.A 10% aliquot was subjected to MNase digestion to confirm successful chromatin assembly.
  • Biotinylation and restriction digestion: Chromatin was biotinylated with iodoacetyl-PEG-2-biotin (IPB), then fixed in 1% formaldehyde at room temperature for 15 minutes.Fixation was followed by quenching with a 2-fold molar excess of 2.5M Glycine.
  • Biotinylation and restriction digestion: Excess IPB and cross-linked glycine were removed by dialysis against 1L of dialysis buffer at 4°C for a minimum of 3 hours.The dialysis buffer contained 10mm Tris-Cl, pH8.0, and 1mM EDTA.
  • Biotinylation and restriction digestion: Chromatin was digested with either MboI or MluCI in 1X CutSmart for 4 hours at 37°C.It was subsequently dialyzed at 4°C for 2 hours and again overnight with fresh buffer to remove enzyme and short, free fragments.
  • Biotinylation and restriction digestion: Dynabead MyOne C1 streptavidin beads were washed, resuspended in PBS + 0.1% Tween-20, added to chromatin, and incubated for 1 hour at room temperature.The beads were then concentrated magnetically, washed, re-concentrated, and resuspended in 100 µL 1X NEBuffer 2.

dNTP fill-in · Ligation · Exonuclease digestion

The protocol fills sticky ends with labeled nucleotides, ligates chromatin aggregates under dilute conditions, and removes biotinylated free ends by exonuclease digestion. Subsequent proteinase K treatment and DNA recovery prepare the ligated DNA for cleanup and elution.

  • dNTP fill-in: Unbound streptavidin sites were blocked with free biotin before the fill-in reaction to prevent labeled dNTP capture.Beads were incubated with free biotin for 15 minutes at RT and then washed twice.
  • dNTP fill-in: Sticky ends were filled in with a-S-dGTP, biotinylated dCTP, other dNTPs, and 25 U Klenow at 25°C for 40 minutes.The reaction volume was 165 µL and was stopped with 7 µL of 0.5M EDTA.
  • dNTP fill-in: After fill-in, beads were washed twice in pre-ligation wash buffer and resuspended in 100µL PLWB.PLWB contained 50mM Tris 7.4, 0.4% Triton X-100, and 0.1mM EDTA.
  • Ligation: Ligation was performed in at least 1mL T4 ligation buffer at 16°C for a minimum of 4 hours.The large volume minimized cross-ligation between different chromatin aggregates.
  • Ligation: The ligation reaction was stopped with 40µL of 0.5M EDTA before beads were concentrated and resuspended in 100µL extraction buffer.The extraction buffer contained 50mM Tris-Cl pH 8.0, 1mM EDTA, and 0.2% SDS.
  • Ligation: DNA underwent overnight digestion with 400ug Proteinase K at 55°C, followed by 2 hours with an additional 200 ug at 55°C.DNA was recovered using SPRI beads at a 2:1 ratio, a column purification kit, or phenol:chloroform extraction.
  • Ligation: Recovered DNA was eluted into Low TE containing 10 mM Tris-Cl pH 8.0 and 0.5 mM EDTA.This followed the available SPRI, column, or phenol:chloroform recovery methods.
  • Exonuclease digestion: 100 U Exonuclease III removed biotinylated free ends during a 40-minute digestion at 37°C, followed by SPRI cleanup and elution into 101 µL low TE.The digestion used Exonuclease III before the final cleanup step.

Shearing and Library Prep · Read mapping · De novo assemblies

The workflow sheared and prepared DNA for sequencing, applied junction-aware read mapping with duplicate removal and quality filtering, and generated human and alligator assemblies from processed short-read datasets.

  • Shearing and Library Prep: DNA was sheared with a Bioruptor for 60 cycles of 30 seconds on and 30 seconds off.The sheared DNA was then filled in with Klenow polymerase and T4 PNK at 20°C for 30 minutes.
  • Shearing and Library Prep: Filled-in DNA was captured on washed C1 beads, incubated for 20 minutes with rocking, washed three times, and resuspended in Low TE.Unbiotinylated fragments were removed before sequencing libraries were generated using established protocols.
  • Read mapping: Sequence reads were truncated at MboI or MluCI junctions before forward and reverse reads were independently aligned with SMALT.Junction sequences were GATCGATC for MboI and AAT-TAATT for MluCI.
  • Read mapping: PCR duplicates were marked with Picard-tools, and analysis retained non-duplicate pairs whose reads both mapped with mapping quality greater than 10.Both reads in each retained pair had to map successfully.
  • De novo assemblies: Human and alligator de novo shotgun assemblies were generated with Meraculous 2.0.3 using publicly available short-insert and mate-pair reads.The alligator mate-pair reads were adapter-trimmed with Trimmomatic.
  • De novo assemblies: Overlapping alligator short-insert reads that had been merged were unmerged back into forward and reverse reads.This processing followed adapter trimming of the alligator mate-pair reads.

Chicago HighRise (HiRISE) Scaffolder · Input pre-processing · Estimation of likelihood model parameters

HiRise pre-processes Chicago data by masking repetitive intervals and excluding highly connected links, then models Chicago-pair likelihoods to guide assembly decisions and estimate model parameters. The likelihood model incorporates contig orientations, gaps, genomic separations, and noise pairs.

  • Input pre-processing: HiRise uses aligned shotgun-read depth to identify repetitive genomic intervals likely to produce misleading Chicago links.A double-threshold strategy flags intervals with elevated overall depth containing at least one exceptionally high-depth base.
  • Input pre-processing: About 0.5% of the assembly was masked using thresholds t1 and t2 selected for practical filtering.The thresholds identify intervals exceeding t1 that contain at least one base exceeding t2.
  • Input pre-processing: HiRise excludes Chicago links within 1 Kbp windows connected to more than four input contigs by at least two links.This removes links from locally highly connected regions during input pre-processing.
  • Estimation of likelihood model parameters: HiRise uses a Chicago-data likelihood model to guide assembly decisions and optimize contig order and orientation within scaffolds.The likelihood describes spanning-pair counts and implied separations between contigs with orientations o ∈ ++, +−, −+, −− and gap length g.
  • Estimation of likelihood model parameters: The separation distribution f(x) combines noise pairs sampling the genome independently with a modeled distribution f′(x).The model is expressed as f(x) = pn/G + (1 −pn)f′(x), with f′(x) represented as a sum of exponentials.
  • Estimation of likelihood model parameters: For limited-contiguity assemblies, HiRise first fixes an estimate of Npn, representing the total number of noise pairs.This provides a robust starting point for estimating N, pn, G, and f′(x).
  • Estimation of likelihood model parameters: HiRise estimates noise-pair parameters from link densities between sampled contig pairs after excluding the highest and lowest 1% of densities.It sets G to the sum of input-contig lengths and uses the density-based estimate for Nn.
  • Estimation of likelihood model parameters: Remaining parameters in Nf(x) are fit by least squares to observed Chicago-pair separation histograms after a multiplicative correction involving G.The correction is applied to smoothed counts at separation x.

Break low-support joins in the input contigs · Contig-contig linking graph construction

The assembly pipeline identifies low-support regions using a likelihood-ratio model, breaks qualifying contig segments, and represents the resulting contigs in a Chicago-linking graph. Scaffolding partitions this graph by link-count thresholds for parallel processing.

  • Break low-support joins in the input contigs: The likelihood model evaluates the log likelihood change from joining the left and right sides at each position i in every starting-assembly contig.The support is expressed as the log likelihood ratio L_i = ln L(g = 0)/L(g = ∞).
  • Break low-support joins in the input contigs: A maximal internal segment is defined as low support when its support falls below threshold t_b.The threshold is applied to the support for the two contigs that would result from breaking at position i.
  • Break low-support joins in the input contigs: Low-support segments within 300bp of one another are merged before candidate breaks are selected.Segments within 1 Kbp of a contig end are excluded from this process.
  • Break low-support joins in the input contigs: Segments are broken at their midpoint, unless longer than 1000 bp, in which case breaks are introduced at both ends.The procedure applies only after merging nearby low-support segments and excluding segments near contig ends.
  • Contig-contig linking graph construction: The Chicago linking data forms a graph whose nodes are the broken starting-assembly contigs.Each edge is labeled with ordered integer pairs representing the positions in the two contigs of reads from a mapped Chicago pair.
  • Contig-contig linking graph construction: Initial scaffolding steps run in parallel on data subsets created by partitioning the graph into connected components.The partitioning excludes edges with fewer than t_L Chicago links.
  • Contig-contig linking graph construction: The lowest integer threshold t_L is chosen so that no connected component comprises more than 5% of the input contigs.This threshold controls which Chicago-link edges remain during graph partitioning.

Seed scaffold construction · Edge filtering

HiRise seeded iterative scaffold construction by filtering the contig-contig graph, extracting high-confidence linear subgraphs, and selecting maximum-likelihood contig orientations. Edge filtering excluded promiscuous contigs using graph-connectivity thresholds calibrated to remove approximately 5% of an upper-tail distribution.

  • Seed scaffold construction: The iterative scaffold-construction phase began by filtering contig-contig graph edges and decomposing the graph into high-confidence linear subgraphs.The filtered graph’s minimum spanning forest was found before linearization.
  • Seed scaffold construction: The filtered contig-contig graph was linearized through three successive rounds removing degree-1 nodes, followed by removal of nodes with degree greater than 2.These operations produced a graph whose connected components had linear topology.
  • Seed scaffold construction: Each resulting connected component defined an ordering of a subset of the input contigs.The components had a linear topology after node removal.
  • Seed scaffold construction: Initial scaffolds were completed by selecting the maximum-likelihood contig orientations for each linear component.Orientation choice was the final step in creating the initial scaffolds.
  • Edge filtering: Before linearization, edges from “promiscuous” contigs were excluded from the contig-contig graph.The edge-filtering procedure was applied before graph linearization.
  • Edge filtering: “Promiscuous” contigs were identified when graph degree divided by contig length in basepairs exceeded tp.The criterion used the degree of the corresponding graph node and the contig length in basepairs.
  • Edge filtering: Edges were also excluded for contigs having links with at least tL links to more than dm other contigs.This was an additional promiscuity criterion in the edge filter.
  • Edge filtering: Thresholds tp and dm were selected to exclude approximately 5% of the upper tail of the relevant distribution.The thresholds governed the promiscuous-contig filtering criteria.

Contig orienting · Merge scaffolds within components

The method determines scaffold orientations by dynamic programming, with links to scaffolds several steps back improving orientation accuracy. Within connected components, candidate joins, insertions, and inversions are evaluated by their total LLR-score change, and score-increasing moves are accepted.

  • Contig orienting: Each input scaffold is assigned a forward or reverse orientation corresponding to the Watson or Crick DNA strand.The optimal orientations are determined separately for each linear scaffold string.
  • Contig orienting: Dynamic programming finds the highest-scoring sequence of scaffold orientation choices in an ordered list.The recurrence tracks orientations of scaffolds up to position i, including a window of preceding scaffolds.
  • Contig orienting: Including links from scaffolds k steps back significantly improves orientation accuracy.Intercalated scaffolds may link on only one side, while links that jump over them provide orientation information for flanking scaffolds.
  • Merge scaffolds within components: Contig ends are classified as free at scaffold termini or buried when internal to a scaffold.For every pair of contig ends within each connected component, a joining LLR is computed using the standard gap size g0.
  • Merge scaffolds within components: Candidate operations are sorted by decreasing score and tested as end-to-end joins, scaffold insertions, or inversions according to end locations.The tested configurations include joining different scaffolds, inserting a scaffold beside a buried end, inverting a same-scaffold interval, and testing four joins for buried ends on different scaffolds.
  • Merge scaffolds within components: A candidate move is accepted when it increases the total LLR score.The total change is computed by summing LLR scores between all affected contig pairs, and the best move is accepted if the change is positive.

Local order and orientation refinement · Iterative joining

The method refines contig order and orientation within scaffolds using a windowed dynamic program, then iteratively joins scaffolds by likelihood-scored end-to-end and intercalating candidates. The refinement considers all local ordering and orientation configurations but is practically limited to small window sizes because w!2w grows steeply.

  • Local order and orientation refinement: A dynamic-programming algorithm slides a window of size w across the ordered and oriented contigs in each scaffold.It refines both local ordering and orientations.
  • Local order and orientation refinement: At each position i, the algorithm considers all w!2w orderings and orientations of contigs within the window.These configurations represent the local order-and-orientation alternatives evaluated by the method.
  • Local order and orientation refinement: For each window position, it stores a score for the optimal ordering and orientation of contigs up to the window’s end.The stored score ends with the current ordering and orientation of the window’s contigs.
  • Local order and orientation refinement: Scores from compatible orderings and orientations at positions i −1 through i −w are used to score extensions with the current configuration.The recurrence compares compatible prior-window configurations when extending the current ordering.
  • Local order and orientation refinement: Because w!2w is a steep function, the refinement method is limited in practice to small values of w.The combinatorial number of configurations constrains the usable window size.
  • Iterative joining: After initial scaffolding within each connected component, the resulting scaffolds are pooled for multiple rounds of end-to-end and intercalating joins.The iterative procedure returns component-level scaffolds to a single pool before joining.
  • Iterative joining: In each round, all scaffold pairs receive parallel likelihood scores for end-to-end and intercalating joins, after which candidate joins are sorted and non-conflicting joins accepted by decreasing likelihood.Acceptance proceeds in decreasing order of likelihood score.

Competing financial interests

The authors report patent applications covering the manuscript’s technology and disclose commercial and advisory ties to Dovetail Genomics.

  • Competing financial interests: The authors have applied for patents on the technology described in the manuscript.Dovetail Genomics LLC was established to commercialize this technology.
  • Competing financial interests: R.E.G. is Founder and Chief Scientific Officer of Dovetail Genomics, while D.H. and D.S.R. are members of its Scientific Advisory Board.
Loading 1502.05331v1…