# -*- coding: utf-8 -*-
"""
教程示例图生成脚本 Part B —— 数据分析示例图（模拟数据）
图06 火山图 | 图07 热图 | 图08 PCA | 图09 GO/KEGG 富集气泡图
图10 GSEA | 图11 PPI 网络 | 图12 KM 生存曲线 | 图13 Cox 森林图 | 图14 ROC/诺模图

重要说明：所有数据均为随机模拟（seed 固定），仅供教学演示图形样式与解读方法，
不代表任何真实生物学结论。读者跑真实数据时用 R 脚本（code/）生成自己的图。
运行：python gen_figures_partB_analysis.py （输出到 ../../figures/）
"""
import os
import numpy as np
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches
from matplotlib.patches import FancyBboxPatch
from scipy import stats

plt.rcParams["font.sans-serif"] = ["Microsoft YaHei", "SimHei"]
plt.rcParams["axes.unicode_minus"] = False

rng = np.random.default_rng(42)
OUT = os.path.join(os.path.dirname(__file__), "..", "..", "figures")
os.makedirs(OUT, exist_ok=True)

C_UP = "#C44E52"; C_DOWN = "#4C72B0"; C_NS = "#BBBBBB"
C_A = "#4C72B0"; C_B = "#DD8452"


def watermark(ax):
    ax.text(0.99, 0.01, "教学示例图（模拟数据）", transform=ax.transAxes, ha="right",
            va="bottom", fontsize=7.5, color="#AAAAAA", style="italic")


def fig06_volcano():
    """火山图：模拟 20000 基因差异表达"""
    n = 20000
    logfc = rng.normal(0, 0.45, n)
    pval = rng.uniform(0, 1, n)
    # 注入显著基因
    sig_up = rng.choice(n, 260, replace=False); sig_dn = rng.choice(n, 220, replace=False)
    logfc[sig_up] = rng.uniform(1.1, 3.4, len(sig_up))
    logfc[sig_dn] = rng.uniform(-3.4, -1.1, len(sig_dn))
    pval[sig_up] = 10 ** rng.uniform(-8, -1.4, len(sig_up))
    pval[sig_dn] = 10 ** rng.uniform(-8, -1.4, len(sig_dn))
    neglog = -np.log10(pval)
    up = (logfc > 1) & (neglog > 1.301)
    dn = (logfc < -1) & (neglog > 1.301)
    fig, ax = plt.subplots(figsize=(8.2, 6.6))
    ax.scatter(logfc[~(up | dn)], neglog[~(up | dn)], s=4, c=C_NS, alpha=0.45, label=f"无显著差异 (n={n - up.sum() - dn.sum()})")
    ax.scatter(logfc[dn], neglog[dn], s=5, c=C_DOWN, alpha=0.75, label=f"下调 (n={dn.sum()})")
    ax.scatter(logfc[up], neglog[up], s=5, c=C_UP, alpha=0.75, label=f"上调 (n={up.sum()})")
    ax.axhline(1.301, color="#888888", ls="--", lw=1)
    ax.axvline(1, color="#888888", ls="--", lw=1); ax.axvline(-1, color="#888888", ls="--", lw=1)
    ax.text(2.6, 1.15, "log2FC = 1", fontsize=8, color="#666666")
    ax.text(0.02, 3.4, "adj.P < 0.05", fontsize=8, color="#666666", rotation=90)
    for gene, x, y in [("GSDMD", 2.9, 5.9), ("GSDME", 2.3, 5.3), ("IL1B", -2.6, 5.1),
                       ("CASP1", 1.7, 4.4), ("NLRP3", -1.9, 4.0)]:
        ax.annotate(gene, (x, y), fontsize=8.5, color="#333333",
                    arrowprops=dict(arrowstyle="-", color="#999999", lw=0.8))
    ax.set_xlabel("log2 (Fold Change)  肿瘤 vs 癌旁", fontsize=11)
    ax.set_ylabel("-log10 (adjusted P)", fontsize=11)
    ax.set_title("差异表达基因火山图（limma, TCGA-LIHC 示例）", fontsize=13, fontweight="bold")
    ax.legend(loc="upper left", fontsize=8.5, frameon=False)
    watermark(ax)
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig06_volcano.png"), dpi=150)
    plt.close(fig)


