Information Theory of DNA Shotgun Sequencing

1.1 Background and Motivation

• The information-theoretic framing treats the genome as a random string and defines the minimum reads for reliable reconstruction as an algorithm-independent benchmark.Coverage alone is not generally sufficient: reads must also be stitchable, which depends on DNA repeats and read noise.Source: §1.1 · pp. 3–4

1.2 Main contributions

• The basic model uses an i.i.d. DNA sequence, fixed-length reads at independent uniform locations, and a noiseless read process.The sequence has length G, and symbols are drawn from a distribution p on {A, C, G, T}.Source: §1.2 · p. 4
• Above the threshold 2/H2(p), coverage is sufficient for reconstruction; below it, reconstruction is impossible in the asymptotic regime.The threshold separates a repeat-limited regime, with many length-L repeats, from a coverage-limited regime, with no length-L repeats, with high probability.Source: §1.2 · pp. 5–6
• For Markov DNA, the threshold result replaces H2(p) with the Renyi entropy rate; noisy reads raise the read-length threshold for impossibility.The paper analyzes a modification of the greedy algorithm for noisy reads, while noting that long repeats and nonuniform read sampling remain important limitations.Source: §1.2 · p. 6

1.3 Related Work

• For equiprobable DNA, reconstruction is possible when L > log G with G/L ln(G/Lǫ) reads, versus Li’s L > 4 log G and O(G/L ln G) reads.The authors say their characterization is necessary and sufficient and gives the correct pre-log constant on the read count.Source: §1.3 · p. 6
• Direct DNA reconstruction is not equivalent to shortest-common-superstring assembly: the shortest superstring is not guaranteed to recover the original sequence and may merge repeats, especially in repeat-rich sequences [15].For, reads are independent rather than drawn from one mother sequence; in this paper’s model, overlapping reads are highly correlated.Source: §1.3 · pp. 7–8

1.4 Notations and Outline

• The paper proceeds from problem formulation and impossibility results to an optimal algorithm, then treats noise, richer DNA statistics, and future work.Formal proofs are placed in the appendices.Source: §1.4 · p. 8

2.1 Formulation and Result

• The DNA is modeled as a length-G i.i.d. sequence over {A, C, G, T} with symbol distribution p, and is treated as circular to avoid boundary effects.The results also apply to non-circular sequences with minor modifications.Source: §2.1 · p. 9
• The task is perfect reconstruction from length-L reads whose starting locations are independent and uniformly sampled from the DNA sequence.The section studies the asymptotic ratio of the minimum reads required to the number required for coverage as L and G grow with L/log G constant.Source: §2.1 · pp. 9–10
• For L_bar < 2/H_2(p), reconstruction is impossible; for L_bar > 2/H_2(p), c_min(L_bar)=1 and a greedy algorithm is asymptotically optimal.Here H_2(p) is the order-2 Renyi entropy; c_min is the minimum normalized coverage depth.Source: §2.1 · p. 10

H2(p): Repeat-limited regime

• The two obstructing patterns are a three-way repeat and two interleaved repeat pairs of length L−1; Arratia et al. identify the latter as typical of reconstruction failure.Figure 4 illustrates how interleaved repeats make the order of sequences x and y ambiguous.Source: H2(p): Repeat-limited regime · p. 11, Figure 4 · p. 11
• The heuristic places the repeat-count transition at L̄ = 2/H2(p): expected repeats become numerous below it, and approximately uniform positions make interleaving likely.The probability that two physically disjoint subsequences match is governed by the order-2 Rényi entropy.Source: H2(p): Repeat-limited regime · pp. 11–12
• Swapping the intervening sequences leaves every read’s occurrence count and likelihood unchanged, so MAP reconstruction error is at least 1/2 when such repeats occur.This holds regardless of the number of reads.Source: H2(p): Repeat-limited regime · p. 12
• For L̄ < 2/H2(p), interleaved repeats occur with probability approaching one, forcing Nmin(ε, G, L) = ∞ for ε < 1/2 regardless of read count.The repeats occur when L = L̄ log G.Source: H2(p): Repeat-limited regime · pp. 11–12

H2(p): Coverage-limited regime

• For L̄ > 2/H2(p), greedy overlap-based merging reaches the coverage lower bound when reads cover the sequence and no length-L repeats occur.Coverage is necessary because missing nucleotides—and gaps that split the sequence into ambiguously ordered pieces—prevent reconstruction; the greedy method merges contigs by overlap score.Source: H2(p): Coverage-limited regime · pp. 13–15

2.4 Performance of Existing Algorithms

• The paper evaluates normalized coverage depth against the information-theoretic limit, not to rank assemblers designed for noisier reads and repeats.Greedy merging was used in Sanger assemblers including phrap, TIGR Assembler, and CAP3.Source: §2.4 · p. 16

2.4.1 Sequential Algorithm

• L̄H2(p)^−1 is the sequential algorithm’s minimum normalized coverage depth when L̄ > 2/H2(p), and its performance is strictly worse than greedy.This strategy grows one contig sequentially by repeatedly extending it with the read having the largest overlap.Source: §2.4.1 · p. 16

