Structural variants (1000 Genomes)
TL;DR: take one whole-gene deletion out of the 1000 Genomes ensemble SV callset, sort the cohort's genotypes at it, then open three of those samples' reads and watch the coverage go to zero, halve, and stay flat. The page ends on a call from the same cohort whose coverage does none of that.
Prerequisites
- nothing to install: the demo instance already carries the callset and the alignments
Where the data comes from
The 1000 Genomes 2022 high-coverage ensemble SV callset (Byrska-Bishop et al., 2022), read against three of its own samples' CRAMs and QuicK-mer2 copy number for the whole cohort.
- the ensemble SV callset, 3202 samples. EBI publishes it and nobody mirrors it, so the demo reads our own byte-for-byte copy1: https://jbrowse.org/demos/1000g/1KGP_3202.Illumina_ensemble_callset.freeze_V1.vcf.gz
- HG00113 (homozygous alt) high-coverage CRAM: https://1000genomes.s3.amazonaws.com/1000G_2504_high_coverage/data/ERR3240129/HG00113.final.cram
- HG00096 (heterozygous) high-coverage CRAM: https://1000genomes.s3.amazonaws.com/1000G_2504_high_coverage/data/ERR3240114/HG00096.final.cram
- HG00097 (homozygous reference) high-coverage CRAM: https://1000genomes.s3.amazonaws.com/1000G_2504_high_coverage/data/ERR3240115/HG00097.final.cram
- QuicK-mer2 copy number for the cohort, the store the
copy-number tutorial also reads. A directory
of chunks rather than a file, so it is the
urian adapter takes: https://jbrowse.org/demos/1000g/qm2_cn_1kb.zarr
The 1000 Genomes SV callset
The 1000 Genomes Project sequenced genomes from 2,504 individuals across 26 populations. The 2022 high-coverage re-analysis produced a comprehensive SV callset (Byrska-Bishop et al., 2022) that includes deletions, insertions, inversions, and translocations with per-sample genotypes across all 3,202 individuals.
HGSV_1821 is a deletion on chr1 whose span contains the whole of RHD, so the
samples called homozygous here carry no copy of that gene. Deleting RHD is the
most common cause of the RhD-negative blood type; inactivating variants and
RHD-CE hybrid genes produce it too.
The call is PASS and common enough in the cohort to fill all three genotype
classes, and read depth settles whether a gene is present twice, once, or not at
all.
The tracks are added with the usual jbrowse add-track workflow. The callset is
bgzip-compressed and tabix-indexed, and the alignment tracks stream the
published high-coverage CRAMs directly from the urls above.
In the track selector, enable the 1KGP 2022 Illumina ensemble SV callset under
1000 Genomes → SV callsets, listed by its file name
1KGP_3202.Illumina_ensemble_callset.freeze_V1.vcf, then navigate to
chr1:25,200,000-25,400,000, the window the figures below use.
Everything below works the same in JBrowse Desktop, which opens your own VCF and BAM files from disk.
Genotypes across the cohort
Switch the track to the Multi-sample variant display (regular) from the track menu. Each sample becomes a row, drawn at the variant's real genomic span, so the deletion is a wide block rather than a tick. Clicking it opens the feature details panel, whose SAMPLES section lists every sample with its genotype, read depth, and other per-sample fields.
Rows arrive in the callset's own order. Right-click the deletion and pick Sort by genotype: rows order by their genotype at that call first, then by how far each row keeps matching the rows to either side of it. The three classes separate into three contiguous bands. The track menu's Clustering → Cluster rows by genotype... keys every row on the whole window and draws the dendrogram it built.
Dark blue is a sample with no copy of RHD, light blue one copy, grey two, and the olive stripe running through the block is a separate nested call.
A matrix cell says a sample carries something at that column, not which call it carries, and a good many records overlap in this window. Loading the same VCF a second time in the ordinary variant display puts each of them on its own row with its id, class and size, so a band in the matrix reads off a named record. Cell coloring by SV type is the other way to ask that question, shown in the multi-variant track guide.
Three lanes read below, over NCBI RefSeq genes:
- the callset as a genotype matrix, one row per sample, sorted by genotype at the RHD deletion
- the same records drawn ordinarily and colored by SV class, so a band in the matrix reads off a named record
- QuicK-mer2 copy number for 2504 individuals, one row each and clustered on this window, where blue is a copy lost against the diploid white and red a copy gained
The olive no-call column is a copy-number gain in the lane beneath it.
The olive stripe is HGSV_1823, a small copy-number record sitting inside the
deletion, and most of the cohort is uncalled for it. The display gives a no-call
its own color, since a no-call is the caller declining to answer. In the record
lane the same column is a copy-number call.
A genotype is a caller's discrete verdict per record, so a record the caller declined leaves a hole; copy number is one continuous quantity per bin per individual, taken from the reads, so the column that is olive above is red below. The copy-number tutorial is where that lane comes from, and it is the same store.
Reading the genotypes off the reads
The genotypes are the caller's answer, and the reads it read are in the demo too, so every row of the block above can be checked against them.
Open three samples' alignments from 1000 Genomes → Alignments, one per genotype: HG00113 called homozygous alt, HG00096 heterozygous, HG00097 homozygous reference. Two settings make them comparable:
- Turn the pileup off from the track menu's Show... submenu, since at this width the individual reads are a solid mass and it is the coverage curve that carries the comparison.
- Pin each lane's axis from the track menu's Score → Set min/max score..., so the three coverage lanes share one scale.
A closer look at the empty span
Look again at the top row and the deleted span is not quite at zero: reads are there, in a sample that carries no RHD at all.
They belong to the gene next door. RHCE sits just to the right of RHD and is nearly identical to it, so when a sample has no RHD for its reads to come from, some RHCE reads land in the empty footprint instead. An aligner with nowhere better to put a read puts it somewhere, and records how sure it was in the read's mapping quality.
A display setting separates that residue from real coverage. Open HG00113's pileup inside the deleted span and set Color by... → Mapping quality: the ramp runs red at MAPQ 0 through orange to yellow at MAPQ 60, so reads the aligner could not confidently place come out red.
Residual coverage inside a called deletion turns up wherever the deleted sequence has a close paralog. Raising the track's mapping quality filter empties the span the same way.
An SV the coverage cannot see
A deletion moves the coverage, which is why every step so far could be read off the curve. Most structural variants leave it alone. The same demo carries a complex call on chromosome 1 in HG02768, and the profile under it looks like anywhere else on the arm.
Put 1:39,658,200-39,661,800 in the location box and open HG02768's alignments
from 1000 Genomes → Alignments. The ensemble callset draws the call above
the reads, and nothing in the coverage marks where it starts or stops.
Which way the pairs point does. Turn on Track menu → Read connections → SV channels (pairs by orientation): the reads split into one band per orientation class, each with its own coverage curve and its own arcs.
- The normal band holds the flat profile, the same reading the homozygous reference sample gave at RHD.
- The two same-strand bands each carry a bundle of arcs standing on one pair of breakpoints. Pairs that point the same way are the inversion signature.
- The outward-pointing band stays near empty, which is where a tandem duplication would have gone.
The call names a duplicated copy in its INFO.CPX_INTERVALS too, and no band
shows it. A copy that lands beside where it came from leaves pair orientation
alone, so this half of the call is back to reading the coverage, at a size where
the profile's own noise is the same shape.
See also
- Low-mappability regions (SMN)
- SV inspector view
- SV visualization
- Multi-sample variant display
- Clustering rows
- Phased trio analysis (1000 Genomes)
- Structural variants (Dog10K)
- CNV across a population (1000 Genomes)
- Structural variants (Cancer GIAB)
- jbrowse-anywidget
Notes
-
The file EBI publishes lives at https://ftp.1000genomes.ebi.ac.uk/vol1/ftp/data_collections/1000G_2504_high_coverage/working/20210124.SV_Illumina_Integration/1KGP_3202.Illumina_ensemble_callset.freeze_V1.vcf.gz and our copy is byte-for-byte that. The CRAM urls above are the Registry of Open Data mirror of the same 1000 Genomes ftp tree, which answers a range request in a fraction of the time EBI takes. ↩
Feedback on this tutorial is welcome: contact us.