#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
复现 Nature Microbiology 2026 (doi:10.1038/s41564-026-02303-9) Fig. 6a:
慢性期基孔肯雅病毒感染小鼠经小分子抗病毒药 SRI-42718 治疗 7 天后，
踝关节组织 CHIKV 基因组拷贝数（RT-qPCR，log10 转换）。
数据: chikv_genomes_fig6a.csv —— 论文补充材料 Source Data (MOESM10, sheet
'Fig 6A-Viral RNA Copies') 的 16 个真实样本原始值，非模拟。
"""
from pathlib import Path
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from matplotlib.font_manager import FontProperties, fontManager
import numpy as np
import pandas as pd
from scipy import stats

# ---------- 字体 ----------
HERE = Path(__file__).resolve().parent
FONT = HERE / "NotoSansSC-Regular.ttf"
fontManager.addfont(str(FONT))
fp = FontProperties(fname=str(FONT))

# ---------- 数据 ----------
df = pd.read_csv(HERE / "chikv_genomes_fig6a.csv")
df["log10_copies"] = np.log10(df["copies_per_ug_rna"])
groups = ["Vehicle", "SRI-42718"]
vals = {g: df.loc[df.group == g, "log10_copies"].to_numpy() for g in groups}
means = {g: v.mean() for g, v in vals.items()}
p = stats.ttest_ind(vals["Vehicle"], vals["SRI-42718"]).pvalue
print(f"Vehicle mean={means['Vehicle']:.4f}, SRI-42718 mean={means['SRI-42718']:.4f}, t-test p={p:.3e}")

# ---------- 作图 ----------
BLACK, WHITE = "#000000", "#FFFFFF"
rng = np.random.default_rng(42)

fig, ax = plt.subplots(figsize=(4.0, 5.2), dpi=300)
markers = {"Vehicle": "o", "SRI-42718": "s"}

for i, g in enumerate(groups):
    x = i + rng.uniform(-0.16, 0.16, size=len(vals[g]))
    ax.scatter(x, vals[g], s=64, marker=markers[g],
               facecolor=WHITE, edgecolor=BLACK, linewidths=1.2, zorder=3)
    ax.hlines(means[g], i - 0.28, i + 0.28, colors=BLACK, linewidths=1.6, zorder=4)

# 显著性括号
y_br, y_txt = 6.92, 7.02
ax.plot([0.12, 0.12, 0.88, 0.88], [y_br - 0.06, y_br, y_br, y_br - 0.06],
        color=BLACK, linewidth=1.2, clip_on=False)
ax.text(0.5, y_txt, "****", ha="center", va="bottom", fontsize=13,
        fontproperties=fp, color=BLACK)

ax.set_xlim(-0.55, 1.55)
ax.set_ylim(4, 7.15)
ax.set_xticks([])
ax.set_yticks([4, 5, 6, 7])
ax.set_ylabel("CHIKV genomes per \u00b5g RNA\n(log$_{10}$)",
              fontproperties=fp, fontsize=12, color=BLACK)
ax.tick_params(axis="y", labelsize=11, colors=BLACK)
for lab in ax.get_yticklabels():
    lab.set_fontproperties(fp)

ax.set_title("CHIKV genomes", fontproperties=fp, fontsize=14, pad=34, color=BLACK)

# 图例(顶部): ○ Vehicle   □ SRI-42718
from matplotlib.lines import Line2D
handles = [
    Line2D([0], [0], marker="o", color="none", markerfacecolor=WHITE,
           markeredgecolor=BLACK, markersize=9, label="Vehicle"),
    Line2D([0], [0], marker="s", color="none", markerfacecolor=WHITE,
           markeredgecolor=BLACK, markersize=9, label="SRI-42718"),
]
leg = ax.legend(handles=handles, loc="upper center", bbox_to_anchor=(0.5, 1.10),
                ncol=2, frameon=False, fontsize=11, handletextpad=0.4,
                columnspacing=1.2)
for t in leg.get_texts():
    t.set_fontproperties(fp)

# 左下角面板标号
fig.text(0.02, 0.955, "a", fontsize=20, fontweight="bold",
         fontproperties=fp, color=BLACK)

for spine in ("top", "right"):
    ax.spines[spine].set_visible(False)
for spine in ("left", "bottom"):
    ax.spines[spine].set_color(BLACK)
    ax.spines[spine].set_linewidth(1.0)

fig.tight_layout(rect=[0.02, 0, 1, 0.94])
fig.savefig(HERE / "repro_python.png", dpi=300)
print("saved repro_python.png")