def fig07_heatmap():
    """热图：模拟 30 基因 × 20 样本"""
    genes = ["GSDMD", "GSDME", "CASP1", "IL1B", "IL18", "NLRP3", "AIM2", "PYCARD",
             "CASP3", "CASP4", "CASP5", "CASP6", "CASP8", "CASP9", "GSDMA", "GSDMB",
             "GSDMC", "GSDMF", "TNF", "NFKB1", "NFKB2", "TIRAP", "MEFV", "PJVK",
             "IL1A", "IL6", "IL18R1", "TNFSF10", "TNFRSF1A", "TNFRSF1B"]
    n_genes, n_s = len(genes), 20
    base = rng.normal(0, 1, (n_genes, n_s))
    for i in range(n_genes):
        if i % 3 == 0:   # 上调
            base[i, 10:] += rng.uniform(1.5, 3.0)
        elif i % 3 == 1:  # 下调
            base[i, 10:] -= rng.uniform(1.5, 3.0)
    from scipy.cluster.hierarchy import linkage, leaves_list
    from scipy.spatial.distance import pdist
    row_link = linkage(pdist(base), method="average")
    row_ord = leaves_list(row_link)
    base = base[row_ord, :]
    genes = [genes[i] for i in row_ord]
    fig, ax = plt.subplots(figsize=(9.5, 7.2))
    im = ax.imshow(base, cmap="RdBu_r", aspect="auto", vmin=-3, vmax=3)
    ax.set_xticks(range(n_s))
    ax.set_xticklabels([f"{'癌旁' if j < 10 else '肿瘤'}{j % 10 + 1}" for j in range(n_s)],
                       rotation=90, fontsize=7)
    ax.set_yticks(range(n_genes))
    ax.set_yticklabels(genes, fontsize=8.5)
    ax.axvline(9.5, color="black", lw=1.2)
    for j in range(10):
        ax.add_patch(mpatches.Rectangle((j - 0.5, -0.5), 1, n_genes, fill=False, edgecolor=C_A, lw=1.6))
    for j in range(10, 20):
        ax.add_patch(mpatches.Rectangle((j - 0.5, -0.5), 1, n_genes, fill=False, edgecolor=C_B, lw=1.6))
    ax.legend(handles=[mpatches.Patch(color=C_A, label="癌旁 (n=10)"),
                       mpatches.Patch(color=C_B, label="肿瘤 (n=10)")],
              loc="upper right", fontsize=9, frameon=True)
    cb = fig.colorbar(im, ax=ax, shrink=0.75, pad=0.02)
    cb.set_label("z-score 标准化表达", fontsize=9)
    ax.set_title("焦亡相关基因表达热图（示例，行聚类）", fontsize=13, fontweight="bold")
    watermark(ax)
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig07_heatmap.png"), dpi=150)
    plt.close(fig)


def fig08_pca():
    """PCA 图：肿瘤 vs 癌旁分离"""
    n1 = n2 = 60
    x = np.concatenate([rng.normal(0, 1, n1), rng.normal(3.2, 1.3, n2)])
    y = np.concatenate([rng.normal(0, 1, n1), rng.normal(1.1, 1.2, n2)])
    grp = np.array(["癌旁"] * n1 + ["肿瘤"] * n2)
    fig, ax = plt.subplots(figsize=(7.6, 6.4))
    ax.scatter(x[grp == "癌旁"], y[grp == "癌旁"], c=C_A, s=30, alpha=0.75, edgecolor="white", label="癌旁 (n=60)")
    ax.scatter(x[grp == "肿瘤"], y[grp == "肿瘤"], c=C_B, s=30, alpha=0.75, edgecolor="white", label="肿瘤 (n=60)")
    for label, col in [("癌旁", C_A), ("肿瘤", C_B)]:
        xs = x[grp == label]; ys = y[grp == label]
        ex, ey = xs.mean(), ys.mean()
        t = stats.t.ppf(0.975, len(xs) - 1)
        ax.add_patch(Ellipse_ci(ax, ex, ey, np.sqrt(xs.var()), np.sqrt(ys.var()), 0, col, 0.12))
    ax.set_xlabel("PC1 (34.2%)", fontsize=11)
    ax.set_ylabel("PC2 (12.8%)", fontsize=11)
    ax.set_title("主成分分析：肿瘤与癌旁样本分离（示例）", fontsize=13, fontweight="bold")
    ax.legend(frameon=False, fontsize=9.5)
    watermark(ax)
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig08_pca.png"), dpi=150)
    plt.close(fig)


