Volcano plot

Effect size against −log10 p-value, with color highlighting points that clear both an effect-size and a significance cutoff.

Use the two sliders to define how large an effect must be and how strong its p-value evidence must be. The guide lines move with the sliders, and points that pass both cutoffs are recolored in the browser.

Data use and provenance

The packaged table is a subset of HapMap variants from the MIT-licensed manhattanly and Plotly datasets projects. The association p-values and effect sizes are simulated, so the plot is for visualization demonstration only. During data loading, the package calculates -log10(p). GenomeSpy then renders the prepared table and classifies the points again whenever a threshold slider moves.

Code

"""Volcano plot.

Effect size against −log10 p-value, with color highlighting points that clear
both an effect-size and a significance cutoff.
"""

import genome_spy as gs
from genome_spy.datasets._hapmap import hapmap_volcano_data
from genome_spy.schema import Scale


# Load association results and starting axis limits and cutoffs.
data, domains = hapmap_volcano_data()
X_DOMAIN = domains["x_domain"]
Y_DOMAIN = domains["y_domain"]
EFFECT_CUTOFF = domains["effect_cutoff"]
NEGLOG_P_CUTOFF = domains["neglog_pvalue_cutoff"]
# Grow points gently when zooming in, without letting them become too large.
ZOOM_LEVEL = gs.Expression("zoomLevel")
POINT_SIZE = gs.expr(gs.expr.min(16 * gs.expr.pow(ZOOM_LEVEL, 0.75), 64))

# Add sliders that update the guide lines and point colors together.
effect_cutoff = gs.param(
    "hapmapEffectCutoff",
    value=EFFECT_CUTOFF,
    bind=gs.binding_range(
        min=0,
        max=X_DOMAIN[1],
        step=0.1,
        name="Absolute effect cutoff: ",
    ),
)
significance_cutoff = gs.param(
    "hapmapSignificanceCutoff",
    value=NEGLOG_P_CUTOFF,
    bind=gs.binding_range(
        min=0,
        max=Y_DOMAIN[1],
        step=0.25,
        name="−log10 p cutoff: ",
    ),
)
# Classify points as protective or risk only when they pass both cutoffs.
ASSOCIATION_EXPRESSION = gs.expr.if_(
    (gs.datum.neglog >= significance_cutoff)
    & (gs.expr.abs(gs.datum.EFFECTSIZE) >= effect_cutoff),
    gs.expr.if_(gs.datum.EFFECTSIZE < 0, "protective", "risk"),
    "n.s.",
)

# Use blue for protective effects, red for risk, and grey for the rest.
association_colors = Scale(
    domain=["protective", "n.s.", "risk"],
    range=["#3e8cb6", "#c9d1d9", "#c53b2c"],
)

# Plot effect size against significance, coloring points by the current cutoffs.
points = (
    gs.Chart()
    .transform_collect()
    .transform_formula(expr=ASSOCIATION_EXPRESSION, as_="association")
    .mark_point(size=POINT_SIZE, filled=True, opacity=0.6)
    .encode(
        x=gs.X("EFFECTSIZE:Q")
        .scale(domain=X_DOMAIN, zoom=True)
        .title("Effect size (beta)"),
        y=gs.Y("neglog:Q")
        .scale(reverse=False, domain=Y_DOMAIN, zoom=True)
        .title("−log10 p"),
        color=gs.Color("association:N")
        .scale(association_colors)
        .legend(title="Association"),
    )
)

# Mark the effect-size cutoff on both sides of zero.
effect_cutoffs = (
    gs.Chart([{"side": -1}, {"side": 1}])
    .transform_collect()
    .transform_formula(expr=gs.datum.side * effect_cutoff, as_="x")
    .mark_rule(strokeDash=[4, 4], size=1, color="#8f98a3")
    .encode(x=gs.X("x:Q").scale(domain=X_DOMAIN, zoom=True).title("Effect size (beta)"))
)

# Move the horizontal significance line with its slider.
significance_rule = (
    gs.Chart([{}])
    .transform_collect()
    .transform_formula(expr=significance_cutoff, as_="y")
    .mark_rule(strokeDash=[4, 4], size=1, color="#8f98a3")
    .encode(
        y=gs.Y("y:Q").scale(reverse=False, domain=Y_DOMAIN, zoom=True).title("−log10 p")
    )
)

# Combine the guide lines and points, then attach the sliders.
chart = (
    (effect_cutoffs + significance_rule + points)
    .properties(
        title="HapMap association volcano",
        data=data,
        description=(
            "Effect size versus significance, with interactive effect-size and "
            "p-value cutoffs controlling the guide lines and point colors."
        ),
    )
    .add_params(effect_cutoff, significance_cutoff)
)