2.4.2 K-mer based Algorithms

• For L̄ > 2/H2(p), K-mer-based algorithms perform strictly worse than greedy because they require successive reads to have normalized overlap of at least 2/H2(p), whereas greedy only requires coverage.K is chosen so that, with high probability, K-mers from physically disjoint DNA regions are distinct.Source: §2.4.2 · pp. 18–19

2.5 Complexity of the Greedy Algorithm

• A sort-and-extend implementation reduces greedy assembly from O(N^2) pairwise comparisons to O(LN log N); Ukkonen’s data structure achieves O(NL).The implementation extracts ℓ-mers at each overlap stage; Turner described a similar-complexity idea using a suffix tree.Source: §2.5 · p. 19

3.1 Formulation and Result

• The model uses a first-order Markov source with transition matrix Q, and the result generalizes to any fixed order m.The order m must not grow with sequence length G.Source: §3.1 · p. 20
• The minimum normalized coverage depth is 1 when normalized read length exceeds 2/H2(Q).The result extends the i.i.d. model’s dependence on order-2 Rényi entropy to Markov sources.Source: §3.1 · p. 20

3.2 Sketch of Proof

• Below the threshold 2/H2(Q), the proof establishes that Markov DNA contains interleaved repeats.This is the lower-threshold argument in the proof sketch.Source: §3.2 · pp. 20–21
• Above the threshold, the greedy algorithm with the i.i.d.-case overlap score achieves minimum normalized coverage depth 1, unaffected asymptotically by physically overlapping reads.This is the upper-threshold argument.Source: §3.2 · p. 21

4.1 Formulation and Result

• The algorithm scores overlaps by applying a maximum a posteriori test to decide whether read suffixes and prefixes came from the same physical source subsequence.The read model perturbs each nucleotide independently according to Q(r|s), and noise is independent across reads.Source: §4.1 · pp. 22–23
• Perfect layout means mapping every read to its true location, but does not guarantee perfect reconstruction; coverage yields correctness at most positions.Under coverage, most positions receive O(log G) reads.Source: §4.1 · p. 22
• The model treats sequencing noise as a discrete memoryless channel, while the optimal assembly algorithm for noisy reads remains an open problem.The passage notes actual technologies may also have increasing noise, insertions, or deletions.Source: §4.1 · p. 22
• Theorem 10 states that the modified greedy algorithm can achieve normalized coverage depth c(¯L) = 1; the supplied excerpt does not specify a threshold condition.Source: §4.1 · p. 23

4.2 Sketch of Proof for Uniform Source and Symmetric Noise

• For a uniform source with symmetric noise, full normalized coverage is achievable when ¯L > 2/I∗(ǫ), a threshold that increases continuously with error rate ǫ.At ǫ = 0, the threshold is 2/H2(p) = 1.Source: §4.2 · pp. 23–24, and α∗satisfies · pp. 24–25, Figure 9 · p. 25

and α∗satisfies

• The proof sketch uses a Hamming-distance overlap cutoff α to balance false overlaps between unrelated reads against missed overlaps caused by excessive noise.Optimizing α equates the two quantities in the bound; the corollary follows by applying Theorem 10.Source: and α∗satisfies · pp. 24–26

Long repeats

• Because the i.i.d. and Markov models capture only short-range correlations, they miss long repeats that bottleneck real-genome assembly; follow-up work [3] extends analysis to arbitrary repeat statistics.Work [3] generalizes reconstruction bounds using repeat statistics and develops a de Bruijn graph assembly pipeline that computes feasibility plots for read length and coverage depth.Source: §5 · p. 27

Double-strandedness

• For i.i.d. and Markov models, sampling both DNA strands leaves asymptotic minimum normalized coverage unchanged; greedy assembly includes reverse complements and stops at the two reconstructed strands.The error probability changes by only a factor of two, and the proof shows that overlapping reads—including reverse-complement overlaps—make a negligible contribution.Source: §5 · pp. 27–28

Read process

• More realistic read models could capture indels, quality scores, correlated or nonuniform errors, sequence-dependent error, biased sampling, and mate-pairs that help resolve long repeats.The paper notes that noise and sampling bias make assembly more difficult; paired-end reads have an approximately known separation.Source: §5 · p. 28

Partial reconstruction

• When perfect reconstruction is infeasible, partial assembly requires defining an optimal tradeoff between the number of contigs and the number of misjoins.The output must allow multiple contigs as well as errors.Source: §5 · p. 28

A Proof of Theorem 1, Part 2

• Although the lemma can be found in, the authors include its proof to ease generalization to the Markov case.Source: A Proof of Theorem 1, Part 2 · p. 29
• Lemma 12 concerns any distinct substrings of length ℓ from the i.i.d. DNA sequence.Source: A Proof of Theorem 1, Part 2 · p. 29
• For distinct length-ℓ substrings with no physical overlap, the probability of identity is e^(-ℓH_2(p)).This is the stated result of Lemma 12.Source: A Proof of Theorem 1, Part 2 · p. 29

