---
title: "Genomics and Sequencing"
book: "University Biology — Year 3"
subject: biology
language: en
chapter: 4
exercises: 12
source: https://one-course.com/books/biology/5/en/chapter/4-genomics-and-sequencing
---

# Chapter 4 — Genomics and Sequencing

The first complete genome of a free-living organism, the $1.8\,\mathrm{Mb}$ of a bacterium, was published in 1995 after a year’s work by a team of forty. The human genome, $3.2$ billion base pairs, took an international consortium thirteen years and about three billion dollars and was declared finished in 2003. Today a bench-top machine [reads](#def-b3-genomics-ngs) a human genome overnight for a few hundred dollars, a hospital sequences the tumour of a patient to choose a drug, and a museum extracts and [reads](#def-b3-genomics-ngs) the genome of a forty-thousand-year-old bone. The technology that made this possible, the mathematics that turns millions of short [reads](#def-b3-genomics-ngs) into a genome, and what the finished genomes have taught about the size, content and history of our DNA are the subject of this chapter. The chapter after it takes up the algorithms that compare the sequences once they are in hand.

## 4.1 Reading DNA

**Definition 4.1 (Sequencing by chain termination).**

*Sanger sequencing* copies a single-stranded template from a primer with DNA polymerase in the presence of the four normal nucleotides and a small proportion of *dideoxynucleotides*, which lack the 3$'$ hydroxyl and so terminate the chain wherever they are incorporated. Each dideoxynucleotide carries a different fluorescent dye. The product is a mixture of fragments, one for every position of the template, each ending in a base of known identity; separated by size in a capillary of gel, they pass a detector in order of length, and the sequence of colours is the sequence of the template. A run [reads](#def-b3-genomics-ngs) $700\text{ to }900\,\mathrm{bp}$ with an error rate below $10^{-3}$; it remains the method for verifying a construct or a single gene.

**Evidence.** Sanger, Nicklen and Coulson published the method in 1977 and used it that year to [read](#def-b3-genomics-ngs) the $5386\,\mathrm{bp}$ of the phage $\phi$X174 — the first complete DNA genome — and in 1981 the $16\,569\,\mathrm{bp}$ of the human mitochondrion. Fleischmann and colleagues (1995) [read](#def-b3-genomics-ngs) the $1.83\,\mathrm{Mb}$ of *Haemophilus influenzae* by breaking the whole genome into random fragments, sequencing $24\,000$ of them and assembling the [reads](#def-b3-genomics-ngs) by computer — the *whole-genome shotgun* strategy that every later project scaled up. ∎

![Left: Frederick Sanger, who devised the first practical methods for reading proteins and DNA and received a Nobel prize for each. Right: an autoradiograph of a sequencing gel from the pre-capillary era, one lane per base, the sequence read upward from the ladder of bands.](https://one-course.com/images/onecourse/chapters/biology-5/b3-genomics/img-449ad1d73146.jpg)

![Left: Frederick Sanger, who devised the first practical methods for reading proteins and DNA and received a Nobel prize for each. Right: an autoradiograph of a sequencing gel from the pre-capillary era, one lane per base, the sequence read upward from the ladder of bands.](https://one-course.com/images/onecourse/chapters/biology-5/b3-genomics/img-874706cdd095.jpg)

*Left: Frederick Sanger, who devised the first practical methods for reading proteins and DNA and received a Nobel prize for each. Right: an autoradiograph of a sequencing gel from the pre-capillary era, one lane per base, the sequence [read](#def-b3-genomics-ngs) upward from the ladder of bands.*

![Chain-termination sequencing. A dideoxy base ends the copy at every position where it is incorporated; the fragments, one per length, are separated by size and the colour of each terminal base is read in order.](https://one-course.com/images/onecourse/chapters/biology-5/b3-genomics/fig-b39964002efe.svg)

*Chain-termination sequencing. A dideoxy base ends the copy at every position where it is incorporated; the fragments, one per length, are separated by size and the colour of each terminal base is [read](#def-b3-genomics-ngs) in order.*

**Definition 4.2 (Massively parallel sequencing).**

Second-generation instruments read millions to billions of fragments at once. In *sequencing by synthesis* the fragments, with adapters ligated to their ends, are bound to a glass flow cell and amplified in place into clusters of identical molecules; the clusters are then extended one base per cycle with fluorescent, reversibly blocked nucleotides, imaged, unblocked, and extended again, so that each cycle adds one base to every cluster’s *read*. Reads are $100\text{ to }300\,\mathrm{bp}$, usually from both ends of a fragment (*paired ends*), with an error rate of about $10^{-3}$ per base, and one run yields up to $10^{12}$ bases. Third-generation *long-read* instruments read single molecules without amplification: by watching one polymerase incorporate fluorescent nucleotides in real time, or by threading the DNA through a protein *nanopore* and recording the ionic current, which each sequence of bases modulates in its own way. Reads of $10\text{ to }100\,\mathrm{kb}$ and more span the repeats that short reads cannot, at a higher raw error rate that consensus reduces.

![Left: a flow cell, the glass slide on which billions of DNA clusters are grown and read one base per cycle. Right: a pocket-sized nanopore sequencer reading single molecules as changes in an ionic current.](https://one-course.com/images/onecourse/chapters/biology-5/b3-genomics/img-1248e19cab5b.jpg)

![Left: a flow cell, the glass slide on which billions of DNA clusters are grown and read one base per cycle. Right: a pocket-sized nanopore sequencer reading single molecules as changes in an ionic current.](https://one-course.com/images/onecourse/chapters/biology-5/b3-genomics/img-7e84822a8975.jpg)

*Left: a flow cell, the glass slide on which billions of DNA clusters are grown and [read](#def-b3-genomics-ngs) one base per cycle. Right: a pocket-sized [nanopore](#def-b3-genomics-ngs) sequencer reading single molecules as changes in an ionic current.*

**Theorem 4.3 (Coverage and gaps in a shotgun project).**

Let $N$ [reads](#def-b3-genomics-ngs) of length $L$ be taken at random positions along a genome of length $G$, and let $c = NL/G$ be the *coverage*, the average number of [reads](#def-b3-genomics-ngs) covering a base. Then a given base is covered by a number of [reads](#def-b3-genomics-ngs) that is Poisson with mean $c$: the fraction of the genome left unsequenced is

$$
P(\text{uncovered}) = e^{-c},
$$

and if two [reads](#def-b3-genomics-ngs) are recognised as overlapping only when they share at least $T$ bases, so that $\theta = T/L$, the expected number of *contigs* (islands of overlapping [reads](#def-b3-genomics-ngs)) is

$$
\text{contigs} = N\,e^{-c(1-\theta)},
$$

and their mean length is $L\,\bigl(e^{c(1-\theta)} - 1\bigr)/c + L\theta$, approximately.

**Proof.** [Read](#def-b3-genomics-ngs) start points fall at random with density $N/G$ per base. A base is covered by the [reads](#def-b3-genomics-ngs) that start in the $L$ bases before it; the number of starts in a window of length $L$ is Poisson with mean $LN/G
= c$, and the probability of none is $e^{-c}$. A [read](#def-b3-genomics-ngs) is the rightmost of its contig if no other [read](#def-b3-genomics-ngs) starts within the $L - T = L(1-\theta)$ bases after its own start (a [read](#def-b3-genomics-ngs) starting later would overlap it by less than $T$ and not be joined); that probability is $e^{-c(1-\theta)}$, and since each contig has exactly one rightmost [read](#def-b3-genomics-ngs), the expected number of contigs is $N e^{-c(1-\theta)}$. The mean contig length follows from $G$ divided by the number of contigs, corrected for the uncovered fraction. ∎

**Example 4.4 (How much is enough).**

At $c = 5$ the unsequenced fraction is $e^{-5} = 0.7\,\%$ — for a $3.2\,\mathrm{Gb}$ genome, twenty million bases in some tens of thousands of gaps. At $c = 10$ it is $4.5\times 10^{-5}$, $150\,\mathrm{kb}$ in total. Human genomes are routinely sequenced at $c = 30$ not because of coverage gaps ($e^{-30} \approx 10^{-13}$) but because each base must be [read](#def-b3-genomics-ngs) several times on each of the two chromosomes to call a heterozygous variant with confidence against an error rate of $10^{-3}$ per [read](#def-b3-genomics-ngs). Bacterial genomes are sequenced at $c = 50$–$100$ for the same reason and because it is cheap. The formula also shows what coverage cannot fix: a repeat longer than a [read](#def-b3-genomics-ngs) is a place where the overlap graph branches, and no amount of short [reads](#def-b3-genomics-ngs) resolves it. That is what long [reads](#def-b3-genomics-ngs) are for.

![The two Lander–Waterman quantities against coverage: the fraction of the genome never read falls as e-c, and the number of contigs (here scaled per read) peaks at low coverage and then falls as the islands merge.](https://one-course.com/images/onecourse/chapters/biology-5/b3-genomics/fig-2095ebc59ed3.svg)

*The two Lander–Waterman quantities against coverage: the fraction of the genome never [read](#def-b3-genomics-ngs) falls as $e^{-c}$, and the number of contigs (here scaled per [read](#def-b3-genomics-ngs)) peaks at low coverage and then falls as the islands merge.*

## 4.2 From reads to a genome

**Definition 4.5 (Assembly and annotation).**

*Assembly* reconstructs a genome from its [reads](#def-b3-genomics-ngs) by overlap: in an *overlap graph* each [read](#def-b3-genomics-ngs) is a node joined to the [reads](#def-b3-genomics-ngs) it overlaps; in a *de Bruijn graph*, used for billions of short [reads](#def-b3-genomics-ngs), every $k$-mer (subsequence of length $k$) is a node and the genome is a path through them. Both are broken by repeats longer than the [read](#def-b3-genomics-ngs), which give branching paths. Contigs are ordered and oriented into *scaffolds* by paired-end [reads](#def-b3-genomics-ngs) and long-range information (long [reads](#def-b3-genomics-ngs), optical maps, chromosome-contact maps), and the scaffolds are placed on chromosomes. Assembly quality is summarised by the *N50*: the contig length such that half of the assembled bases lie in contigs at least that long. *Annotation* then finds the genes: in bacteria, open reading frames longer than chance; in eukaryotes, by combining sequence signals (splice sites, promoters, codon bias), homology to known proteins, and transcripts sequenced from RNA. The result, for a species, is a *reference genome* to which every later [read](#def-b3-genomics-ngs) from that species is aligned rather than assembled anew.

**Method 4.6 (From a sample to variants).**

For a resequencing study of an individual against a reference: (1) extract DNA, fragment it to $300\text{ to }500\,\mathrm{bp}$, ligate adapters (the *library*); (2) sequence to the coverage required (30$\times$ for a human genome, 100$\times$ for an exome, which captures the $1.5\,\%$ of the genome that codes for protein); (3) align each [read](#def-b3-genomics-ngs) to the reference, tolerating mismatches; (4) at each position count the bases in the [reads](#def-b3-genomics-ngs): a position where about half the [reads](#def-b3-genomics-ngs) disagree with the reference is a heterozygous *variant*, where nearly all disagree a homozygous one, where a few disagree an error; (5) filter by depth, base quality and strand balance; (6) annotate each variant with its effect on any gene (synonymous, missense, nonsense, splice, frameshift) and its frequency in population databases; (7) for a diagnosis, keep the rare variants predicted damaging in genes consistent with the phenotype, and confirm by [Sanger sequencing](#def-b3-genomics-sanger).

## 4.3 What genomes look like

**Proposition 4.7 (Genome size and gene number).**

[Genome sizes](#prop-b3-genomics-sizes) span a factor of $10^{5}$ among eukaryotes — $12\,\mathrm{Mb}$ for yeast, $100\,\mathrm{Mb}$ for the worm, $140\,\mathrm{Mb}$ for the fly, $3.2\,\mathrm{Gb}$ for a human, $16\,\mathrm{Gb}$ for an onion, $130\,\mathrm{Gb}$ for a lungfish, $150\,\mathrm{Gb}$ for the lily *Paris japonica* — while gene numbers span barely a factor of $10$: $6000$ in yeast, $20\,000$ in the worm, $14\,000$ in the fly, about $20\,000$ protein-coding genes in a human, $40\,000$ in rice. This is the *C-value paradox*: [genome size](#prop-b3-genomics-sizes) does not measure complexity, and neither does gene number. What varies is the non-coding content — introns, transposable elements, satellite repeats — and the number of proteins a genome can produce is multiplied far beyond its gene count by alternative splicing ([Chapter 2](https://one-course.com/books/biology/5/en/chapter/2-non-coding-rnas-and-post-transcriptional-regulation#ch-b3-rna-regulation)) and by regulation, which is where the complexity of an organism is mostly written. Bacterial genomes, by contrast, are compact and their size does track gene number: about one gene per kilobase, $88\,\%$ coding in *E. coli*.

**Example 4.8 (The human genome by content).**

Of the $3.2\,\mathrm{Gb}$: protein-coding exons, $1.5\,\%$; introns and untranslated regions of genes, about $35\,\%$; transposable elements and their fossils, about $45\,\%$ — LINE-1 retrotransposons $17\,\%$, Alu elements $10\,\%$ (over a million copies of a $300\,\mathrm{bp}$ sequence), endogenous retroviruses $8\,\%$; segmental duplications $5\,\%$; simple repeats and satellites, including the centromeres, about $5\,\%$; the rest unique non-coding sequence, where the regulatory elements live. Some $100\text{ to }200$ genes still encode active LINE-1 machinery, and new insertions occur about once in every twenty births. About $8\,\%$ of the genome is under detectable purifying selection — far more than the exons — and most of that is regulatory.

![The human genome by content, to scale. The protein-coding exons (red) are a sliver; nearly half of the genome descends from transposable elements.](https://one-course.com/images/onecourse/chapters/biology-5/b3-genomics/fig-dedb3915680f.svg)

*The human genome by content, to scale. The protein-coding exons (red) are a sliver; nearly half of the genome descends from transposable elements.*

**Definition 4.9 (Comparative genomics).**

Genes in two species descended from one gene in their common ancestor are *orthologues*; genes within a species descended from a duplication are *paralogues*. Blocks of chromosome in which gene order is conserved between species are *syntenic*; synteny lets a gene be located in one genome from its position in another and reveals the rearrangements that separate two karyotypes (about $1000$ between human and mouse). Sequences conserved across distant species that code for no protein — *conserved non-coding elements*, some ultraconserved to the base over $200\,\mathrm{bp}$ between human and fish — are mostly enhancers of developmental genes. *Whole-genome duplications* have marked the history of lineages: two rounds at the origin of vertebrates (the four Hox clusters of mammals against the one of invertebrates), one in the ancestor of salmonids, one in the yeast lineage, several in flowering plants; the duplicate genes are mostly lost over tens of millions of years, and the survivors diverge in function.

**Example 4.10 (Human and chimpanzee).**

Aligned single-copy sequence differs between human and chimpanzee by $1.2\,\%$ of bases — some $35$ million substitutions — and by insertions and deletions that together amount to another $3\,\%$ of each genome. Two humans differ at about one base in a thousand, some $4\text{ to }5$ million sites, plus a few thousand structural variants; two chimpanzees, whose population has been larger for longer, at rather more. A child carries about $70$ new mutations absent from both parents, four fifths of them from the father, and the number rises by about two per year of paternal age — the arithmetic of [Chapter 3](https://one-course.com/books/biology/5/en/chapter/3-genome-stability-dna-damage-repair-and-recombination#ch-b3-dna-repair) applied to the many divisions of spermatogenesis.

## 4.4 Genomes and traits

**Definition 4.11 (Genome-wide association).**

A *single-nucleotide polymorphism* (SNP) is a position at which two bases each occur in at least $1\,\%$ of a population; about ten million are common in humans. A *genome-wide association study* (GWAS) genotypes hundreds of thousands of SNPs in thousands of people with a disease and thousands without, and asks at each SNP whether one allele is more frequent in the cases. Because a million tests are made, a result counts only below a significance threshold of about $5\times 10^{-8}$ ($0.05$ divided by the million effectively independent tests). Associated SNPs mark a region, not a causal variant, since neighbouring alleles travel together (linkage disequilibrium) over tens of kilobases. For most common diseases the associated alleles each shift risk by a few percent and lie mostly in regulatory sequence; their combined effect, summed over thousands of SNPs as a *polygenic score*, explains a fraction of the heritability and predicts risk about as well as family history.

**Remark 4.12 (What genomics has and has not delivered).**

Sequencing has been decisive where a single gene has a large effect: several thousand Mendelian diseases have their gene, and sequencing the exome of a child with an undiagnosed disorder yields a diagnosis in about a third of cases. Cancer genomes reveal which drivers a tumour carries and which drugs may work ([Chapter 11](https://one-course.com/books/biology/5/en/chapter/11-cancer-biology#ch-b3-cancer-biology)); pathogen genomes trace outbreaks strain by strain ([Chapter 12](https://one-course.com/books/biology/5/en/chapter/12-bacteriology-growth-physiology-and-genetics#ch-b3-bacteriology)); ancient genomes have rewritten human prehistory, showing that people outside Africa carry about $2\,\%$ Neanderthal DNA. For common diseases — diabetes, heart disease, schizophrenia — genomics has delivered thousands of small effects and few mechanisms, and the promise of prediction from a genome for a healthy person remains modest. The genome is a parts list; the biology of how the parts interact is not [read](#def-b3-genomics-ngs) off it.

## 4.5 Exercises

**Exercise 4.1 ★.**

Explain why a [dideoxynucleotide](#def-b3-genomics-sanger) terminates a growing DNA chain, and why a Sanger reaction must contain both normal and dideoxy forms of each nucleotide.

**Solution of Exercise 4.1.**

A [dideoxynucleotide](#def-b3-genomics-sanger) has no 3$'$ hydroxyl, so no phosphodiester bond can be made to the next nucleotide and the chain stops. With only dideoxy forms every chain would stop at the first position; with only normal forms none would stop. The mixture makes termination a random event at each position, so that the products form a complete ladder, one length for every base of the template.

**Exercise 4.2 ★.**

A run produces $4\times 10^{8}$ paired [reads](#def-b3-genomics-ngs) of $2\times 150\,\mathrm{bp}$. What coverage does that give of a $3.2\,\mathrm{Gb}$ human genome? Of a $5\,\mathrm{Mb}$ bacterial genome?

**Solution of Exercise 4.2.**

$4\times 10^{8}\times 300\,\mathrm{bp} = 1.2\times 10^{11}$ bases. Human: $1.2\times 10^{11}/3.2\times 10^{9} \approx 38\times$. Bacterium: $1.2\times 10^{11}/5\times 10^{6} = 24\,000\times$.

**Exercise 4.3 ★.**

Define contig, [scaffold](#def-b3-genomics-assembly) and [N50](#def-b3-genomics-assembly). An [assembly](#def-b3-genomics-assembly) of $100\,\mathrm{Mb}$ has ten contigs of $8\,\mathrm{Mb}$ and $2000$ of $10\,\mathrm{kb}$. What is its [N50](#def-b3-genomics-assembly)?

**Solution of Exercise 4.3.**

A contig is a contiguous sequence assembled from overlapping [reads](#def-b3-genomics-ngs); a [scaffold](#def-b3-genomics-assembly) is an ordered, oriented set of contigs with gaps of estimated size between them; the [N50](#def-b3-genomics-assembly) is the contig length at which the sorted contigs reach half the assembled bases. Here the ten $8\,\mathrm{Mb}$ contigs alone contain $80\,\mathrm{Mb}$, past the $50\,\mathrm{Mb}$ half-way point (the seventh reaches $56\,\mathrm{Mb}$): [N50](#def-b3-genomics-assembly) $= 8\,\mathrm{Mb}$.

**Exercise 4.4 ★.**

State the C-value paradox with two examples, and say what mostly accounts for the excess DNA of large genomes.

**Solution of Exercise 4.4.**

[Genome size](#prop-b3-genomics-sizes) does not track organismal complexity: an onion has five times the DNA of a human, a lungfish forty times; yeast and the worm differ tenfold in DNA with similar gene numbers. The excess is non-coding: transposable elements and their remains, introns, satellite repeats.

**Exercise 4.5 ★★.**

Using [Theorem 4.3](#thm-b3-genomics-lander-waterman), find the coverage needed to leave at most one base in a million unsequenced, and the expected number of contigs for a $5\,\mathrm{Mb}$ genome [read](#def-b3-genomics-ngs) at $c = 8$ with $150\,\mathrm{bp}$ [reads](#def-b3-genomics-ngs) and a minimum overlap of $30\,\mathrm{bp}$.

**Solution of Exercise 4.5.**

$e^{-c} = 10^{-6}$ gives $c = 6\ln 10 = 13.8$. For $c = 8$: $N =
cG/L = 8\times 5\times 10^{6}/150 \approx 267\,000$ [reads](#def-b3-genomics-ngs), $\theta
= 30/150 = 0.2$, contigs $= N e^{-6.4} = 267\,000\times 0.00166
\approx 440$.

**Exercise 4.6 ★★.**

A heterozygous variant is covered by $30$ [reads](#def-b3-genomics-ngs). Assuming each [read](#def-b3-genomics-ngs) shows either allele with probability $1/2$, what is the probability that fewer than $8$ [reads](#def-b3-genomics-ngs) show the variant allele (so that it might be mistaken for errors)? Use a normal approximation with mean $15$ and standard deviation $\sqrt{7.5}$. Why is $30\times$ the standard?

**Solution of Exercise 4.6.**

$P(X < 8) = P(X \le 7) \approx P\bigl(Z < (7.5 - 15)/2.74\bigr) =
P(Z < -2.74) \approx 0.003$. At $30\times$ both alleles of a heterozygote are almost always seen many times, coverage is uneven (GC-rich regions get fewer [reads](#def-b3-genomics-ngs), so some sites see half the average), and enough [reads](#def-b3-genomics-ngs) remain to separate a real allele from the $10^{-3}$ errors.

**Exercise 4.7 ★★.**

A genome has a repeat of $6\,\mathrm{kb}$ present in $500$ copies. Explain why an [assembly](#def-b3-genomics-assembly) from $150\,\mathrm{bp}$ [reads](#def-b3-genomics-ngs) collapses it, what the resulting graph looks like, and what [read](#def-b3-genomics-ngs) length would resolve it.

**Solution of Exercise 4.7.**

Every $150\,\mathrm{bp}$ [read](#def-b3-genomics-ngs) from inside the repeat is identical whichever copy it came from, so the assembler builds one $6\,\mathrm{kb}$ node with $500$ paths entering and $500$ leaving; it cannot tell which entry pairs with which exit, and the [assembly](#def-b3-genomics-assembly) breaks into $500$ gaps, the repeat present once with a coverage $500$ times the mean. [Reads](#def-b3-genomics-ngs) longer than the repeat plus unique flanks on both sides — $8\,\mathrm{kb}$ or more — span each copy and resolve it.

**Exercise 4.8 ★★.**

Two humans differ at one base in a thousand. How many differences in their exomes ($48\,\mathrm{Mb}$ of coding sequence)? If two thirds of coding differences are synonymous or benign and the rest alter a protein, how many protein-altering variants does one person carry relative to another?

**Solution of Exercise 4.8.**

$4.8\times 10^{7}/1000 = 48\,000$ coding differences; one third, about $16\,000$, alter the protein.

**Exercise 4.9 ★★.**

A GWAS tests $10^{6}$ SNPs at threshold $p < 5\times 10^{-8}$. How many false positives are expected if no SNP is truly associated? A SNP allele raises the risk of a disease with frequency $2\,\%$ to $2.3\,\%$. Explain why such an effect is undetectable in a study of a thousand people and useless for an individual, yet may still point to a mechanism.

**Solution of Exercise 4.9.**

Expected false positives $10^{6}\times 5\times 10^{-8} = 0.05$. A relative risk of $1.15$ on a $2\,\%$ disease changes a few cases per thousand; a thousand people contain about twenty cases, far too few to see it (power grows with the number of cases and the square of the effect). For an individual, $0.3$ percentage points is meaningless. But the allele marks a gene whose modest change in activity alters disease risk — the gene, and its pathway, may be a drug target whose complete inhibition has a large effect.

**Exercise 4.10 ★★★.**

Bacterial genomes are about $88\,\%$ coding; the human genome $1.5\,\%$. Give three hypotheses — population genetic, structural and regulatory — for the difference, and for each a genomic observation that supports or undermines it.

**Solution of Exercise 4.10.**

Population genetic: in bacteria, with enormous populations, selection removes even slightly costly insertions; in mammals, with small effective populations, drift lets mildly deleterious non-coding DNA accumulate — supported by the inverse correlation across lineages between [genome size](#prop-b3-genomics-sizes) and population size, undermined by exceptions. Structural: eukaryotic genomes are invaded by transposons that bacteria, with a bias toward deletion and no meiotic refuge for selfish elements, purge — supported by the correlation of [genome size](#prop-b3-genomics-sizes) with transposon content. Regulatory: complex development needs more regulatory DNA — supported by the [conserved non-coding elements](#def-b3-genomics-comparative) around developmental genes, undermined by the fact that only $8\,\%$ of the human genome shows any selection at all, so most non-coding DNA is not regulatory.

**Exercise 4.11 ★★★.**

Derive the Lander–Waterman contig count for [reads](#def-b3-genomics-ngs) of two different lengths mixed: $N_{1}$ short [reads](#def-b3-genomics-ngs) of length $L_{1}$ and $N_{2}$ long [reads](#def-b3-genomics-ngs) of length $L_{2}$, overlaps counted with the same minimum $T$. (Treat each [read](#def-b3-genomics-ngs) as the rightmost of its contig if no [read](#def-b3-genomics-ngs) of either kind starts within its own length minus $T$.) Show that a few long [reads](#def-b3-genomics-ngs) reduce the contig count more than the same number of bases in short [reads](#def-b3-genomics-ngs).

**Solution of Exercise 4.11.**

Starts of both kinds fall with total density $\rho = (N_{1} +
N_{2})/G$. A [read](#def-b3-genomics-ngs) of length $L_{i}$ is rightmost in its contig if no [read](#def-b3-genomics-ngs) of either kind starts in the $L_{i} - T$ bases after it: probability $e^{-\rho(L_{i} - T)}$. Hence contigs $= N_{1}
e^{-\rho(L_{1} - T)} + N_{2} e^{-\rho(L_{2} - T)}$. A long [read](#def-b3-genomics-ngs) is rightmost with an exponentially smaller probability than a short one, and each long [read](#def-b3-genomics-ngs) also removes the chance of a gap over its whole length; the same number of bases as short [reads](#def-b3-genomics-ngs) adds many start points but each protected window is short, so the long [reads](#def-b3-genomics-ngs) win.

**Exercise 4.12 ★★★.**

The four Hox clusters of mammals are argued to descend from one cluster by two rounds of [whole-genome duplication](#def-b3-genomics-comparative) at the base of the vertebrates. State what pattern of [paralogues](#def-b3-genomics-comparative) across the genome, and what pattern in the genomes of lampreys and of the invertebrate chordates, that hypothesis predicts, and how one would distinguish it from two independent duplications of the cluster alone.

**Solution of Exercise 4.12.**

Two [whole-genome duplications](#def-b3-genomics-comparative) predict that the fourfold pattern is genome-wide: quartets of paralogous chromosome segments (paralogons) carrying the same gene families in the same order, with the duplications dated to the same time for all families; a single cluster in amphioxus and *Ciona*, which diverged before the events; and in lampreys, which diverged around them, an intermediate or independently derived state. Independent duplications of the Hox cluster alone predict paralogons for Hox only, neighbours with different histories, and duplication dates that differ between families. The genome-wide paralogons and the shared dating, as observed, favour [whole-genome duplication](#def-b3-genomics-comparative).

## 4.6 Problem: Two Genomes on One Machine

**Problem 4.1.**

Weekend problem — a bacterium and a human sequenced on the same run: the reads shared out, the coverage and the gaps computed, the contigs counted, the heterozygous variants called against the error rate, the content of the human genome weighed, and an association study sized, ending on the contig count at low coverage, the coverage that makes a variant call safe and the number of new mutations in a child

Data: one run yields $6\times 10^{8}$ [reads](#def-b3-genomics-ngs) of $150\,\mathrm{bp}$. A bacterial genome is $4.6\,\mathrm{Mb}$; the human genome $3.2\,\mathrm{Gb}$ (haploid). Minimum overlap $T = 30\,\mathrm{bp}$. Error rate $10^{-3}$ per base per [read](#def-b3-genomics-ngs). Human coding sequence $48\,\mathrm{Mb}$; the genome is $45\,\%$ transposon-derived, with $1.1$ million Alu copies of $300\,\mathrm{bp}$. Two humans differ at one site in $1000$. Mutation rate $1.2\times 10^{-8}$ per base per generation.

**Part I — The bacterium.**

1. The bacterium is given $0.2\,\%$ of the run. How many [reads](#def-b3-genomics-ngs) , and what coverage?
2. What fraction of its genome is unsequenced? How many bases is that?
3. Compute $\theta$ and the expected number of contigs.
4. An earlier pilot used only $30\,000$ [reads](#def-b3-genomics-ngs) . Coverage, uncovered fraction, and contig number?
5. The genome contains $7$ copies of a $5\,\mathrm{kb}$ rRNA operon. Explain what happens to them in the [assembly](#def-b3-genomics-assembly) and how many contig ends this alone produces.
6. The bacterium has about one gene per kilobase. How many genes, and, if $88\,\%$ of the genome is coding, what is the mean gene length?

**Part II — The human.**

7. The rest of the run goes to one human genome. Coverage of the diploid genome per haploid copy (that is, total bases divided by $3.2\,\mathrm{Gb}$ )?
8. At a heterozygous site each of the two alleles is covered by about half the [reads](#def-b3-genomics-ngs) . With the coverage of question 7, what is the expected number of [reads](#def-b3-genomics-ngs) showing each allele?
9. An error shows a wrong base in a [read](#def-b3-genomics-ngs) with probability $10^{-3}$ . At a homozygous site covered by $c$ [reads](#def-b3-genomics-ngs) , what is the expected number of [reads](#def-b3-genomics-ngs) showing any particular wrong base ( $10^{-3}/3$ each)? Why does a rule “call a variant if at least $3$ [reads](#def-b3-genomics-ngs) and at least $20\,\%$ of [reads](#def-b3-genomics-ngs) show it” reject errors at this coverage?
10. How many heterozygous sites does one human carry (half of the differences between two random genomes lie in each, roughly: take one site in $1500$ )? How many in coding sequence?
11. The [reads](#def-b3-genomics-ngs) are aligned to a reference. Explain why [reads](#def-b3-genomics-ngs) from Alu elements often align to the wrong copy and what consequence this has for calling variants inside them.
12. A child is sequenced with both parents. How many new mutations do you expect? How many [reads](#def-b3-genomics-ngs) in the child, at this coverage, must show a variant absent from *both* parents for it to be believed, and why is the parents’ coverage as important as the child’s?
13. A clinical laboratory captures the exome ( $48\,\mathrm{Mb}$ ) and sequences it at $100\times$ . How many [reads](#def-b3-genomics-ngs) does that take, and why is it cheaper per patient than a genome at $28\times$ although its coverage is higher?

**Part III — The content of a genome.**

14. What fraction of the human genome is Alu, and what fraction of the [reads](#def-b3-genomics-ngs) of question 7 come from Alu elements?
15. The reference is $3.2\,\mathrm{Gb}$ but a diploid cell holds $6.4\,\mathrm{Gb}$ ; and a $16\,\mathrm{Gb}$ onion has about $60\,000$ genes. What is the coding fraction of the onion genome, if its genes average $1.5\,\mathrm{kb}$ of coding sequence?
16. Human and chimpanzee differ at $1.2\,\%$ of aligned bases. How many substitutions is that in $2.9\,\mathrm{Gb}$ of alignable sequence? Split evenly between the two lineages over $6.5$ million years, what mutation rate per base per year does this imply, and per generation of $25$ years? Compare with the direct measurement given in the data.
17. Two vertebrate genome duplications should give up to four copies of each ancestral gene. Humans have about $20\,000$ genes and the invertebrate chordate *Ciona* about $16\,000$ . What fraction of the duplicates has been lost, on the hypothesis that the ancestor had $16\,000$ ?
18. Explain why gene number is a poor measure of organismal complexity, with the alternative splicing count of [Chapter 2](https://one-course.com/books/biology/5/en/chapter/2-non-coding-rnas-and-post-transcriptional-regulation#ch-b3-rna-regulation) as one argument.
19. $8\,\%$ of the genome is under purifying selection but only $1.5\,\%$ codes for protein. What is the rest likely to be, and how would you test one candidate element?

**Part IV — An association study.**

20. $10^{6}$ SNPs are tested. State the Bonferroni threshold for a family-wise error of $0.05$ , and the expected number of false positives at $p < 10^{-5}$ .
21. A disease has prevalence $1\,\%$ . A risk allele of frequency $0.3$ raises the odds by $5\,\%$ . Compute the disease risk of a carrier of two copies, of one copy and of none, taking the non-carrier risk as $0.9\,\%$ . (Multiply the odds.)
22. Two hundred such alleles are found. Explain what a [polygenic score](#def-b3-genomics-gwas) is, and why the score of a person in the top $1\,\%$ may correspond to a several-fold risk although each allele does almost nothing.
23. Explain why a SNP that reaches significance is usually not the causal variant, and what experiment would identify the causal one in the region.
24. An ancient genome from a $40\,000$ -year-old bone is sequenced at $c = 0.5$ . What fraction of it is [read](#def-b3-genomics-ngs) ? Why can this still answer questions about population history that a modern genome alone cannot?
25. Summarise: the number of contigs in the pilot [assembly](#def-b3-genomics-assembly) (question 4), the coverage per haploid genome that gives about $15$ [reads](#def-b3-genomics-ngs) per allele (question 7–8), and the number of new mutations in a child (question 12).

**Solution of Problem 4.1.**

**1.** $0.002\times 6\times 10^{8} = 1.2\times 10^{6}$ [reads](#def-b3-genomics-ngs); $c = 1.2\times 10^{6}\times 150/4.6\times 10^{6} \approx 39$. **2.** $e^{-39} \approx 10^{-17}$: effectively no base uncovered. **3.** $\theta = 0.2$; contigs $= 1.2\times 10^{6} e^{-31} \approx
0$: one contig in theory, gaps only at repeats. **4.** $c = 30\,000\times 150/4.6\times 10^{6} = 0.98$; uncovered $e^{-0.98} = 0.38$, about $1.7\,\mathrm{Mb}$; contigs $= 30\,000
\times e^{-0.78} \approx 13\,700$. **5.** The seven operons give identical [reads](#def-b3-genomics-ngs) and collapse into one $5\,\mathrm{kb}$ contig, entered by seven unique left flanks and left by seven right flanks whose pairing is unknown: $14$ contig ends, seven gaps. **6.** About $4600$ genes; coding $0.88\times 4.6\times 10^{6}
= 4.0\times 10^{6}$ bp, mean gene $880\,\mathrm{bp}$. **7.** $0.998\times 6\times 10^{8}\times 150/3.2\times 10^{9}
\approx 28\times$. **8.** About $14$ [reads](#def-b3-genomics-ngs) per allele. **9.** $28\times 10^{-3}/3 \approx 0.009$ [reads](#def-b3-genomics-ngs) showing a given wrong base — the chance of three or more is about $10^{-7}$ — and $20\,\%$ of $28$ is $5.6$ [reads](#def-b3-genomics-ngs); an error passes neither test, while a real heterozygous allele is expected in $14$. **10.** $3.2\times 10^{9}/1500 \approx 2.1\times 10^{6}$ heterozygous sites; in coding sequence $4.8\times 10^{7}/1500 \approx
32\,000$. **11.** A $150\,\mathrm{bp}$ [read](#def-b3-genomics-ngs) from an Alu matches many of the $1.1$ million copies almost equally well, so it is placed on the wrong copy or given a low mapping confidence. Differences between copies then look like heterozygous variants, and true variants are hidden among them: variant calls inside Alu are usually filtered out, and the elements are blind spots of [short-read](#def-b3-genomics-ngs) sequencing. **12.** $1.2\times 10^{-8}\times 6.4\times 10^{9} \approx 77$ new mutations. About $14$ [reads](#def-b3-genomics-ngs) should show the variant in the child; but if a parent’s site is covered by only five [reads](#def-b3-genomics-ngs), a heterozygous allele is missed with probability $2^{-5} = 3\,\%$ and an inherited variant is miscalled as new — so the parents must be sequenced as deeply as the child. **13.** $4.8\times 10^{7}\times 100/150 = 3.2\times 10^{7}$ [reads](#def-b3-genomics-ngs), twenty times fewer than the genome: the sequencing cost falls twentyfold, more than the cost of the capture step. **14.** $1.1\times 10^{6}\times 300 = 3.3\times 10^{8}$ bp, about $10\,\%$ of the genome; $10\,\%$ of the [reads](#def-b3-genomics-ngs), some $6\times
10^{7}$. **15.** $60\,000\times 1500 = 9\times 10^{7}$ bp; $9\times
10^{7}/1.6\times 10^{10} = 0.56\,\%$. **16.** $0.012\times 2.9\times 10^{9} = 3.5\times 10^{7}$ substitutions, $1.7\times 10^{7}$ per lineage; rate $1.7\times
10^{7}/(2.9\times 10^{9}\times 6.5\times 10^{6}) = 9\times 10^{-10}$ per base per year, $2.3\times 10^{-8}$ per generation of 25 years — about twice the pedigree rate of $1.2\times 10^{-8}$, a known discrepancy that suggests longer past generation times or an older split. **17.** Four copies of $16\,000$ would be $64\,000$; $20\,000$ remain, so of the $48\,000$ extra copies only $4000$ survive: $92\,\%$ of the duplicates were lost. **18.** Worm and human have the same number of genes; complexity lies in how they are used — alternative splicing multiplies one gene into thousands of proteins ($38\,016$ for *Dscam*), regulation combines transcription factors and enhancers in cell-type-specific ways, [non-coding RNAs](https://one-course.com/books/biology/5/en/chapter/2-non-coding-rnas-and-post-transcriptional-regulation#def-b3-rna-regulation-ncrna) add layers, and proteins interact in networks. **19.** Regulatory elements (enhancers, promoters, insulators), [non-coding RNA](https://one-course.com/books/biology/5/en/chapter/2-non-coding-rnas-and-post-transcriptional-regulation#def-b3-rna-regulation-ncrna) genes, splicing and localisation signals, origins and structural sequences. Test an enhancer candidate by placing it before a reporter gene in a transgenic embryo and looking at where the reporter is expressed, then deleting the element in the genome with CRISPR and measuring the neighbouring genes. **20.** $0.05/10^{6} = 5\times 10^{-8}$; at $p < 10^{-5}$, $10^{6}
\times 10^{-5} = 10$ false positives expected. **21.** Non-carrier odds $0.009/0.991 = 0.00908$; one copy odds $0.00954$, risk $0.945\,\%$; two copies odds $0.01001$, risk $0.99\,\%$. **22.** The score sums, over the $200$ alleles, the number of risk copies a person carries weighted by each allele’s log odds ratio. The count has mean $200\times 2\times 0.3 = 120$ and standard deviation $\sqrt{200\times 2\times 0.3\times 0.7} \approx 9$; the top $1\,\%$ carries about $21$ alleles more than average, and $1.05^{21} \approx
2.8$: nearly threefold the odds of an average person, from alleles that individually change risk by $5\,\%$. **23.** Alleles within a block of tens of kilobases are inherited together (linkage disequilibrium), so any of them shows the same association; the genotyped SNP is merely the one on the array. To find the causal variant: sequence the region in cases and controls, narrow the set statistically, then test each candidate — reporter assays for enhancer activity of each allele, editing each allele in cells and measuring the nearby genes’ expression, colocalisation with expression-quantitative-trait signals. **24.** $1 - e^{-0.5} = 39\,\%$ of the genome is [read](#def-b3-genomics-ngs). It gives a direct sample of a population at a known time in the past — allele frequencies before later migrations and admixture, the length of Neanderthal segments (long, because few generations of recombination had broken them) and hence a date for the admixture — which modern genomes can only infer through models. **25.** About $13\,700$ contigs in the pilot; $28\times$ per haploid genome, about $14$ [reads](#def-b3-genomics-ngs) per allele; some $77$ new mutations in a child.