Skip to content

count_bam() silently returns zero counts for numeric chromosome names (non-human genomes) #268

Description

@stephenrong

Bug

count_bam() in workflow/scripts/neighborhoods.py silently returns zero read counts
for a subset of chromosomes when using non-human Ensembl genomes with numeric chromosome
names (e.g. 1, 2, ..., 20).

Cause

The BED file is read with:

bed_regions = pd.read_table(bed_file, header=None)

Pandas infers the dtype of column 0 independently per chunk (default chunk size 2^17 =
131072 rows). If a chunk contains only numeric chromosome names, they are parsed as
integers. However, pysam.AlignmentFile.references returns chromosome names as strings.
The check row.chr in read_chrs then compares e.g. integer 1 to the set {"1", "2", ...}, which evaluates to False, and the count is reported as zero.

Chunks that happen to contain non-numeric chromosome names (e.g. X, Y, or scaffold
accessions like ML143117.1) force pandas to use string dtype for the entire chunk, so
those chromosomes are counted correctly.

Which chromosomes are affected depends on the BED file size and sort order — it is
determined by which chromosomes fall in all-numeric chunks. In our case (Macaca mulatta,
Mmul_10), the first chunk covered chr1 through partway into chr11 (~131k rows). Chr1–10
had entirely zero counts, chr11 transitioned from zero to non-zero at row 131072, and
chr12–20/X/Y were unaffected.

The bug is silent — no error is raised, and zero counts propagate into EnhancerList.txt
and GeneList.txt as zero activity, producing incorrect ABC scores.

Fix

Force string dtype when reading the BED file:

bed_regions = pd.read_table(bed_file, header=None, dtype={0: str})

This is consistent with how read_bed() elsewhere in the same file handles chromosome
names via pd.Categorical.

Environment

  • ABC-Enhancer-Gene-Prediction (current main branch)
  • pandas 1.x / 2.x (any version with chunked type inference)
  • Non-human Ensembl genome (e.g. Macaca mulatta Mmul_10, chromosomes named 1–20, X, Y)
  • Human genomes with chr prefixes are unaffected since all chromosome names are strings

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions