Skip to content

mecfs_bio.build_system.task.locus_plot_helpers

Matplotlib helpers shared by the locus plots: Manhattan panels, recombination-rate tracks from a genetic map, and gene-track panels.

-log10 p is the negative base-10 logarithm of a variant's association p-value. The recombination rate is in cM/Mb (centimorgans per megabase) and comes from a genetic map table with the GMAP_* columns, in whatever genome build the map was made for.

Functions:

  • load_recomb –

    Locus-windowed recombination rate (cM/Mb) from a genetic map task.

  • plot_association –

    Scatter -log10 p against position and return the index of the lead variant.

  • plot_gene_track_panel –

    Gene tracks of the genes on chrom within [start, end], with a GRCh38 x-axis label.

  • plot_recomb –

load_recomb 

load_recomb(
    fetch: Fetch,
    task: Task,
    chrom: int,
    bp_min: int,
    bp_max: int,
) -> pl.DataFrame

Locus-windowed recombination rate (cM/Mb) from a genetic map task.

Source code in mecfs_bio/build_system/task/locus_plot_helpers.py 
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
def load_recomb(
    fetch: Fetch, task: Task, chrom: int, bp_min: int, bp_max: int
) -> pl.DataFrame:
    """Locus-windowed recombination rate (cM/Mb) from a genetic map task."""
    return (
        scan_dataframe_asset(fetch(task.asset_id), task.meta)
        .filter(
            (nw.col(GWASLAB_CHROM_COL) == chrom)
            & (nw.col(GMAP_POS_COL) >= bp_min)
            & (nw.col(GMAP_POS_COL) <= bp_max)
        )
        .select(GWASLAB_CHROM_COL, GMAP_POS_COL, GMAP_RATE_COL)
        .collect()
        .to_polars()
        .sort(GMAP_POS_COL)
    )

plot_association 

plot_association(
    ax: Axes,
    pos: ndarray,
    mlog10p: ndarray,
    label: str,
    significance_line: bool,
) -> int

Scatter -log10 p against position and return the index of the lead variant.

The lead is not drawn differently: with no LD coloring there is nothing to relate it to.

Source code in mecfs_bio/build_system/task/locus_plot_helpers.py 
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
def plot_association(
    ax: Axes,
    pos: np.ndarray,
    mlog10p: np.ndarray,
    label: str,
    significance_line: bool,
) -> int:
    """Scatter -log10 p against position and return the index of the lead variant.

    The lead is not drawn differently: with no LD coloring there is nothing to relate it to.
    """
    ax.scatter(pos, mlog10p, s=8, color="tab:blue").set_rasterized(True)
    lead = int(np.nanargmax(mlog10p))
    if significance_line:
        ax.axhline(
            y=GWAS_SIGNIFICANCE_MLOG10P, color="grey", linestyle="--", linewidth=1
        )
    ax.set_ylabel(f"{label}\n-log10 p")
    return lead

plot_gene_track_panel 

plot_gene_track_panel(
    ax: Axes,
    genes: DataFrame,
    chrom: int,
    start: int,
    end: int,
) -> None

Gene tracks of the genes on chrom within [start, end], with a GRCh38 x-axis label.

Source code in mecfs_bio/build_system/task/locus_plot_helpers.py 
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
def plot_gene_track_panel(
    ax: Axes, genes: pl.DataFrame, chrom: int, start: int, end: int
) -> None:
    """Gene tracks of the genes on chrom within [start, end], with a GRCh38 x-axis label."""
    chrom_name = GWASLAB_CHROM_NAME_FOR_CODE.get(chrom, None) or str(chrom)
    plot_gene_tracks(
        ax=ax,
        gene_df=genes.filter(pl.col(GENE_INFO_CHROM_COL).cast(pl.String) == chrom_name),
        start_bp=start,
        end_bp=end,
        gene_start_col=GENE_INFO_START_COL,
        gene_end_col=GENE_INFO_END_COL,
        gene_name_col=GENE_INFO_NAME_COL,
        gene_strand_col=GENE_INFO_STRAND_COL,
    )
    ax.set_ylabel("genes")
    ax.set_xlabel(f"GRCh38 chr{chrom_name} position (bp)")

plot_recomb 

plot_recomb(ax: Axes, recomb: DataFrame) -> None
Source code in mecfs_bio/build_system/task/locus_plot_helpers.py 
52
53
54
55
56
57
58
59
60
def plot_recomb(ax: Axes, recomb: pl.DataFrame) -> None:
    bp = recomb[GMAP_POS_COL].to_numpy().astype(float)
    rate = recomb[GMAP_RATE_COL].to_numpy().astype(float)
    if len(bp) >= 1:
        (line,) = ax.plot(bp, rate, color="tab:red", alpha=0.4, linewidth=1)
        line.set_rasterized(True)
    # The recomb line is red, but its axis label/ticks stay black to match the
    # other panels' axis labels (the color only carries meaning on the line).
    ax.set_ylabel("cM/Mb")