Genomic scores

Genomic scores associate values with genomic positions, alleles, or intervals. GAIn distinguishes three types of genomic scores according to what each record describes:

Type	What a record represents	Examples
Position scores	A genomic position or a run of positions	Conservation scores such as `phastCons100way`
Allele scores	An allele or substitution at a genomic position	Variant-level scores such as `CADD`
Fragment scores	A genomic interval with associated attributes	CNV collections, ATAC fragments, and region annotations

All three types share the GenomicScore base class. It provides common functionality for opening and closing score resources, inspecting score definitions, retrieving records, and filtering data. Each specialized class also provides methods appropriate to the genomic features it represents. This chapter covers the shared Python interface. The following chapters describe the three score types individually.

Accessing a genomic score

The build_score_from_resource_id() function builds a genomic score object from a resource ID and a repository:

import gain.grr as grr

repo = grr.build_genomic_resource_repository()
score = grr.build_score_from_resource_id("hg38/scores/phastCons100way", repo)

If you already have a GenomicResource object, use build_score_from_resource() instead:

resource = repo.get_resource("hg38/scores/phastCons100way")
score = grr.build_score_from_resource(resource)

Both functions are available through gain.grr. They determine the score type from the resource definition and return the corresponding PositionScore, AlleleScore, or FragmentScore object. Each type also has its own builder functions, described in the corresponding chapter.

Opening and closing a genomic score

A GenomicScore object is initially closed. It uses an underlying data source, such as a tabix-indexed file, a bigWig, or a VCF, that must be opened before reading scores. Call open() to open the object and close() when you finish. The open() method returns the score object, so it can be chained to a builder function.

A with block closes the score automatically, including when an error occurs inside the block:

with grr.build_score_from_resource_id(
    "hg38/scores/phastCons100way", repo
).open() as score:
    print(score.get_scores_at_position("chr21", 5_030_000))

The with statement does not open the score automatically, so the .open() call is necessary. Attempting to read a closed score raises ValueError: genomic score <id> is not open. Requiring an explicit open() call makes the cost of initializing a score clear and helps prevent file handles from being left open. You can check whether a score is open using is_open(). This behavior applies to all three genomic score types.

Score definitions

A genomic score resource can provide one or more scores. The resource definition describes these in a scores: block, as explained in the Scores section of the GRR documentation.

In Python, each declared score is represented by a GenomicScoreDef object. It contains information such as the score ID, value type, rules for converting values read from the data file, and the associated histogram configuration. Score definitions can be inspected through score.score_definitions or retrieved individually with get_score_definition():

for score_def in score.score_definitions.values():
    print(score_def.score_id, score_def.value_type)

one = score.get_score_definition("phastCons100way")

parse_value() converts an individual raw value to its declared type. parse_array() performs the same conversion on an array of values. The array method supports efficient reads of large regions without constructing a separate Python object for every record.

Chromosome names and lengths

get_all_chromosomes() returns the contigs represented in the score’s data, in the order stored in the table. has_chromosome() checks whether a particular contig is present.

Chromosome lengths can come from several sources. GAIn may use the reference genome associated with the score, the header of a bigWig file, an estimate from a tabix index, or the extent of records in an in-memory table. These sources are identified by ChromLengthSource. A length obtained from a reference genome or bigWig header can describe the complete chromosome, whereas a tabix estimate provides an upper bound derived from the indexed data rather than the true chromosome length.

By default, GAIn selects the best available source. You can inspect the available sources with chrom_length_sources, request a specific source, and determine which source supplied a result:

with grr.build_score_from_resource_id(
    "hg38/scores/phastCons100way", repo
).open() as score:
    print(score.chrom_length_sources)
    print(score.get_chrom_length("chr21"))
    print(score.get_chrom_length_source("chr21"))
    print(score.get_chrom_length("chr21", source="tabix_estimate"))
    print(score.get_all_chrom_lengths(source="reference_genome"))

get_all_chrom_lengths() returns lengths for the contigs covered by the requested source. The reference-genome lengths reach the score through the statistics/chrom_lengths.json file created during repository repair. If a resource was repaired before that file was introduced, it may have only table-derived lengths until grr_manage repo-repair is run again. Requesting the unavailable reference_genome source then raises a ValueError.