def Ellipse_ci(ax, cx, cy, sx, sy, angle, color, alpha):
    from matplotlib.patches import Ellipse as E
    return E((cx, cy), width=3.0 * sx, height=3.0 * sy, angle=angle,
             facecolor=color, alpha=alpha, edgecolor=color, lw=1.4)


def fig09_go_kegg():
    """GO/KEGG 富集气泡图（水平，模拟）"""
    rows = [
        ("炎症小体复合物组装", "BP", 0.32, 18, 1e-9, C_UP),
        ("白细胞介素-1β 产生", "BP", 0.28, 15, 1e-8, C_UP),
        ("细胞因子介导的信号通路", "BP", 0.41, 26, 1e-7, C_UP),
        ("程序性细胞死亡", "BP", 0.45, 22, 1e-6, C_UP),
        ("NLRP3 炎症小体复合物", "CC", 0.18, 8, 1e-6, C_GREEN := "#55A868"),
        ("细胞膜成孔复合物", "CC", 0.15, 6, 1e-5, C_GREEN),
        ("半胱天冬酶激活", "MF", 0.12, 7, 1e-5, C_ORANGE := "#DD8452"),
        ("IL-1 受体结合", "MF", 0.10, 5, 1e-4, C_ORANGE),
        ("NOD 样受体信号通路", "KEGG", 0.22, 12, 1e-7, C_PURPLE := "#8172B3"),
        ("细胞因子-细胞因子受体互作", "KEGG", 0.26, 14, 1e-6, C_PURPLE),
        ("凋亡通路", "KEGG", 0.19, 10, 1e-5, C_PURPLE),
    ]
    names, ont, gr, cnt, pv, cols = zip(*rows)
    fig, ax = plt.subplots(figsize=(9.8, 6.6))
    ypos = np.arange(len(names))
    ax.scatter(gr, ypos, s=np.array(cnt) * 14, c=cols, alpha=0.85, edgecolor="white", linewidth=1.2)
    for i, p in enumerate(pv):
        ax.text(gr[i] + 0.012, ypos[i] + 0.18, f"p={p:.0e}", fontsize=7, color="#555555")
    ax.set_yticks(ypos)
    ax.set_yticklabels(names, fontsize=9.5)
    ax.set_xlabel("Gene Ratio", fontsize=11)
    ax.set_xlim(0.05, 0.55)
    for i, (start, end, label, color) in enumerate([(0, 4, "BP 生物过程", C_UP), (4, 6, "CC 细胞组分", C_GREEN),
                                                    (6, 8, "MF 分子功能", C_ORANGE), (8, 11, "KEGG 通路", C_PURPLE)]):
        ax.axhspan(start - 0.5, end - 0.5, color=color, alpha=0.06, zorder=0)
        ax.text(0.52, (start + end - 1) / 2, label, fontsize=7.5, color=color, va="center", ha="right")
    h = [mpatches.Patch(color=C_UP, label="BP"), mpatches.Patch(color=C_GREEN, label="CC"),
         mpatches.Patch(color=C_ORANGE, label="MF"), mpatches.Patch(color=C_PURPLE, label="KEGG")]
    ax.legend(handles=h, loc="lower right", fontsize=8.5, frameon=False, ncol=2)
    ax.set_title("GO/KEGG 功能富集分析（示例，点大小=基因数）", fontsize=13, fontweight="bold")
    watermark(ax)
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig09_go_kegg.png"), dpi=150)
    plt.close(fig)


