{"spec_id":"ma-differential-expression","library":"bokeh","language":"python","code":"\"\"\" anyplot.ai\nma-differential-expression: MA Plot for Differential Expression\nLibrary: bokeh 3.9.1 | Python 3.13.14\nQuality: 89/100 | Updated: 2026-06-21\n\"\"\"\n\n# Remove the script's own directory from sys.path so \"bokeh.py\" doesn't shadow the package\nimport os as _os\nimport sys as _sys\n\n\n_sys.path = [p for p in _sys.path if _os.path.abspath(p) != _os.path.dirname(_os.path.abspath(__file__))]\n\nimport os\nimport time\nfrom pathlib import Path\n\nimport numpy as np\nfrom bokeh.io import output_file, save\nfrom bokeh.models import BasicTicker, ColorBar, ColumnDataSource, HoverTool, Label, LinearColorMapper, Range1d, Span\nfrom bokeh.plotting import figure\nfrom bokeh.resources import CDN\nfrom bokeh.transform import transform\nfrom selenium import webdriver\nfrom selenium.webdriver.chrome.options import Options\n\n\nTHEME = os.getenv(\"ANYPLOT_THEME\", \"light\")\n\n# Theme-adaptive chrome tokens\nPAGE_BG = \"#FAF8F1\" if THEME == \"light\" else \"#1A1A17\"\nELEVATED_BG = \"#FFFDF6\" if THEME == \"light\" else \"#242420\"\nINK = \"#1A1A17\" if THEME == \"light\" else \"#F0EFE8\"\nINK_SOFT = \"#4A4A44\" if THEME == \"light\" else \"#B8B7B0\"\nINK_MUTED = \"#6B6A63\" if THEME == \"light\" else \"#A8A79F\"\n\n\n# Imprint sequential colormap (green → blue) for significance strength\ndef _lerp_hex(c0, c1, t):\n    r0, g0, b0 = (int(c0[i : i + 2], 16) for i in (1, 3, 5))\n    r1, g1, b1 = (int(c1[i : i + 2], 16) for i in (1, 3, 5))\n    r = int(round(r0 + (r1 - r0) * t))\n    g = int(round(g0 + (g1 - g0) * t))\n    b = int(round(b0 + (b1 - b0) * t))\n    return f\"#{r:02X}{g:02X}{b:02X}\"\n\n\nIMPRINT_SEQ256 = [_lerp_hex(\"#009E73\", \"#4467A3\", t / 255.0) for t in range(256)]\n\n# Data: Simulated RNA-seq differential expression results\nnp.random.seed(42)\nn_genes = 15000\n\ngene_prefixes = [\n    \"BRCA\",\n    \"TP\",\n    \"MYC\",\n    \"EGFR\",\n    \"KRAS\",\n    \"PTEN\",\n    \"AKT\",\n    \"MAPK\",\n    \"STAT\",\n    \"JAK\",\n    \"CDK\",\n    \"RB\",\n    \"VEGF\",\n    \"HIF\",\n    \"TNF\",\n    \"IL\",\n    \"NOTCH\",\n    \"WNT\",\n    \"SHH\",\n    \"FGF\",\n    \"SOX\",\n    \"PAX\",\n    \"HOX\",\n    \"GATA\",\n    \"FOXP\",\n    \"RUNX\",\n    \"ETS\",\n    \"MMP\",\n    \"COL\",\n    \"FN\",\n]\ngene_names = [f\"{gene_prefixes[i % len(gene_prefixes)]}{i}\" for i in range(n_genes)]\n\nmean_expression = np.random.gamma(shape=2.5, scale=2.5, size=n_genes)\nlog_fold_change = np.random.normal(0, 0.4, n_genes)\n\nn_de = int(n_genes * 0.08)\nde_indices = np.random.choice(n_genes, n_de, replace=False)\nlog_fold_change[de_indices] = np.random.choice([-1, 1], n_de) * np.random.uniform(1.0, 4.0, n_de)\n\np_values = np.ones(n_genes)\np_values[de_indices] = 10 ** (-np.random.uniform(2, 10, n_de))\np_values[~np.isin(np.arange(n_genes), de_indices)] = np.random.uniform(0.01, 1.0, n_genes - n_de)\n\nsignificant = p_values < 0.05\nneg_log10_p = -np.log10(np.clip(p_values, 1e-15, 1.0))\n\nsig_mask = significant\nnonsig_mask = ~significant\n\nsource_nonsig = ColumnDataSource(\n    data={\n        \"x\": mean_expression[nonsig_mask],\n        \"y\": log_fold_change[nonsig_mask],\n        \"gene\": [gene_names[i] for i in np.where(nonsig_mask)[0]],\n        \"pval\": p_values[nonsig_mask],\n    }\n)\n\nsource_sig = ColumnDataSource(\n    data={\n        \"x\": mean_expression[sig_mask],\n        \"y\": log_fold_change[sig_mask],\n        \"gene\": [gene_names[i] for i in np.where(sig_mask)[0]],\n        \"pval\": p_values[sig_mask],\n        \"neg_log10_p\": neg_log10_p[sig_mask],\n    }\n)\n\nx_limit = np.percentile(mean_expression, 99)\n\nplot = figure(\n    width=3200,\n    height=1800,\n    title=\"ma-differential-expression · python · bokeh · anyplot.ai\",\n    x_axis_label=\"Mean Expression (log₂)\",\n    y_axis_label=\"Log₂ Fold Change (M)\",\n    toolbar_location=None,\n    x_range=Range1d(-0.5, x_limit + 0.5),\n    min_border_bottom=160,\n    min_border_left=180,\n    min_border_top=110,\n    min_border_right=80,\n)\n\n# Imprint sequential colormap for significance strength (single-polarity continuous)\ncolor_mapper = LinearColorMapper(\n    palette=IMPRINT_SEQ256,\n    low=neg_log10_p[sig_mask].min() if sig_mask.sum() > 0 else 0,\n    high=neg_log10_p[sig_mask].max() if sig_mask.sum() > 0 else 10,\n)\n\n# Non-significant genes (muted, behind)\nr_nonsig = plot.scatter(\n    x=\"x\", y=\"y\", source=source_nonsig, size=5, color=INK_MUTED, alpha=0.12, legend_label=\"Not Significant\"\n)\n\n# Significant genes with color mapped by -log10(p) using Imprint sequential cmap\nr_sig = plot.scatter(\n    x=\"x\",\n    y=\"y\",\n    source=source_sig,\n    size=9,\n    color=transform(\"neg_log10_p\", color_mapper),\n    alpha=0.7,\n    legend_label=f\"Significant (n={sig_mask.sum():,})\",\n)\n\nhover = HoverTool(\n    renderers=[r_sig],\n    tooltips=[(\"Gene\", \"@gene\"), (\"Mean Expr\", \"@x{0.2f}\"), (\"Log₂ FC\", \"@y{0.2f}\"), (\"p-value\", \"@pval{0.2e}\")],\n)\nplot.add_tools(hover)\n\n# Reference lines: M=0 (solid) and M=±1 (dashed fold-change thresholds)\nzero_line = Span(location=0, dimension=\"width\", line_color=INK_SOFT, line_width=2, line_alpha=0.8)\nplot.add_layout(zero_line)\n\nupper_fc = Span(location=1, dimension=\"width\", line_color=INK_SOFT, line_width=1.5, line_dash=\"dashed\", line_alpha=0.6)\nplot.add_layout(upper_fc)\n\nlower_fc = Span(location=-1, dimension=\"width\", line_color=INK_SOFT, line_width=1.5, line_dash=\"dashed\", line_alpha=0.6)\nplot.add_layout(lower_fc)\n\n# LOESS approximation via binned moving average\nsort_idx = np.argsort(mean_expression)\nx_sorted = mean_expression[sort_idx]\ny_sorted = log_fold_change[sort_idx]\n\nvis_mask = x_sorted <= x_limit\nx_vis = x_sorted[vis_mask]\ny_vis = y_sorted[vis_mask]\n\nn_bins = 80\nbin_edges = np.linspace(x_vis.min(), x_vis.max(), n_bins + 1)\nbin_centers = []\nbin_means = []\nfor i in range(n_bins):\n    in_bin = (x_vis >= bin_edges[i]) & (x_vis < bin_edges[i + 1])\n    if in_bin.sum() > 5:\n        bin_centers.append((bin_edges[i] + bin_edges[i + 1]) / 2)\n        bin_means.append(np.mean(y_vis[in_bin]))\n\nbin_centers = np.array(bin_centers)\nbin_means = np.array(bin_means)\n\nwindow = 11\npad = window // 2\npadded = np.pad(bin_means, pad, mode=\"edge\")\ny_smooth = np.convolve(padded, np.ones(window) / window, mode=\"valid\")\n\nsmooth_source = ColumnDataSource(data={\"x\": bin_centers, \"y\": y_smooth})\n# Imprint lavender (#C475FD, position 2) — distinct from the green→blue sig colormap\nplot.line(x=\"x\", y=\"y\", source=smooth_source, line_width=3.5, color=\"#C475FD\", alpha=0.85, legend_label=\"LOESS Trend\")\n\n# Fold-change threshold labels\nlabel_x = x_limit * 0.88\n\nupper_label = Label(\n    x=label_x,\n    y=1,\n    text=\"FC = 2\",\n    text_font_size=\"26pt\",\n    text_color=INK_SOFT,\n    text_baseline=\"bottom\",\n    y_offset=5,\n    text_font_style=\"italic\",\n)\nplot.add_layout(upper_label)\n\nlower_label = Label(\n    x=label_x,\n    y=-1,\n    text=\"FC = −2\",\n    text_font_size=\"26pt\",\n    text_color=INK_SOFT,\n    text_baseline=\"top\",\n    y_offset=-5,\n    text_font_style=\"italic\",\n)\nplot.add_layout(lower_label)\n\n# Annotate top DE genes by |fold change| × −log10(p)\nsig_indices = np.where(sig_mask)[0]\nif len(sig_indices) > 0:\n    de_score = np.abs(log_fold_change[sig_indices]) * neg_log10_p[sig_indices]\n    top_n = 5\n    top_local = np.argsort(de_score)[-top_n:]\n    top_global = sig_indices[top_local]\n\n    y_offsets = [6, 30, -18, 50, -40]\n    for rank, idx in enumerate(top_global):\n        gx = mean_expression[idx]\n        gy = log_fold_change[idx]\n        if gx <= x_limit:\n            gene_label = Label(\n                x=gx,\n                y=gy,\n                text=f\" {gene_names[idx]}\",\n                text_font_size=\"24pt\",\n                text_color=INK,\n                text_font_style=\"bold\",\n                x_offset=10,\n                y_offset=y_offsets[rank % len(y_offsets)],\n            )\n            plot.add_layout(gene_label)\n\n# Summary annotation: up/downregulated counts\nn_up = int(np.sum(significant & (log_fold_change > 1)))\nn_down = int(np.sum(significant & (log_fold_change < -1)))\nsummary_label = Label(\n    x=70,\n    y=70,\n    x_units=\"screen\",\n    y_units=\"screen\",\n    text=f\"▲ {n_up} upregulated  ·  ▼ {n_down} downregulated  (|FC| > 2, p < 0.05)\",\n    text_font_size=\"24pt\",\n    text_color=INK_SOFT,\n    text_font_style=\"italic\",\n)\nplot.add_layout(summary_label)\n\n# ColorBar for significance gradient\ncolor_bar = ColorBar(\n    color_mapper=color_mapper,\n    ticker=BasicTicker(desired_num_ticks=6),\n    label_standoff=14,\n    major_label_text_font_size=\"26pt\",\n    major_label_text_color=INK_SOFT,\n    title=\"-log₁₀(p-value)\",\n    title_text_font_size=\"28pt\",\n    title_text_font_style=\"italic\",\n    title_text_color=INK,\n    title_standoff=12,\n    width=30,\n    location=(0, 0),\n    padding=20,\n    background_fill_color=ELEVATED_BG,\n    background_fill_alpha=1,\n)\nplot.add_layout(color_bar, \"right\")\n\n# Typography — canonical bokeh sizes for 3200×1800 canvas\nplot.title.text_font_size = \"50pt\"\nplot.title.text_color = INK\nplot.xaxis.axis_label_text_font_size = \"42pt\"\nplot.yaxis.axis_label_text_font_size = \"42pt\"\nplot.xaxis.major_label_text_font_size = \"34pt\"\nplot.yaxis.major_label_text_font_size = \"34pt\"\nplot.xaxis.axis_label_text_color = INK\nplot.yaxis.axis_label_text_color = INK\nplot.xaxis.major_label_text_color = INK_SOFT\nplot.yaxis.major_label_text_color = INK_SOFT\nplot.xaxis.axis_line_color = INK_SOFT\nplot.yaxis.axis_line_color = INK_SOFT\nplot.xaxis.axis_line_width = 1.5\nplot.yaxis.axis_line_width = 1.5\nplot.xaxis.minor_tick_line_color = None\nplot.yaxis.minor_tick_line_color = None\nplot.xaxis.major_tick_line_color = INK_SOFT\nplot.yaxis.major_tick_line_color = INK_SOFT\n\n# Grid — y-axis only for cleaner look\nplot.xgrid.grid_line_color = None\nplot.ygrid.grid_line_color = INK\nplot.ygrid.grid_line_alpha = 0.15\n\n# Background and borders\nplot.background_fill_color = PAGE_BG\nplot.border_fill_color = PAGE_BG\nplot.outline_line_color = None\n\n# Legend\nplot.legend.label_text_font_size = \"28pt\"\nplot.legend.label_text_color = INK_SOFT\nplot.legend.location = \"top_left\"\nplot.legend.background_fill_color = ELEVATED_BG\nplot.legend.border_line_color = INK_SOFT\nplot.legend.border_line_width = 1\nplot.legend.padding = 12\nplot.legend.spacing = 8\nplot.legend.margin = 15\n\n# Output paths relative to this script's directory (so files land correctly\n# regardless of the working directory the script is invoked from)\nOUT_DIR = Path(__file__).parent\nhtml_path = OUT_DIR / f\"plot-{THEME}.html\"\npng_path = OUT_DIR / f\"plot-{THEME}.png\"\n\n# Save interactive HTML\noutput_file(str(html_path))\nsave(plot, resources=CDN, title=\"MA Plot for Differential Expression\")\n\n# Screenshot via headless Chrome (export_png unavailable on this host)\n# CDP viewport override ensures the screenshot is exactly W×H pixels, regardless\n# of any virtual browser-chrome offset that --window-size alone doesn't eliminate.\nW, H = 3200, 1800\nopts = Options()\nfor arg in (\n    \"--headless=new\",\n    \"--no-sandbox\",\n    \"--disable-dev-shm-usage\",\n    \"--disable-gpu\",\n    f\"--window-size={W},{H}\",\n    \"--hide-scrollbars\",\n    \"--force-device-scale-factor=1\",\n):\n    opts.add_argument(arg)\ndriver = webdriver.Chrome(options=opts)\ndriver.set_window_size(W, H)\ndriver.execute_cdp_cmd(\n    \"Emulation.setDeviceMetricsOverride\", {\"width\": W, \"height\": H, \"deviceScaleFactor\": 1, \"mobile\": False}\n)\ndriver.get(f\"file://{html_path.resolve()}\")\ntime.sleep(3)\ndriver.save_screenshot(str(png_path))\ndriver.quit()\n"}