A ValueError is also raised when a chromosome-length operation cannot provide an answer, for example when the score is closed, the contig is absent from the score, the requested source has no length for it, or no available source can determine its length.

Reading score records and arrays

fetch_records() yields the records in a genomic region, one at a time. This is the general interface for reading the underlying score data.

For score types and backends that support it, fetch_region_value_arrays() retrieves entire columns as arrays instead of returning a separate record for each row. This is useful when working with large numbers of score values. Before using it, check supports_region_value_arrays() for the scores you want to retrieve. Support depends on the requested value types as well as the underlying data backend.

The position, allele, and fragment score chapters describe methods for retrieving score values at specific genomic locations.

Filtering scores

compile_filter() converts a boolean expression involving the score’s columns into a ScoreFilter object. This filter can be passed to read methods through the score_filter argument, so that records are selected as they are read.

For example, the following code selects records with a phastCons100way score greater than 0.9:

with grr.build_score_from_resource_id(
    "hg38/scores/phastCons100way", repo
).open() as score:
    keep = score.compile_filter("phastCons100way > 0.9")
    for record in score.fetch_records(
        "chr21", 5_030_000, 5_040_000, score_filter=keep
    ):
        print(record)

Each filter is associated with the score object that created it. A filter compiled for one score cannot be used with another score object, even when their columns have the same names.

Aggregating scores across regions

Allele and fragment scores provide methods for aggregating the values found in a genomic region. These methods are get_allele_scores_in_region_agg() and get_fragment_scores_overlapping_region_agg(). They process the matching records in a single traversal and return an AlleleAggregate or FragmentAggregate object, respectively.

Each aggregate contains the reduced values and information about the records used to calculate them. For allele scores, this includes the matched allele keys, when the call asks for them with the allele_keys argument. For fragment scores, it includes the number of fragments. Keeping these results together ensures that the values and accompanying information describe the same set of records.

The values are returned in the order of the aggregation queries, rather than indexed by score ID. For example, requesting the same score with two different aggregators produces two values, one for each query. The allele and fragment score chapters describe their aggregation methods in more detail.

Score histograms

A genomic score may have a histogram describing the distribution of its values. The get_score_histogram() method retrieves the histogram associated with a score:

hist = score.get_score_histogram("phastCons100way")

The result may be a numeric, categorical, or null histogram. These types are described in Reading a number histogram, Reading a categorical histogram, and Reading a null histogram. Check the histogram type before using attributes specific to that type.

When a resource defines only one score, the score ID can be omitted from get_score_histogram(), get_histogram_config(), get_score_range(), and the two methods that retrieve histogram image URLs. These methods then use that resource’s single score automatically. The dedicated Histograms chapter explains histogram configurations, contents, and retrieval in more detail.

API

GenomicScore is defined in the module gain.genomic_resources.genomic_scores.base and is imported from the package gain.genomic_resources.genomic_scores. GenomicScoreDef is defined in the module gain.genomic_resources.score_def. ScoreFilter is defined in the module gain.genomic_resources.score_filter. GenomicScore is also exported by gain.grr, as grr.GenomicScore. The generic builders build_score_from_resource_id and build_score_from_resource are defined in the module gain.genomic_resources.genomic_scores.builders and are imported from the package gain.genomic_resources.genomic_scores. They are also exported by gain.grr, as grr.build_score_from_resource_id and grr.build_score_from_resource.