def fig10_gsea():
    """GSEA 富集图（模拟）"""
    n = 20000
    rank = rng.normal(0, 1, n)
    rank[:400] += np.linspace(2.2, 0.6, 400)  # 基因集富集于前列
    es = np.cumsum(rank / np.abs(rank).sum()) * 1000
    fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(8.2, 6.2), sharex=True,
                                   gridspec_kw={"height_ratios": [1, 0.45], "hspace": 0.08})
    x = np.arange(len(rank))
    ax1.fill_between(x, 0, es, where=es >= 0, color=C_UP, alpha=0.75)
    ax1.fill_between(x, 0, es, where=es < 0, color=C_DOWN, alpha=0.75)
    ax1.axhline(0, color="black", lw=0.8)
    ax1.set_ylabel("Enrichment score", fontsize=10)
    ax1.set_title("GSEA：炎症小体相关基因集在肿瘤中富集（示例）", fontsize=13, fontweight="bold")
    ax1.text(0.97, 0.8, "NES = 2.31\nFDR < 0.001", transform=ax1.transAxes, ha="right",
             fontsize=9, color=C_UP, fontweight="bold")
    # 底部条带
    for i in range(n):
        c = C_UP if rank[i] > 0 else "#CCCCCC"
        ax2.add_patch(mpatches.Rectangle((i, 0), 1, 0.9, color=c, linewidth=0))
    ax2.set_yticks([])
    ax2.set_xlabel("按 log2FC 排序的基因（前=上调）", fontsize=10)
    ax2.set_ylim(0, 1)
    watermark(ax1)
    fig.savefig(os.path.join(OUT, "fig10_gsea.png"), dpi=150)
    plt.close(fig)


def fig11_ppi():
    """PPI 网络（networkx 力导向）"""
    import networkx as nx
    G = nx.Graph()
    hubs = ["GSDMD", "GSDME", "CASP1", "IL1B", "IL18", "NLRP3", "AIM2", "PYCARD",
            "CASP3", "TNF", "NFKB1", "IL6"]
    nodes = hubs + ["CASP4", "CASP5", "CASP8", "GSDMB", "IL18R1", "NFKB2", "TIRAP",
                    "MEFV", "IL1A", "TNFSF10", "TNFRSF1A", "IL18RAP", "BAX", "BCL2",
                    "GSDMC", "GSDMA", "CASP6", "CASP9", "CASP12", "PJVK"]
    G.add_nodes_from(nodes)
    rng2 = np.random.default_rng(7)
    for i, a in enumerate(nodes):
        for b in nodes[i + 1:]:
            if a in hubs and b in hubs and rng2.random() < 0.55:
                G.add_edge(a, b)
            elif rng2.random() < 0.04:
                G.add_edge(a, b)
    pos = nx.spring_layout(G, seed=11, k=0.55)
    hub_colors = {n: C_UP for n in hubs}
    other_c = "#B0C4DE"
    sizes = [520 if n in hubs else 150 for n in G.nodes]
    fig, ax = plt.subplots(figsize=(9.5, 7.4))
    nx.draw_networkx_edges(G, pos, ax=ax, edge_color="#CCCCCC", width=1.0, alpha=0.6)
    nx.draw_networkx_nodes(G, pos, ax=ax, node_size=sizes,
                           node_color=[hub_colors.get(n, other_c) for n in G.nodes],
                           edgecolors="white", linewidths=1.5, alpha=0.9)
    for n in hubs:
        d = G.degree(n)
        ax.annotate(f"{n}", (pos[n][0], pos[n][1]), fontsize=8.5, fontweight="bold", color="#8B1A1A",
                    ha="center", va="center", bbox=dict(boxstyle="round,pad=0.15", fc="white", ec=C_UP, lw=0.8))
    for n in G.nodes:
        if n not in hubs:
            ax.annotate(n, (pos[n][0], pos[n][1]), fontsize=6.5, color="#444444",
                        ha="center", va="center")
    ax.axis("off")
    ax.legend(handles=[mpatches.Patch(color=C_UP, label="hub 基因（MCC 前 12）"),
                       mpatches.Patch(color=other_c, label="其余互作蛋白")],
              loc="lower right", fontsize=9, frameon=False)
    ax.set_title("蛋白互作网络 PPI（STRING 数据 · cytoHubba MCC 排序，示例）",
                 fontsize=13, fontweight="bold")
    watermark(ax)
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig11_ppi.png"), dpi=150)
    plt.close(fig)


