Skip to content

Commit

Permalink
Improved read cloud display for long reads with inversions (#3707)
Browse files Browse the repository at this point in the history
  • Loading branch information
cmdcolin committed May 18, 2023
1 parent fdfd4fa commit d0b284c
Show file tree
Hide file tree
Showing 26 changed files with 913 additions and 922 deletions.
Expand Up @@ -3,7 +3,7 @@ import { createTheme, ThemeProvider } from '@mui/material/styles'
import React from 'react'

const react = jest.requireActual('@testing-library/react')
const render = args => {
const render = (args: React.ReactNode) => {
return react.render(
<ThemeProvider theme={createTheme()}>{args}</ThemeProvider>,
)
Expand Down
19 changes: 12 additions & 7 deletions plugins/alignments/src/LinearReadArcsDisplay/drawFeats.ts
Expand Up @@ -4,13 +4,13 @@ import { Assembly } from '@jbrowse/core/assemblyManager/assembly'

// locals
import {
getOrientationColor,
getInsertSizeColor,
getInsertSizeAndOrientationColor,
getPairedOrientationColor,
getPairedInsertSizeColor,
getPairedInsertSizeAndOrientationColor,
} from '../shared/color'
import { featurizeSA } from '../MismatchParser'
import { LinearReadArcsDisplayModel } from './model'
import { hasPairedReads } from './util'
import { hasPairedReads } from '../shared/util'

type LGV = LinearGenomeViewModel

Expand Down Expand Up @@ -108,11 +108,16 @@ export default function drawFeats(
} else {
if (hasPaired) {
if (type === 'insertSizeAndOrientation') {
ctx.strokeStyle = getInsertSizeAndOrientationColor(k1, k2, stats)[0]
ctx.strokeStyle = getPairedInsertSizeAndOrientationColor(
k1,
k2,
stats,
)[0]
} else if (type === 'orientation') {
ctx.strokeStyle = getOrientationColor(k1)[0]
ctx.strokeStyle = getPairedOrientationColor(k1)[0]
} else if (type === 'insertSize') {
ctx.strokeStyle = getInsertSizeColor(k1, k2, stats)?.[0] || 'grey'
ctx.strokeStyle =
getPairedInsertSizeColor(k1, k2, stats)?.[0] || 'grey'
} else if (type === 'gradient') {
ctx.strokeStyle = `hsl(${Math.log10(absrad) * 10},50%,50%)`
}
Expand Down
191 changes: 9 additions & 182 deletions plugins/alignments/src/LinearReadCloudDisplay/drawFeats.ts
@@ -1,29 +1,14 @@
import { getConf } from '@jbrowse/core/configuration'
import { getContainingView, getSession, max, min } from '@jbrowse/core/util'
import { getContainingView, getSession } from '@jbrowse/core/util'
import { LinearGenomeViewModel } from '@jbrowse/plugin-linear-genome-view'

// locals
import {
getOrientationColor,
getInsertSizeColor,
getInsertSizeAndOrientationColor,
} from '../shared/color'
import { ChainStats, ReducedFeature } from '../shared/fetchChains'
import { LinearReadCloudDisplayModel } from './model'
import { fillRectCtx, strokeRectCtx } from './util'
import { hasPairedReads } from '../shared/util'
import { drawPairChains } from './drawPairChains'
import { drawLongReadChains } from './drawLongReadChains'

type LGV = LinearGenomeViewModel

interface ChainCoord {
distance: number
r1s: number
r1e: number
r2s: number
r2e: number
v0: ReducedFeature
v1: ReducedFeature
}

export default function drawFeats(
self: LinearReadCloudDisplayModel,
ctx: CanvasRenderingContext2D,
Expand All @@ -33,176 +18,18 @@ export default function drawFeats(
return
}
const { assemblyManager } = getSession(self)
const featureHeight = getConf(self, 'featureHeight')
const displayHeight = self.height
const view = getContainingView(self) as LGV
const assemblyName = view.assemblyNames[0]
const asm = assemblyManager.get(assemblyName)
if (!asm) {
return
}

const { chains, stats } = chainData
const type = self.colorBy?.type || 'insertSizeAndOrientation'

function fill({
r1s,
r1e,
r2s,
r2e,
distance,
v0,
v1,
}: {
r1s: number
r1e: number
r2s: number
r2e: number
distance: number
v0: ReducedFeature
v1: ReducedFeature
}) {
const w1 = Math.max(r1e - r1s, 2)
const w2 = Math.max(r2e - r2s, 2)
// eslint-disable-next-line @typescript-eslint/no-non-null-assertion
const [fill, stroke] = getColor({ type, v0, v1, stats: stats! })

const top = (Math.log(distance) - minD) * scaler
const halfHeight = featureHeight / 2 - 0.5
const w = r2s - r1e
fillRectCtx(r1e - view.offsetPx, top + halfHeight, w, 1, ctx, 'black')

strokeRectCtx(r1s - view.offsetPx, top, w1, featureHeight, ctx, stroke)
strokeRectCtx(r2s - view.offsetPx, top, w2, featureHeight, ctx, stroke)
fillRectCtx(r1s - view.offsetPx, top, w1, featureHeight, ctx, fill)
fillRectCtx(r2s - view.offsetPx, top, w2, featureHeight, ctx, fill)
}

const coords: ChainCoord[] = []
for (let i = 0; i < chains.length; i++) {
const chain = chains[i]
// if we're looking at a paired read (flag 1) then assume it is just
// two reads (some small cases may defy this assumption such as
// secondary alignments but this may be uncommon)
if (chain[0].flags & 1) {
if (chain.length > 1) {
const v0 = chain[0]
const v1 = chain[1]
const ra1 = asm.getCanonicalRefName(v0.refName) || v0.refName
const ra2 = asm.getCanonicalRefName(v1.refName) || v1.refName
const r1s = view.bpToPx({ refName: ra1, coord: v0.start })?.offsetPx
const r1e = view.bpToPx({ refName: ra1, coord: v0.end })?.offsetPx
const r2s = view.bpToPx({ refName: ra2, coord: v1.start })?.offsetPx
const r2e = view.bpToPx({ refName: ra2, coord: v1.end })?.offsetPx

let distance = 0

if (
r1s !== undefined &&
r1e !== undefined &&
r2s !== undefined &&
r2e !== undefined
) {
if (v0.refName === v1.refName) {
const s = Math.min(v0.start, v1.start)
const e = Math.max(v0.end, v1.end)
distance = Math.abs(e - s)
}
coords.push({
r1s,
r1e,
r2s,
r2e,
v0,
v1,
distance,
})
}
} else if (self.drawSingletons) {
const v0 = chain[0]

const ra1 = asm.getCanonicalRefName(v0.refName) || v0.refName
const r1s = view.bpToPx({ refName: ra1, coord: v0.start })?.offsetPx
const r1e = view.bpToPx({ refName: ra1, coord: v0.end })?.offsetPx
if (r1s !== undefined && r1e !== undefined) {
const w1 = Math.max(r1e - r1s, 2)
fillRectCtx(r1s - view.offsetPx, 0, w1, featureHeight, ctx, '#f00')
strokeRectCtx(r1s - view.offsetPx, 0, w1, featureHeight, ctx, '#a00')
}
}
} else {
// if we're not looking at pairs, then it could be a
// multiply-split-long read, so traverse chain
for (let i = 1; i < chain.length; i++) {
const v0 = chain[i - 1]
const v1 = chain[i]
const ra1 = asm.getCanonicalRefName(v0.refName) || v0.refName
const ra2 = asm.getCanonicalRefName(v1.refName) || v1.refName
const r1s = view.bpToPx({ refName: ra1, coord: v0.start })?.offsetPx
const r1e = view.bpToPx({ refName: ra1, coord: v0.end })?.offsetPx
const r2s = view.bpToPx({ refName: ra2, coord: v1.start })?.offsetPx
const r2e = view.bpToPx({ refName: ra2, coord: v1.end })?.offsetPx

let distance = 0

if (
r1s !== undefined &&
r1e !== undefined &&
r2s !== undefined &&
r2e !== undefined
) {
if (v0.refName === v1.refName) {
const s = Math.min(v0.start, v1.start)
const e = Math.max(v0.end, v1.end)
distance = Math.abs(e - s)
}
coords.push({
r1s,
r1e,
r2s,
r2e,
v0,
v1,
distance,
})
}
}
}
}

const maxD = Math.log(max(coords.map(c => c.distance)))
const minD = Math.max(Math.log(min(coords.map(c => c.distance))) - 1, 0)
const scaler = (displayHeight - 20) / (maxD - minD)

for (let i = 0; i < coords.length; i++) {
fill(coords[i])
}
}
const hasPaired = hasPairedReads(chainData)

function getColor({
type,
v0,
v1,
stats,
}: {
type: string
v0: ReducedFeature
v1: ReducedFeature
stats: ChainStats
}) {
if (type === 'insertSizeAndOrientation') {
return getInsertSizeAndOrientationColor(v0, v1, stats)
} else if (type === 'orientation') {
return getOrientationColor(v0)
} else if (type === 'insertSize') {
return getInsertSizeColor(v0, v1, stats) || 'grey'
} else if (type === 'gradient') {
const s = Math.min(v0.start, v1.start)
const e = Math.max(v0.end, v1.end)
return [
`hsl(${Math.log10(Math.abs(e - s)) * 10},50%,50%)`,
`hsl(${Math.log10(Math.abs(e - s)) * 10},50%,30%)`,
]
if (hasPaired) {
drawPairChains({ self, view, asm, ctx, chainData })
} else {
drawLongReadChains({ self, view, asm, ctx, chainData })
}
return []
}
@@ -0,0 +1,85 @@
import { getConf } from '@jbrowse/core/configuration'
import { max, min } from '@jbrowse/core/util'
import { Assembly } from '@jbrowse/core/assemblyManager/assembly'
import { LinearGenomeViewModel } from '@jbrowse/plugin-linear-genome-view'

// locals
import { ChainData } from '../shared/fetchChains'
import { LinearReadCloudDisplayModel } from './model'
import { fillColor, strokeColor } from '../shared/color'
import { fillRectCtx, strokeRectCtx } from './util'

export function drawLongReadChains({
ctx,
self,
chainData,
view,
asm,
}: {
ctx: CanvasRenderingContext2D
self: LinearReadCloudDisplayModel
chainData: ChainData
view: LinearGenomeViewModel
asm: Assembly
}) {
const distances: number[] = []
const minXs: number[] = []
const { chains } = chainData
const { height } = self
const featureHeight = getConf(self, 'featureHeight')

// get bounds on the 'distances' (pixel span that a particular split long
// read 'chain' would have in view)
for (const chain of chains) {
let minX = Number.MAX_VALUE
let maxX = Number.MIN_VALUE
for (const elt of chain) {
const refName = asm.getCanonicalRefName(elt.refName) || elt.refName
const rs = view.bpToPx({ refName, coord: elt.start })?.offsetPx
const re = view.bpToPx({ refName, coord: elt.end })?.offsetPx
if (rs !== undefined && re !== undefined) {
minX = Math.min(minX, rs)
maxX = Math.max(maxX, re)
}
}
const distance = Math.abs(maxX - minX)
distances.push(distance)
minXs.push(minX)
}

const maxD = Math.log(max(distances))
const minD = Math.max(Math.log(min(distances)) - 1, 0)
const scaler = (height - 20) / (maxD - minD)
const halfHeight = featureHeight / 2 - 0.5

// draw split long read 'chains' as connected entities
for (let i = 0; i < chains.length; i++) {
const chain = chains[i]
const w = distances[i]
const top = (Math.log(w) - minD) * scaler
const min = minXs[i]
fillRectCtx(min - view.offsetPx, top + halfHeight, w, 1, ctx, 'black')
const c1 = chain[0]
let primaryStrand: undefined | number
if (!(c1.flags & 2048)) {
primaryStrand = c1.strand
} else {
const res = c1.SA?.split(';')[0].split(',')[2]
primaryStrand = res === '-' ? -1 : 1
}
for (const v0 of chain) {
const ra = asm.getCanonicalRefName(v0.refName) || v0.refName
const rs = view.bpToPx({ refName: ra, coord: v0.start })?.offsetPx
const re = view.bpToPx({ refName: ra, coord: v0.end })?.offsetPx
if (rs !== undefined && re !== undefined) {
const w = Math.max(re - rs, 2)
const l = rs - view.offsetPx
const effectiveStrand = v0.strand * primaryStrand
const c =
effectiveStrand === -1 ? 'color_rev_strand' : 'color_fwd_strand'
strokeRectCtx(l, top, w, featureHeight, ctx, strokeColor[c])
fillRectCtx(l, top, w, featureHeight, ctx, fillColor[c])
}
}
}
}

0 comments on commit d0b284c

Please sign in to comment.