• For overlapping strings sharing k bases, equality requires the longer sequence to have a repeating UV pattern, with component lengths determined by dividing ℓ by ℓ−k.The proof bounds the matching probability using the probability of observing this constrained structure; it also notes that overlapping reads can have greater confusion probability, so both overlapping and non-overlapping cases matter.Source: §2 · pp. 29, 30, 33
• If every read is correctly merged to its successor at the correct physical offset and reads cover the sequence, the greedy algorithm outputs the DNA sequence.The error event is contained in the union of incorrect merging and lack of coverage.Source: §2 · p. 30
• With the chosen N, P(E1) and P(E2) vanish, so reconstruction error vanishes and the minimum-read requirement is bounded.E1 is the event of an incorrect merge; E2 is the event that the DNA sequence is not covered by reads.Source: §2 · p. 34

B Proof of Theorem 6

• The source’s stationary distribution is denoted p = (p1, p2, p3, p4)t.Source: B Proof of Theorem 6 · p. 34
• Because Q̄ has positive entries, its largest eigenvalue is real and positive, and its corresponding eigenvector π has positive components.The passage then introduces a useful inequality without specifying it.Source: B Proof of Theorem 6 · pp. 34–35

B.1 Proof of Lemma 7

• The proof follows Arratia et al.’s Poisson-approximation approach to bounding recoverability from the L-spectrum, using the theorem of Arratia et al..Arratia et al. identify interleaved repeat pairs as the dominant cause of non-recoverability.Source: §B.1 · p. 35
• Given k repeats, the probability of having no interleaved repeat pairs is calculated using Catalan numbers.The proof then approximates this event through the repeat-count indicators and a Chen-Stein bound.Source: §B.1 · pp. 36–37
• The proof partitions the sequence into length-L blocks and counts repeats that cannot extend left, avoiding clumping that would undermine the Poisson approximation.The indicators track pairs of blocks matching over positions 2 through L while differing at the first position.Source: §B.1 · p. 36
• Because the right-hand side of (32) approaches zero, with probability 1 − o(1) the sequence contains two pairs of interleaved repeats.Source: §B.1 · p. 38

1. If the strings have no physical overlap, the probability that they are identical is bounded

• The supplied passages provide no details of the no-physical-overlap case or its bound.This limitation is specific to the excerpt; it does not establish that the paper omits that proof.Source: §B.2 · p. 39

B.2 Proof of Lemma 8

• The proof follows the i.i.d. model argument, replacing Lemma 12 with Lemma 14.Source: §B.2 · p. 39
• Lemma 14 concerns any distinct length-ℓ substrings of the Markov DNA sequence.Source: §B.2 · p. 39

2. If the strings have physical overlap, the probability that they are identical is bounded

• For overlapping strings, the proof assumes X = S[1, ℓ] and Y = S[ℓ−k + 1, 2ℓ−k] for k in {1, …, ℓ−1}.Source: §2 · p. 39
• It decomposes S[1, 2ℓ−k] into q blocks of length ℓ−k and a remainder V of length r.Here q and r are the quotient and remainder when 2ℓ−k is divided by ℓ−k.Source: §2 · p. 39
• Under that decomposition, X = Y exactly when every block U_i and the remainder V match the corresponding prefix of S.The excerpt specifies the matching prefixes as S[1, ℓ−k] for each U_i and S[1, r] for V.Source: §2 · pp. 39–40
• The overlap bound’s displayed inequality is absent, but its proof invokes Cauchy–Schwarz, ∑P_i^2 ≤ 1, (25), and k < ℓ with ρ_max( Q̄) ≤ 1.Source: §2 · p. 40

C Proof of Theorem 10

• The overlap score is the longest suffix-prefix match passing MAP rule (11), whose null hypothesis says the reads share a physical source subsequence.The null hypothesis H0 represents reads from the same physical source subsequence.Source: C Proof of Theorem 10 · p. 40

H0 and H1 and their corresponding MAP rule (11),

• Lemma 16 bounds false acceptance of H0 for distinct reads both when they lack physical overlap and when they overlap.The supplied passages identify these two cases but do not provide the bound expressions.Source: H0 and H1 and their corresponding MAP… · p. 41
• For physically overlapping reads, the false-alarm proof partitions dependent χj terms into equal-sized sets whose terms are independent, enabling use of Theorem 15.The proof constructs the partition so the underlying pairs of sequence symbols within each set are disjoint.Source: H0 and H1 and their corresponding MAP… · pp. 41–42
• Lemma 17 bounds mis-detection: for distinct reads from the same physical location, it bounds the probability that H1 is accepted.This addresses the possibility that noise prevents physically overlapping reads from passing the overlap criterion.Source: H0 and H1 and their corresponding MAP… · p. 42
• Reconstruction error is decomposed into erroneous merging and incomplete coverage, with the latter vanishing exponentially when R > L̄.The erroneous-merge event is further split by the stage at which the first merging error occurs.Source: H0 and H1 and their corresponding MAP… · p. 42
• The proof shows both error events vanish, then optimizes over θ to obtain Theorem 10’s result.It establishes P(E1) = o(1) for erroneous merges and P(E2) = o(1) for incomplete coverage for any θ before optimizing.Source: H0 and H1 and their corresponding MAP… · p. 44