def fig12_km():
    """KM 生存曲线（模拟，指数分布生成）"""
    t_hi = rng.exponential(18, 100); t_lo = rng.exponential(34, 100)
    c_hi = rng.uniform(6, 60, 100); c_lo = rng.uniform(6, 60, 100)
    obs_hi = np.minimum(t_hi, c_hi); ev_hi = (t_hi <= c_hi).astype(int)
    obs_lo = np.minimum(t_lo, c_lo); ev_lo = (t_lo <= c_lo).astype(int)

    def km_est(times, events):
        order = np.argsort(times)
        ts, es = times[order], events[order]
        surv = [1.0]; st = [0.0]
        for i, t in enumerate(ts):
            n_at = (ts >= t).sum()
            d = es[i]
            if n_at > 0 and d > 0:
                surv.append(surv[-1] * (1 - d / n_at))
                st.append(t)
        return np.array(st), np.array(surv)

    st1, sv1 = km_est(obs_hi, ev_hi)
    st2, sv2 = km_est(obs_lo, ev_lo)
    fig, ax = plt.subplots(figsize=(8.2, 6.2))
    ax.step(st1, sv1, where="post", color=C_UP, lw=2.2, label="高表达组 (n=100)")
    ax.step(st2, sv2, where="post", color=C_B, lw=2.2, label="低表达组 (n=100)")
    ax.scatter(obs_hi[ev_hi == 0], np.interp(obs_hi[ev_hi == 0], st1, sv1), s=6, color=C_UP, alpha=0.6, zorder=5)
    ax.scatter(obs_lo[ev_lo == 0], np.interp(obs_lo[ev_lo == 0], st2, sv2), s=6, color=C_B, alpha=0.6, zorder=5)
    ax.set_xlabel("生存时间（月）", fontsize=11)
    ax.set_ylabel("总生存率", fontsize=11)
    ax.set_ylim(0, 1.02)
    ax.set_title("GSDMD 高/低表达组 Kaplan-Meier 生存曲线（示例）", fontsize=13, fontweight="bold")
    ax.text(0.97, 0.12, "log-rank P < 0.001\nHR = 2.31 (95%CI 1.52–3.51)",
            transform=ax.transAxes, ha="right", fontsize=10, color=C_UP, fontweight="bold")
    ax.legend(frameon=False, fontsize=9.5, loc="upper right")
    watermark(ax)
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig12_km.png"), dpi=150)
    plt.close(fig)


def fig13_forest():
    """Cox 森林图（模拟多因素结果）"""
    vars_ = ["年龄 (≥60 vs <60)", "性别 (男 vs 女)", "TNM 分期 (III/IV vs I/II)",
             "GSDMD 高表达", "IL1B 高表达", "CASP1 高表达"]
    hr = [1.62, 1.05, 2.44, 2.15, 1.38, 1.71]
    lo = [1.08, 0.72, 1.63, 1.44, 0.94, 1.19]
    hi = [2.43, 1.53, 3.66, 3.21, 2.02, 2.46]
    pv = [0.019, 0.778, "<0.001", "<0.001", 0.102, 0.004]
    y = np.arange(len(vars_))[::-1]
    fig, ax = plt.subplots(figsize=(9.6, 5.6))
    for yi, h, l, u in zip(y, hr, lo, hi):
        ax.plot([l, u], [yi, yi], color="#555555", lw=2.2, solid_capstyle="round")
        ax.scatter(h, yi, s=55, c=C_UP if (l > 1 or u < 1) else C_NS, zorder=5, edgecolor="white", linewidth=1)
    ax.axvline(1, color="black", ls="--", lw=1.2)
    ax.set_yticks(y)
    ax.set_yticklabels(vars_, fontsize=10)
    ax.set_xscale("log")
    ax.set_xlim(0.4, 5)
    ax.set_xticks([0.5, 1, 2, 4])
    ax.set_xticklabels(["0.5", "1", "2", "4"])
    ax.set_xlabel("风险比 HR (95% CI, 对数轴)", fontsize=11)
    ax.set_title("多因素 Cox 回归森林图（示例）", fontsize=13, fontweight="bold")
    for yi, p in zip(y, pv):
        ax.text(5.0, yi, f"P = {p}", fontsize=9, va="center", color="#555555")
    ax.text(0.42, -1.05, "森林图：点=HR，横线=95%CI；跨过 1 的竖虚线表示无显著差异", fontsize=8.5, color="#666666")
    watermark(ax)
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig13_forest.png"), dpi=150)
    plt.close(fig)


