Press n or j to go to the next uncovered block, b, p or k for the previous block.
| 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 807 808 809 810 811 812 813 814 815 816 817 818 819 820 821 822 823 824 825 826 827 828 829 830 831 832 833 834 835 836 837 | 2x 19x 19x 19x 19x 19x 19x 19x 5x 5x 19x 6x 25x 25x 25x 6x 19x 3x 7x 7x 7x 3x 19x 4x 4x 19x 1x 19x 19x 2x 269600749x 449890x 449890x 269150859x 269150859x 2x 269150896x 269150896x 169871040x 269150896x 97455635x 269150896x 925686x 269150896x 898535x 269150896x 269150896x 269150896x 2x 41x 41x 27x 41x 9x 41x 4x 41x 1x 41x 41x 41x 2x 65x 65x 65x 58x 58x 58x 65x 28x 28x 28x 28x 65x 65x 2x 373x 373x 373x 3607x 3607x 3607x 3607x 3607x 3607x 3607x 3607x 645886x 645886x 645886x 645886x 645886x 642279x 642279x 642279x 642279x 3607x 2248x 3607x 373x 373x 2x 384x 384x 384x 372x 372x 372x 372x 384x 343x 384x 410453x 410453x 343x 343x 384x 384x 2x 2x 8x 8x 2x 4x 4x 4x 2x 2x 11x 11x 11x 11x 11x 7x 7x 7x 7x 7x 4x 4x 4x 4x 4x 4x 11x 4x 4x 11x 11x 2x 8x 8x 8x 8x 7x 8x 8x 8x 8x 8x 2x 262x 262x 262x 262x 262x 179x 179x 123x 123x 179x 139x 139x 221x 83x 83x 262x 262x 325x 325x 325x 325x 325x 325x 325x 325x 325x 325x 338368823x 338368823x 338368823x 325x 325x 2x 2x 2x 2x 2x 2x 2x 2x 7x 7x 7x 4x 7x 2x 2x 15x 15x 15x 15x 15x 15x 15x 15x 15x 14876x 14876x 15x 15x 15x 14876x 14876x 14876x 15x 15x 8x 8x 8x 8x 8x 6x 6x 8x 246x 246x 246x 6x 6x 3x 3x 3x 2x 2x 1083x 1083x 1083x 1083x 1083x 2x 1083x 1083x 1083x 1083x 1083x 1083x 2x 104x 104x 2x 1748x 1748x 2x 2x 2x 2x 28x 28x 2x 89x 89x 2x 21x 21x 2x 2x 596x 596x 596x 596x 596x 2x 2426x 2426x 2x 2308x 2308x 2308x 2308x 2x 596x 596x 596x 2x 47571x 47571x 2x 560x 560x 2x 300x 300x 300x 2x 1x 1x 2x 6x 6x 6x 2x 1x 1x 1x 1x 1x 1x 1x 2x 2x 2x 2x 2x 2x 2x 1x 1x 1x 1x 2x 2x 1x 1x 2x 2x 2x 2x 2x 2x 2x 2x 2x 1x 1x 1x 2x 1x 1x 1x 1x 1x 2x 1x 1x 1x 1x 1x 2x | /**
* Value transforms, the auto-`vmax` rule, and the one coordinate scale.
*
* The transform and `vmax` halves of this module are a line-by-line mirror of
* `metacooler/ops/show.py` — `_transform_matrix` and the four lines of `show()`
* that pick `vmax`. Spec §7.1 makes that parity a requirement rather than a
* nicety: a figure in the paper and the same view in the browser have to be the
* same picture, and `vmax` is the single number that decides how a contact map
* looks.
*
* The mirror covers raw counts, which is everything the CLI can draw: `show()`
* reads `matrix(balance=False)`. Balanced values are the viewer's own, they are
* fractions rather than counts, and both the zero clamp and the anchored `vmin`
* are wrong for them — see {@link displayValue} and {@link autoRangeOf}, the two
* places where this file deliberately stops mirroring.
*
* The {@link Scale} half is the coordinate model. Spec §8: heatmap and tracks
* share *one* `Scale`, owned by the view state, so a bigWig track and the
* matrix column under it can never drift apart by a pixel.
*/
import type { Extent, Transform } from '../types'
// ---------------------------------------------------------------------------
// Value transform (mirrors `_transform_matrix`, ops/show.py:617)
// ---------------------------------------------------------------------------
/**
* Apply the display transform, exactly as `_transform_matrix` does.
*
* Two details are load-bearing and easy to "improve" into a divergence:
*
* - **`log10` maps non-positive values to `0`, not `-Infinity`.** Python is
* `np.where(arr > 0, np.log10(arr), 0.0)`. Zero is by far the most common
* value in a sparse contact map, so getting this wrong changes the whole
* background of the image.
* - **`NaN > 0` is false**, so under `log10` a NaN becomes `0` — in both
* languages. It survives `log1p` and `none` untouched. That asymmetry is
* inherited from numpy on purpose: an unbalanced row shows as blank under
* `log1p` and as background under `log10` in the CLI too.
*
* The input is copied, never mutated: the caller keeps the raw counts for
* tooltips and for recomputing after a transform change.
*
* @throws RangeError on an unknown transform, with the same wording as the
* Python `ValueError`.
*/
export function transformValues(v: Float32Array, transform: Transform): Float32Array
export function transformValues(v: Float64Array, transform: Transform): Float64Array
export function transformValues(
v: Float32Array | Float64Array,
transform: Transform
): Float32Array | Float64Array {
const n = v.length
const out = v instanceof Float64Array ? new Float64Array(n) : new Float32Array(n)
switch (transform) {
case 'log1p':
for (let i = 0; i < n; i++) out[i] = Math.log1p(v[i])
return out
case 'log10':
for (let i = 0; i < n; i++) {
const x = v[i]
out[i] = x > 0 ? Math.log10(x) : 0
}
return out
case 'log2':
// Deliberately *not* `log10`'s clamp. That one is numpy's, and it is right
// for counts, where zero is the bottom of the picture; `log2` exists for
// ratios, where zero means "none of this in the numerator" and belongs at
// the low end of the scale rather than at its middle. Folding it to `0`
// would draw every empty cell as the neutral colour, which on a ratio map
// reads as "the two files agree here" — the single most misleading thing
// this module could produce.
// `NaN` survives, where `log10` sends it to `0`. That asymmetry is as
// deliberate as the one above: under `log10` a NaN is an unbalanced bin
// and numpy draws it as background, but a NaN under `log2` is a ratio that
// could not be formed — two empty cells, or a bin either file discarded —
// and the raster leaves those transparent rather than claiming agreement.
for (let i = 0; i < n; i++) {
const x = v[i]
out[i] = Number.isNaN(x) ? Number.NaN : x > 0 ? Math.log2(x) : Number.NEGATIVE_INFINITY
}
return out
case 'none':
// A loop rather than `out.set(v)`: `out` is a union of two typed arrays
// and TypeScript refuses to resolve `set` across the union.
for (let i = 0; i < n; i++) out[i] = v[i]
return out
default:
throw new RangeError(`Unknown transform: '${String(transform)}'`)
}
}
/**
* The value a **pixel** is coloured by, which is {@link transformValue} plus one
* exception.
*
* `log10`'s "non-positive becomes `0`" is numpy's rule and it is right for raw
* counts, where every stored value is at least `1` and `0` is therefore the
* bottom of the picture. It is exactly backwards for a *balanced* matrix, whose
* values are fractions: `log10` of those is negative, so the clamp puts the
* empty background *above* every contact in the file and paints it black.
*
* Balancing is a viewer feature — `metacooler show` reads `matrix(balance=False)`
* and has no `--balance` — so there is no parity to keep here. Zero becomes
* `-Infinity`, which {@link rasterize} draws at the low end of the colormap:
* matplotlib's `LogNorm` "under" colour, and byte-identical to the numpy rule
* for a count matrix, where `0` *is* the low end.
*/
export function displayValue(x: number, transform: Transform, balanced = false): number {
if (balanced && transform === 'log10') {
// NaN is not zero: an unbalanced bin has no value, and the raster leaves it
// transparent rather than drawing it at either end of the scale.
return Number.isNaN(x) ? Number.NaN : x > 0 ? Math.log10(x) : Number.NEGATIVE_INFINITY
}
return transformValue(x, transform)
}
/**
* The scalar form of {@link transformValues}, for the handful of single values
* that need it (colourbar tick placement, a tooltip's position on the bar).
*
* Kept as a separate three-line switch rather than being called from the array
* loop: the loop hoists its branch out of a million iterations, and that matters
* more than the duplication. The two must stay in step — the `log10` zero clamp
* is the line that would silently diverge.
*/
export function transformValue(x: number, transform: Transform): number {
switch (transform) {
case 'log1p':
return Math.log1p(x)
case 'log10':
return x > 0 ? Math.log10(x) : 0
case 'log2':
return Number.isNaN(x) ? Number.NaN : x > 0 ? Math.log2(x) : Number.NEGATIVE_INFINITY
case 'none':
return x
default:
throw new RangeError(`Unknown transform: '${String(transform)}'`)
}
}
/**
* Undo {@link transformValues} for a single value — what a colourbar tick label
* needs, because a reader reads counts, not log counts.
*
* `log10` is only a left inverse above zero: everything at or below `1` was
* clamped to `0` on the way out, so `inverseTransform(0, 'log10')` returns `1`
* and cannot distinguish that from a true count of `1`. Label the bottom of a
* `log10` colourbar as "≤ 1" rather than pretending otherwise.
*
* `log2` is a true inverse over the whole range it produces, because it clamps
* nothing: it is `-Infinity` that stands for zero, and `2 ** -Infinity` is `0`.
*/
export function inverseTransform(x: number, transform: Transform): number {
switch (transform) {
case 'log1p':
return Math.expm1(x)
case 'log10':
return 10 ** x
case 'log2':
return 2 ** x
case 'none':
return x
default:
throw new RangeError(`Unknown transform: '${String(transform)}'`)
}
}
// ---------------------------------------------------------------------------
// Auto vmax (mirrors show(), ops/show.py:826-828)
// ---------------------------------------------------------------------------
/**
* Linear-interpolated percentile, reproducing `np.percentile(a, q)` with its
* default `method="linear"` (Hyndman & Fan type 7).
*
* The virtual index is `(n - 1) * q / 100` and the result interpolates between
* the two neighbouring order statistics. A nearest-rank shortcut agrees only
* when the index lands on an integer, which for the 99th percentile happens
* only when `n - 1` is a multiple of 100 — i.e. almost never, and the
* disagreement is a visible change in image brightness.
*
* The two-sided lerp is numpy's `_lerp`, not plain `a + t * (b - a)`: above
* `t = 0.5` numpy anchors on `b` instead, which keeps the last ulp when `a` and
* `b` are far apart. Copied so the two implementations agree to the bit.
*
* @param sorted Values in ascending order. Not copied, not modified.
* @param q Percentile in `[0, 100]`.
*/
export function percentileLinear(sorted: ArrayLike<number>, q: number): number {
const n = sorted.length
if (n === 0) return Number.NaN
if (n === 1) return sorted[0]
const pos = ((n - 1) * q) / 100
const lo = Math.floor(pos)
const hi = Math.ceil(pos)
if (lo === hi) return sorted[lo]
const a = sorted[lo]
const b = sorted[hi]
const t = pos - lo
const diff = b - a
return t < 0.5 ? a + diff * t : b - diff * (1 - t)
}
/**
* The `k`-th smallest element, in place, in expected linear time.
*
* Quickselect with a median-of-three pivot and a Hoare partition. On return the
* array is *partially* ordered around `k`: everything before it is no greater,
* everything after no smaller. That is the whole invariant {@link percentileOf}
* needs, and it is what a full sort spends `log n` extra passes establishing
* for every element nobody is going to look at.
*
* `values` is permuted. Every caller here owns a scratch copy built one line
* earlier, so that is free rather than rude.
*
* Assumes no `NaN`: the collectors above filter on predicates that exclude it,
* and a `NaN` pivot would make the partition loops run off the ends.
*/
export function selectNth(values: Float32Array, k: number): number {
let lo = 0
let hi = values.length - 1
while (lo < hi) {
const mid = (lo + hi) >> 1
const a = values[lo]
const b = values[mid]
const c = values[hi]
// Median of three, by value: any of the three is a member of the range, so
// the partition below is guaranteed to make progress.
const pivot = a < b ? (b < c ? b : a < c ? c : a) : a < c ? a : b < c ? c : b
let i = lo - 1
let j = hi + 1
for (;;) {
do i++
while (values[i] < pivot)
do j--
while (values[j] > pivot)
if (i >= j) break
const swap = values[i]
values[i] = values[j]
values[j] = swap
}
if (k <= j) hi = j
else lo = j + 1
}
return values[k]
}
/**
* {@link percentileLinear} without the sort — the same number, from an unsorted
* array, in O(n).
*
* The two straddling order statistics are all the interpolation reads, and
* selection finds them exactly: `selectNth` places the lower one and leaves
* every larger element after it, so the upper one is the minimum of the tail.
* The lerp that follows is byte-for-byte the one above, so this is a change of
* *cost* and not of value — which matters, because the number it produces sets
* the brightness of every pixel and is compared against numpy in the goldens.
*
* @param values Scratch array, permuted in place.
*/
export function percentileOf(values: Float32Array, q: number): number {
const n = values.length
if (n === 0) return Number.NaN
if (n === 1) return values[0]
const pos = ((n - 1) * q) / 100
const lo = Math.floor(pos)
const hi = Math.ceil(pos)
const a = selectNth(values, lo)
if (lo === hi) return a
let b = values[lo + 1]
for (let i = lo + 2; i < n; i++) {
if (values[i] < b) b = values[i]
}
const t = pos - lo
const diff = b - a
return t < 0.5 ? a + diff * t : b - diff * (1 - t)
}
/** The percentile `metacooler show` uses when `--vmax` is not given. */
export const AUTO_VMAX_PERCENTILE = 99
/**
* `vmax` for a matrix that has already been through {@link transformValues}:
* the 99th percentile of its **non-zero** values.
*
* ```python
* nonzero = mat[mat > 0]
* vmax = float(np.percentile(nonzero, 99)) if nonzero.size > 0 else 1.0
* ```
*
* `mat > 0` is a strict comparison, so it drops zeros *and* NaNs — matching
* `v > 0` here, since `NaN > 0` is false in JavaScript too. Without the
* non-zero filter a sparse MAG matrix, which is mostly zero, would put its 99th
* percentile at or near zero and the image would be uniformly saturated.
*
* `vmin` is always 0; there is no `autoVmin`.
*
* @returns `1.0` when nothing is positive, exactly as Python's fallback.
*/
export function autoVmax(transformed: ArrayLike<number>): number {
return autoVmaxImpl(transformed)
}
/**
* The same percentile, from the **raw** values, applying `transform` as it
* goes — so no `n`-cell transformed array has to exist to hold them.
*
* Only the positive values are collected, and a contact map is sparse: for a
* whole yeast genome at 1 kb that is the difference between 585 MB and the
* size of the pixel list. Float32 because the values are Float32 either way —
* they come from a Float32 matrix and `transformValues` rounded them to Float32
* before, so the order and the two straddling values are identical.
*
* Pass raw counts, never transformed ones: transforming twice is silent and
* darkens the whole image.
*/
export function autoVmaxOf(values: ArrayLike<number>, transform: Transform): number {
const nz = collectDisplay(values, transform, false, (_raw, x) => x > 0)
return nz.length === 0 ? 1.0 : percentileOf(nz, AUTO_VMAX_PERCENTILE)
}
/** The low percentile of the fallback rule in {@link autoRangeOf}. */
export const AUTO_VMIN_PERCENTILE = 1
/** A resolved pair of colour-scale bounds, in transformed units. */
export interface ValueRange {
vmin: number
vmax: number
}
/**
* Both bounds at once — the auto rule the panel uses when neither has been set
* by hand.
*
* **Raw counts** get `show()`'s rule, unchanged and to the bit: `vmin = 0` and
* `vmax` at the 99th percentile of the positive transformed values. Zero is a
* true floor there — a count cannot go below it — so anchoring the scale to it
* is what makes two views of the same file comparable.
*
* **Balanced values** cannot use it. They are fractions, so `log10` of them is
* negative almost everywhere: on a real 10 kb E. coli map, 0.4% of the cells
* transform positive (the diagonal, which balances above 1) and the other 99.6%
* land below a `vmin` of zero — one white panel with a dark diagonal, which is
* the bug this branch exists to fix. Zero is not a floor for these values but an
* arbitrary point in the middle of nothing, so the range becomes the 1st to 99th
* percentile of the cells that carry a contact at all.
*
* Empty cells are excluded from that percentile on purpose. A metagenomic block
* is mostly zeros; including them would put both percentiles on the background
* and produce the same blank panel by a different route.
*
* @param balanced Whether the values are balanced fractions rather than counts.
* It decides both the branch here and how `log10` treats zero
* ({@link displayValue}).
*/
export function autoRangeOf(
values: ArrayLike<number>,
transform: Transform,
balanced = false,
): ValueRange {
if (!balanced) {
const positive = collectDisplay(values, transform, false, (_raw, x) => x > 0)
if (positive.length > 0) {
return { vmin: 0, vmax: percentileOf(positive, AUTO_VMAX_PERCENTILE) }
}
}
const contacts = collectDisplay(
values,
transform,
balanced,
(raw, x) => raw > 0 && Number.isFinite(x),
)
if (contacts.length === 0) return { vmin: 0, vmax: 1.0 }
const vmin = percentileOf(contacts, AUTO_VMIN_PERCENTILE)
const vmax = percentileOf(contacts, AUTO_VMAX_PERCENTILE)
// A constant matrix has no range to speak of; matplotlib collapses one to the
// colormap's low end, and `vmin + 1` is the span that reproduces that.
return vmax > vmin ? { vmin, vmax } : { vmin, vmax: vmin + 1 }
}
/**
* The auto rule for a quantity read against a **middle** rather than against
* zero: a ratio of two files, or of a cell to its distance-decay baseline.
*
* Everything above is built on the assumption that a value has a floor and a
* bright end, which is true of counts and false of a ratio. Two things have to
* change together, and doing only one of them is worse than doing neither:
*
* 1. **The range is symmetric about the neutral value.** A colour scale running
* from the 1st to the 99th percentile of a log ratio is almost never centred
* on zero, so the colour that means "these agree" lands somewhere arbitrary
* in the ramp — and on a diverging colormap it lands on a hue, so a map of
* two identical files comes out uniformly blue. Half-width from the *larger*
* of the two tails, so neither direction is clipped harder than the other.
* 2. **The colormap is diverging** — `render/colormaps.ts`, and the panel's job
* rather than this function's.
*
* The half-width is the 99th percentile of the distance from neutral rather than
* the maximum, for the reason §7.1 clips at all: a handful of cells where one
* file has counts and the other has none are `±Infinity`, and a range that
* accommodated them would compress everything real into the middle two colours.
*
* @param centre The neutral value **in transformed units** — `transformValue(1,
* 'log2')` is `0` for a ratio, and `0` for a difference under `none`.
*/
export function autoSymmetricRange(
display: ArrayLike<number>,
centre = 0,
): ValueRange {
// A neutral point that is not on the scale cannot be the middle of it. This
// is reachable: `log2` is offered on every map now, and `transformValue(0,
// 'log2')` — a difference map's neutral — is `-Infinity`, which would make
// every distance infinite and the range `[-Infinity, -Infinity]`.
if (!Number.isFinite(centre)) return autoRangeOfDisplay(display, 'none')
const spread = collectWhere(display, (x) => Number.isFinite(x))
for (let i = 0; i < spread.length; i++) spread[i] = Math.abs(spread[i] - centre)
// A block with nothing finite in it, and a block whose every cell sits exactly
// on neutral, are both ranges of zero width. `1` is the same fallback span
// `autoRangeOf` uses for a constant matrix, and it puts the whole picture on
// the colormap's middle colour — which is what a map of two identical files
// should look like.
const half = spread.length === 0 ? 1 : percentileOf(spread, AUTO_VMAX_PERCENTILE)
const width = half > 0 ? half : 1
return { vmin: centre - width, vmax: centre + width }
}
/**
* {@link autoRangeOf} over values that have **already** been through
* {@link displayValue} — the pooled raster, rather than the matrix behind it.
*
* The two agree cell for cell when nothing was pooled, because the predicates
* below are the display-space spelling of the raw-space ones:
*
* - `display > 0` is `raw > 0` under `log1p` and `none`, and under `log10` it
* is `raw > 1`, which is what the raw-space rule selects too (`log10` clamps
* everything at or below 1 to `0`);
* - "carries a contact" is `raw > 0 && isFinite(display)`; every transform is
* monotone in `raw`, so that is `isFinite(display) && display > display(0)`
* — and `display(0)` is `-Infinity` exactly in the balanced `log10` case the
* fallback exists for.
*
* Where pooling *has* happened the two deliberately differ, and that is the
* point: the percentile has to describe the pixels being drawn. Scaling
* 146 million cells and then drawing the maximum of every 3 × 3 block put
* ~8.6% of the picture above the 99th percentile instead of 1%, so a whole
* yeast genome came out visibly hotter in the browser than the same file does
* from `metacooler show` (spec §7.1).
*/
export function autoRangeOfDisplay(
display: ArrayLike<number>,
transform: Transform,
balanced = false,
): ValueRange {
if (!balanced) {
const positive = collectWhere(display, (x) => x > 0)
if (positive.length > 0) {
return { vmin: 0, vmax: percentileOf(positive, AUTO_VMAX_PERCENTILE) }
}
}
const floor = displayValue(0, transform, balanced)
const contacts = collectWhere(display, (x) => Number.isFinite(x) && x > floor)
if (contacts.length === 0) return { vmin: 0, vmax: 1.0 }
// Two selections rather than two lookups into one sorted array: each is O(n)
// and the second runs over an array the first already partitioned, which is
// still cheaper than the `n log n` a sort costs to answer both.
const vmin = percentileOf(contacts, AUTO_VMIN_PERCENTILE)
const vmax = percentileOf(contacts, AUTO_VMAX_PERCENTILE)
return vmax > vmin ? { vmin, vmax } : { vmin, vmax: vmin + 1 }
}
/**
* The kept values, unordered. Two passes so the array is sized once.
*
* Deliberately *not* sorted. A 4,096² pooled raster of a dense overview keeps
* ~15 M values, and sorting them to read two order statistics measured at
* 1.8 s — a two-second freeze per colormap change, for a pair of numbers
* {@link percentileOf} extracts in a linear pass.
*/
function collectWhere(
values: ArrayLike<number>,
keep: (x: number) => boolean,
): Float32Array {
const n = values.length
let count = 0
for (let i = 0; i < n; i++) if (keep(values[i])) count++
const out = new Float32Array(count)
let k = 0
for (let i = 0; i < n; i++) {
const x = values[i]
if (keep(x)) out[k++] = x
}
return out
}
/**
* The full extent of the drawable values, which is what a clip *slider* spans —
* as opposed to {@link autoRangeOf}, which is where its two handles start.
*
* `null` when nothing is drawable: an all-NaN block, or an all-zero one under
* `log10` of a balanced matrix.
*/
export function valueExtentOf(
values: ArrayLike<number>,
transform: Transform,
balanced = false,
): ValueRange | null {
let lo = Number.POSITIVE_INFINITY
let hi = Number.NEGATIVE_INFINITY
for (let i = 0; i < values.length; i++) {
const x = displayValue(values[i], transform, balanced)
if (!Number.isFinite(x)) continue
if (x < lo) lo = x
if (x > hi) hi = x
}
return lo <= hi ? { vmin: lo, vmax: hi } : null
}
/**
* The display values whose `(raw, display)` pair passes `keep`, unordered.
*
* Float32 because the values are Float32 either way — they come from a Float32
* matrix — so the order and the two straddling values are identical to what a
* Float64 copy would give.
*/
function collectDisplay(
values: ArrayLike<number>,
transform: Transform,
balanced: boolean,
keep: (raw: number, display: number) => boolean,
): Float32Array {
const n = values.length
let count = 0
for (let i = 0; i < n; i++) {
if (keep(values[i], displayValue(values[i], transform, balanced))) count++
}
const out = new Float32Array(count)
let k = 0
for (let i = 0; i < n; i++) {
const x = displayValue(values[i], transform, balanced)
if (keep(values[i], x)) out[k++] = x
}
return out
}
function autoVmaxImpl(transformed: ArrayLike<number>): number {
const n = transformed.length
let count = 0
for (let i = 0; i < n; i++) if (transformed[i] > 0) count++
if (count === 0) return 1.0
// Float64 for the sort and the interpolation even when the input is Float32:
// the percentile is a single number that sets the whole image's brightness,
// so it is worth carrying at double precision once the values are extracted.
// Float32 rather than Float64, matching its two siblings: the input is a
// Float32 matrix either way, so the order and the two straddling values are
// identical, and the buffer is half the size against the same budget.
const nz = new Float32Array(count)
let k = 0
for (let i = 0; i < n; i++) {
const x = transformed[i]
if (x > 0) nz[k++] = x
}
return percentileOf(nz, AUTO_VMAX_PERCENTILE)
}
// ---------------------------------------------------------------------------
// Coordinate scale (spec §8)
// ---------------------------------------------------------------------------
/**
* A base-pair position on one contig, as {@link Scale}'s genomic layer reports
* it. `contig` is an index into `/contigs`, not a name — names are the UI's
* business and the reader already holds the table.
*/
export interface GenomicPosition {
contig: number
/** 0-based base-pair offset within that contig. */
position: number
}
/** A bin's contig and its base-pair span, as `metacooler/coords.ts` computes it. */
export interface BinLocus {
contig: number
/** 0-based, inclusive. */
start: number
/** Exclusive. Short for the last bin of a contig — cooler's own rule. */
end: number
}
/**
* The fixed-bin coordinate arithmetic, injected rather than imported.
*
* Implemented by `metacooler/coords.ts` over `indexes/chrom_offset` and
* `contigs/length` alone — `/bins` is never read (spec §3.3). Declaring it as
* an interface here keeps `render/` free of any dependency on the reader, so
* the heatmap can be unit-tested with three lines of stub arithmetic and so
* neither module can grow a second, subtly different copy of the mapping.
*/
export interface GenomicMapper {
/** Uniform bin width in base pairs. */
readonly binsize: number
/** Absolute bin index → its contig and base-pair span. */
binToLocus(bin: number): BinLocus
/** Contig index and base-pair position → absolute bin index. */
locusToBin(contig: number, position: number): number
}
function clamp(x: number, lo: number, hi: number): number {
return x < lo ? lo : x > hi ? hi : x
}
/**
* The single mapping between matrix coordinates and screen pixels.
*
* "Matrix coordinates" are **fractional bin indices along one axis of the
* rendered block**: `0` is the left (or top) edge of the block's first cell,
* `1` the boundary between its first and second cells. Integers therefore fall
* on cell boundaries, which is what separator lines are drawn at, and
* `Math.floor` of a coordinate is the cell containing it.
*
* Pixels are **plot-area-relative CSS pixels**, not canvas pixels: the heatmap
* owns a margin for tick labels and a track sits underneath the same plot area,
* so anchoring the scale to the plot rather than the canvas is what lets both
* consume the same object unmodified.
*
* One `Scale` per axis. The view state owns them and hands them to the heatmap
* (`HeatmapView.setScales`) and to every track; nobody else constructs one
* (spec §8).
*/
export class Scale {
private full0: number
private full1: number
private d0: number
private d1: number
private px: number
private mapper: GenomicMapper | null = null
/**
* @param full The block's whole extent in matrix coordinates, normally
* `[0, nbins]`. Also the domain a `reset()` returns to.
* @param pixels Length of the plot area along this axis, in CSS pixels.
*/
constructor(full: Extent, pixels: number) {
this.full0 = full[0]
this.full1 = full[1]
this.d0 = full[0]
this.d1 = full[1]
this.px = pixels
}
/** The whole block, independent of the current zoom. */
get full(): Extent {
return [this.full0, this.full1]
}
/** The currently visible sub-range, in matrix coordinates. */
get domain(): Extent {
return [this.d0, this.d1]
}
/** Width of the visible range in matrix coordinates. Always `> 0`. */
get span(): number {
return this.d1 - this.d0
}
/** Plot-area length along this axis, in CSS pixels. */
get pixels(): number {
return this.px
}
/** CSS pixels per matrix cell — the heatmap's zoom factor. */
get pixelsPerUnit(): number {
return this.px / (this.d1 - this.d0)
}
/** Matrix cells per CSS pixel. */
get unitsPerPixel(): number {
return (this.d1 - this.d0) / this.px
}
/** The genomic layer, if one was attached. */
get genomic(): GenomicMapper | null {
return this.mapper
}
/**
* Replace the block extent, resetting the visible domain to all of it. Called
* when the view switches to a different genome or aggregation level.
*/
setFull(full: Extent): void {
if (!(full[1] > full[0])) {
throw new RangeError(`Scale extent must be non-empty, got [${full[0]}, ${full[1]}]`)
}
this.full0 = full[0]
this.full1 = full[1]
this.reset()
}
/** Resize the plot area. The visible domain is preserved, so zoom changes. */
setPixels(pixels: number): void {
this.px = pixels
}
/** Set the visible sub-range directly. */
setDomain(domain: Extent): void {
if (!(domain[1] > domain[0])) {
throw new RangeError(
`Scale domain must be non-empty, got [${domain[0]}, ${domain[1]}]`
)
}
this.d0 = domain[0]
this.d1 = domain[1]
}
/** Show the whole block again. */
reset(): void {
this.d0 = this.full0
this.d1 = this.full1
}
/** Matrix coordinate → plot-area CSS pixel. */
toPixel(unit: number): number {
return ((unit - this.d0) * this.px) / (this.d1 - this.d0)
}
/** Plot-area CSS pixel → matrix coordinate. */
toUnit(pixel: number): number {
return this.d0 + (pixel * (this.d1 - this.d0)) / this.px
}
/**
* Cell index under a pixel, or `-1` when the pixel falls outside the block.
* `full` is half-open, so a pixel exactly on the trailing edge is outside.
*/
cellAt(pixel: number): number {
const u = Math.floor(this.toUnit(pixel))
return u >= this.full0 && u < this.full1 ? u : -1
}
/** Scroll by a screen distance. Positive `dx` drags the content rightwards. */
panByPixels(dx: number): void {
this.panByUnits(-dx * this.unitsPerPixel)
}
/** Scroll by a distance in matrix coordinates. */
panByUnits(du: number): void {
this.d0 += du
this.d1 += du
}
/**
* Zoom about a fixed screen point: the matrix coordinate under `pixel` stays
* under `pixel`. `factor > 1` zooms in.
*/
zoomAt(pixel: number, factor: number): void {
if (!(factor > 0)) throw new RangeError(`Zoom factor must be positive, got ${factor}`)
const anchor = this.toUnit(pixel)
const span = (this.d1 - this.d0) / factor
const frac = this.px === 0 ? 0 : pixel / this.px
this.d0 = anchor - frac * span
this.d1 = this.d0 + span
}
/**
* Attach the fixed-bin arithmetic so this axis can also speak base pairs.
* Only meaningful for a bin-level axis: at contig or genome level a "cell" is
* a whole contig or genome and has no single coordinate.
*/
attachGenomic(mapper: GenomicMapper | null): void {
this.mapper = mapper
}
/**
* Matrix coordinate → contig and base pair, interpolating inside the bin so a
* track can be drawn at sub-bin precision.
*
* @returns `null` when no {@link GenomicMapper} is attached or the coordinate
* is outside the block.
*/
unitToLocus(unit: number): GenomicPosition | null {
if (this.mapper === null) return null
const bin = Math.floor(unit)
if (bin < this.full0 || bin >= this.full1) return null
const locus = this.mapper.binToLocus(bin)
const width = locus.end - locus.start
const frac = clamp(unit - bin, 0, 1)
return { contig: locus.contig, position: locus.start + frac * width }
}
/** Plot-area CSS pixel → contig and base pair. */
pixelToLocus(pixel: number): GenomicPosition | null {
return this.unitToLocus(this.toUnit(pixel))
}
/**
* Contig and base pair → matrix coordinate, again with sub-bin precision so a
* bigWig interval does not snap to the bin grid.
*
* @returns `null` when no {@link GenomicMapper} is attached.
*/
locusToUnit(contig: number, position: number): number | null {
if (this.mapper === null) return null
const bin = this.mapper.locusToBin(contig, position)
const locus = this.mapper.binToLocus(bin)
const width = locus.end - locus.start
const frac = width > 0 ? clamp((position - locus.start) / width, 0, 1) : 0
return bin + frac
}
/** Contig and base pair → plot-area CSS pixel. */
locusToPixel(contig: number, position: number): number | null {
const unit = this.locusToUnit(contig, position)
return unit === null ? null : this.toPixel(unit)
}
/** An independent copy, sharing the genomic mapper. */
clone(): Scale {
const copy = new Scale([this.full0, this.full1], this.px)
copy.setDomain([this.d0, this.d1])
copy.attachGenomic(this.mapper)
return copy
}
/** Structural equality of extent, domain and viewport — for change detection. */
equals(other: Scale): boolean {
return (
this.full0 === other.full0 &&
this.full1 === other.full1 &&
this.d0 === other.d0 &&
this.d1 === other.d1 &&
this.px === other.px
)
}
}
|