Flagging

How a statistic becomes a Pass, Warning, Fail or Unknown.

testing

Turn a per-class statistic into a Pass, Warning, Fail or Unknown flag.

Every check asks one question: can a single feature tell the two classes apart? The feature's values are read as if they were a classifier's scores, and the AU-ROC is the answer - 0.5 for a feature that says nothing about the class, 1.0 for one that settles it by itself. A fixed boundary makes that a flag: Pass up to 0.6, Warning up to 0.7, Fail above that.

Two floors sit in front of the scoring, and a check has to clear both. Nothing is scored on fewer than MIN_SEQUENCES_PER_CLASS sequences per class - counted per position for the per-position checks, since a position is compared only on the sequences long enough to have it. And a per-position check additionally stops where fewer than DEFAULT_MIN_COVERAGE of a class reaches the position. A check that misses either floor reports Unknown, which says the comparison was not made, not that it came out clean.

This module is the arithmetic the checks share. The checks themselves - which feature each one compares, and how the report shows it - are one module each in genomic_benchmarks_qc.checks, and flag_significant_differences runs all of them on one pair of classes for evaluate-classes. Why the boundaries and the floors are where they are is on the How a flag is decided page.

flag_significant_differences

flag_significant_differences(stats1, stats2)

Run every check on one pair of classes.

The checks are the ones genomic_benchmarks_qc.checks.classes lists, in its order.

Parameters:

Name	Type	Description	Default
stats1	SequenceStatistics	Computed statistics for the first class.	required
stats2	SequenceStatistics	Computed statistics for the second class.	required

Returns:

Type	Description
dict	Tuple of (summary_statuses, failed_by_feature). summary_statuses is
dict	every flag and metric in report order: the headline checks first, then
tuple[dict, dict]	their sub-checks. failed_by_feature is what the plots shade - only the
tuple[dict, dict]	sub-checks that came out Warning or Fail, keyed by feature and then by
tuple[dict, dict]	what was flagged: { 'Per sequence nucleotide content': {'A': 'Warning', 'G': 'Fail'}, 'Per sequence dinucleotide content': {'GG': 'Fail'}, 'Per position nucleotide content': {'A': {52: 'Warning'}, 'G': {66: 'Fail', 70: 'Fail'}}, 'Per position reversed nucleotide content': {}, }

Source code in src/genomic_benchmarks_qc/utils/testing.py

def flag_significant_differences(stats1: 'SequenceStatistics',
                                 stats2: 'SequenceStatistics') -> tuple[dict, dict]:
    """Run every check on one pair of classes.

    The checks are the ones `genomic_benchmarks_qc.checks.classes` lists, in its
    order.

    Args:
        stats1: Computed statistics for the first class.
        stats2: Computed statistics for the second class.

    Returns:
        Tuple of `(summary_statuses, failed_by_feature)`. `summary_statuses` is
        every flag and metric in report order: the headline checks first, then
        their sub-checks. `failed_by_feature` is what the plots shade - only the
        sub-checks that came out Warning or Fail, keyed by feature and then by
        what was flagged:

            {
                'Per sequence nucleotide content': {'A': 'Warning', 'G': 'Fail'},
                'Per sequence dinucleotide content': {'GG': 'Fail'},
                'Per position nucleotide content':
                    {'A': {52: 'Warning'}, 'G': {66: 'Fail', 70: 'Fail'}},
                'Per position reversed nucleotide content': {},
            }
    """
    # Imported here, for the reason in `direct_feature_model`.
    from genomic_benchmarks_qc.checks import classes as class_checks
    from genomic_benchmarks_qc.checks import run_checks

    results, failed_by_feature = run_checks(class_checks.CLASS_CHECKS, stats1, stats2)
    _warn_about_unscored_checks(stats1, stats2, results, class_checks.CLASS_CHECKS)
    return results, failed_by_feature

direct_feature_model

direct_feature_model(stats1, stats2)

Score every feature of one class against the same feature of another.

The checks that score a feature - the ones in the registry with a floor, which leaves out those decided by a rule - in registry order.

Parameters:

Name	Type	Description	Default
stats1	SequenceStatistics	Computed statistics for the first class.	required
stats2	SequenceStatistics	Computed statistics for the second class.	required

Returns:

Type	Description
dict	One entry per check and per sub-check, keyed by name - `'Sequence
dict	lengths','Per sequence nucleotide content - A','Per position
dict	nucleotide content - A position 52'` - each holding its AU-ROC, AU-PR,
dict	Accuracy and Flag. A check made of sub-checks also gets a headline entry
dict	under its own name, holding the worst of them.

Source code in src/genomic_benchmarks_qc/utils/testing.py

def direct_feature_model(stats1: 'SequenceStatistics',
                         stats2: 'SequenceStatistics') -> dict:
    """Score every feature of one class against the same feature of another.

    The checks that score a feature - the ones in the registry with a `floor`,
    which leaves out those decided by a rule - in registry order.

    Args:
        stats1: Computed statistics for the first class.
        stats2: Computed statistics for the second class.

    Returns:
        One entry per check and per sub-check, keyed by name - `'Sequence
        lengths'`, `'Per sequence nucleotide content - A'`, `'Per position
        nucleotide content - A position 52'` - each holding its AU-ROC, AU-PR,
        Accuracy and Flag. A check made of sub-checks also gets a headline entry
        under its own name, holding the worst of them.
    """
    # Imported here: every check is built on the scorers in this module, so the
    # registry imports this one, and importing it back at the top would be a
    # cycle.
    from genomic_benchmarks_qc.checks import classes as class_checks

    results = {}
    for check in class_checks.CLASS_CHECKS:
        if check.floor is not None:
            results.update(check.score(stats1, stats2).rows)
    return results

position_windows

position_windows(stats1, stats2)

The per-position windows a comparison of two classes runs in.

Each class resolves its own windows from its own sequence lengths, and a position belongs to a window only where it belongs to it in both: a position that half of one class does not reach cannot be flagged on the strength of the other class reaching it.

Returns:

Type	Description
int	Tuple of (end_position, scored_end_position), 1-based and inclusive - the
int	last position reported on and the last position allowed to set a flag,
tuple[int, int]	which is also the last position drawn.

Source code in src/genomic_benchmarks_qc/utils/testing.py

def position_windows(stats1, stats2) -> tuple[int, int]:
    """The per-position windows a comparison of two classes runs in.

    Each class resolves its own windows from its own sequence lengths, and a
    position belongs to a window only where it belongs to it in both: a position
    that half of one class does not reach cannot be flagged on the strength of
    the other class reaching it.

    Returns:
        Tuple of (end_position, scored_end_position), 1-based and inclusive - the
        last position reported on and the last position allowed to set a flag,
        which is also the last position drawn.
    """
    return (min(stats1.end_position, stats2.end_position),
            min(stats1.scored_end_position, stats2.scored_end_position))

MIN_SEQUENCES_PER_CLASS module-attribute

MIN_SEQUENCES_PER_CLASS = 250