def fig14_roc_nomogram():
    """ROC 曲线 + 诺模图（模拟）"""
    fig, axes = plt.subplots(1, 2, figsize=(11.6, 5.4), gridspec_kw={"width_ratios": [1, 1.4]})
    ax = axes[0]
    for label, a, col in [("1 年 AUC=0.86", 0.10, C_UP), ("3 年 AUC=0.82", 0.25, C_ORANGE := "#DD8452"),
                          ("5 年 AUC=0.79", 0.40, C_BLUE := "#4C72B0")]:
        tpr = np.clip(1 - np.exp(-np.linspace(0, 4, 100) * (1 - a) / 1.2), 0, 1)
        fpr = np.linspace(0, 1, 100)
        ax.plot(fpr, 1 - np.exp(-np.linspace(0, 3.5, 100) * (1 - a)), color=col, lw=2.2, label=label)
    ax.plot([0, 1], [0, 1], color="#BBBBBB", ls="--", lw=1.2)
    ax.set_xlabel("1 - 特异度 (FPR)", fontsize=10.5)
    ax.set_ylabel("敏感度 (TPR)", fontsize=10.5)
    ax.set_title("(a) 风险模型 ROC 曲线（示例）", fontsize=12, fontweight="bold")
    ax.legend(frameon=False, fontsize=9)
    watermark(ax)

    ax2 = axes[1]
    ax2.set_xlim(0, 10); ax2.set_ylim(0, 10); ax2.axis("off")
    ax2.set_title("(b) 诺模图 Nomogram（示例）", fontsize=12, fontweight="bold")
    scale = np.arange(0, 101, 10)
    ax2.plot([0.6, 9.4], [9.3, 9.3], color="black", lw=1.4)
    for v in scale:
        ax2.plot([0.6 + v / 100 * 8.8, 0.6 + v / 100 * 8.8], [9.15, 9.3], color="black", lw=1)
        ax2.text(0.6 + v / 100 * 8.8, 9.02, str(v), ha="center", fontsize=7)
    ax2.text(0.5, 9.5, "分数", fontsize=9, fontweight="bold")
    rows = [
        ("TNM 分期", [0, 20, 40, 60, 80, 100], "I  II  III  IV"),
        ("GSDMD 表达", [0, 30, 55, 100], "低  中  高"),
        ("CASP1 表达", [0, 25, 60, 100], "低  中  高"),
        ("IL1B 表达", [0, 35, 70, 100], "低  中  高"),
    ]
    y0 = 8.0
    for name, ticks, labels in rows:
        ax2.plot([0.6, 9.4], [y0, y0], color="black", lw=1.2)
        for v in ticks:
            x = 0.6 + v / 100 * 8.8
            ax2.plot([x, x], [y0 - 0.13, y0], color="black", lw=1)
            ax2.text(x, y0 - 0.32, str(v), ha="center", fontsize=6.8)
        ax2.text(0.5, y0, name, fontsize=8.5, fontweight="bold", va="center")
        ax2.text(9.5, y0, labels, fontsize=6.8, va="center", color="#555555")
        y0 -= 1.35
    ax2.plot([0.6, 9.4], [y0, y0], color="black", lw=1.2)
    ax2.text(0.5, y0, "总分数", fontsize=8.5, fontweight="bold", va="center")
    ax2.text(0.5, y0 - 1.1, "3 年生存概率", fontsize=8.5, fontweight="bold", va="center")
    ax2.plot([0.6, 9.4], [y0 - 0.35, y0 - 0.35], color="black", lw=1.2)
    ax2.text(1.2, y0 - 0.35, "0.95", fontsize=7, va="center")
    ax2.text(4.5, y0 - 0.35, "0.5", fontsize=7, va="center")
    ax2.text(8.6, y0 - 0.35, "0.05", fontsize=7, va="center")
    fig.tight_layout()
    fig.savefig(os.path.join(OUT, "fig14_roc_nomogram.png"), dpi=150)
    plt.close(fig)


if __name__ == "__main__":
    fig06_volcano()
    fig07_heatmap()
    fig08_pca()
    fig09_go_kegg()
    fig10_gsea()
    fig11_ppi()
    fig12_km()
    fig13_forest()
    fig14_roc_nomogram()
    print("Part B done:", sorted(os.listdir(OUT)))