`build_score_from_resource`(resource)	Build the score of a score resource, whatever its kind.
`build_score_from_resource_id`(resource_id[, grr])	Build the score of the score resource with a given id.
`GenomicScore`(resource)	Base class of the position, allele and fragment scores.
`GenomicScore.chrom_length_source`	Return the source of chromosome lengths in the score's own table.
`GenomicScore.chrom_length_sources`	The sources that the score can answer a length from, best first.
`GenomicScore.close`()	Close the underlying table and mark the score not open.
`GenomicScore.compile_filter`(expression)	Compile a boolean expression into a filter over this score.
`GenomicScore.fetch_records`(chrom, pos_begin, ...)	Yield the records of a region, optionally filtered.
`GenomicScore.fetch_region_segments_scores`(chrom)	Yield `(begin, end, values)` for each record that touches a region.
`GenomicScore.fetch_region_value_arrays`(...)	Fetch a region as column arrays, without a record for each row.
`GenomicScore.get_all_chrom_lengths`([source])	Return the length of every chromosome.
`GenomicScore.get_all_chromosomes`()	Return the chromosome names that the score's table holds.
`GenomicScore.get_chrom_length`(chrom[, source])	Return the length of a chromosome.
`GenomicScore.get_chrom_length_source`(chrom)	Return the source that answers `get_chrom_length()`.
`GenomicScore.get_config`()	Return the configuration of the score.
`GenomicScore.get_score_value_from_record`(...)	Read one score off a record of this score's table.
`GenomicScore.has_chromosome`(chrom)	Report whether the score's table carries a chromosome.
`GenomicScore.is_open`()	Report whether the score is open.
`GenomicScore.open`()	Open the score for reading and return it.
`GenomicScore.supports_region_value_arrays`(scores)	Report whether `fetch_region_value_arrays()` serves these scores.

`GenomicScoreDef`(score_id, value_type, desc, ...)	The definition of one score of a genomic score resource.
`GenomicScoreDef.parse_array`(cells)	Turn a whole column of raw cells into values, vectorized.
`GenomicScoreDef.parse_value`(value)	Turn one raw cell into the value of this score.

`ScoreFilter`(score, expression, predicate)	A compiled record filter, bound to the score that compiled it.
`ScoreFilter.require_owner`(score)	Refuse the records of a score that did not compile this filter.
`ScoreFilter.select`(score, records)	Yield the records of `score` that this filter accepts.

build_score_from_resource(resource: GenomicResource) → GenomicScore[source]

Build the score of a score resource, whatever its kind.

The function reads the type of the resource and calls the factory of that kind. Use the factory of the kind when the type is known. Use this function when the type is not.

Parameters:: resource (GenomicResource) – A resource of type position_score, allele_score or fragment_score. The legacy type cnv_collection also names a fragment score.
Returns:: A new PositionScore, AlleleScore or FragmentScore. Each call returns a new score that the caller owns, so closing it affects nothing else. The score is closed: call open() before the first read.
Return type:: GenomicScore
Raises:: ValueError – If the resource is not of a score type, or if it declares a retired type such as np_score.

build_score_from_resource_id(resource_id: str, grr: GenomicResourceRepo | None = None) → GenomicScore[source]

Build the score of the score resource with a given id.

The function finds the resource in the repository and calls build_score_from_resource().

Parameters:

resource_id (str) – The id of a score resource in the repository.
grr (GenomicResourceRepo | None, default None) – The repository that resolves resource_id. When it is None, the function builds the default repository with build_genomic_resource_repository().

Returns:

A new PositionScore, AlleleScore or FragmentScore. The score is closed: call open() before the first read.

Return type:

GenomicScore

Raises:

FileNotFoundError – If the repository has no resource with this id.
ValueError – If the resource is not of a score type, or if it declares a retired type such as np_score.

class GenomicScore(resource: GenomicResource)[source]

Base class of the position, allele and fragment scores.

A genomic score gives access to the values of a genomic resource that attaches scores to places in the genome. The class owns what the three kinds share: the score definitions, the open and close lifecycle, the filter compiler, and the record and column-array reads. A score is built closed, so call open() before the first read.

The data of the resource is a tabix file, a VCF, a bigWig or an in-memory table. A score is not safe for concurrent use.

Variables:

resource (GenomicResource) – The genomic resource of the score.
resource_id (str) – The id of the resource.
config (dict) – The validated and normalized configuration.
table (GenomicPositionTable) – The table that holds the data.
score_definitions (dict[str, GenomicScoreDef]) – The definitions of the scores, keyed by score id.
table_loaded (bool) – True while the table is open.
DEFAULT_AGGREGATORS (dict[str, str | None]) – The aggregator that each kind of score value gets when its config states none. The value types of a score are the keys. Each kind of score sets it.