好的,我把“路线A:分类建模择时”的代码按功能分段给出,并为每段配上简明解释。你可以将这些段落依次粘贴到已完成清洗后的 df
变量之后直接运行。
0. 依赖导入(一次即可)
作用:导入所需库;后续各段使用。
import numpy as np
import pandas as pd
from sklearn.model_selection import StratifiedKFold, GridSearchCV
from sklearn.preprocessing import PolynomialFeatures, StandardScaler
from sklearn.pipeline import Pipeline
from sklearn.linear_model import LogisticRegression
from sklearn.calibration import CalibratedClassifierCV
from sklearn.metrics import roc_auc_score, brier_score_loss
import matplotlib.pyplot as plt
1. 构造标签与选择特征
作用:
- 将浓度阈值0.04转为“是否达标”的二分类标签。
- 选择用于建模的特征列(存在才会使用)。
df_cls = df.copy()
df_cls['达标'] = (df_cls['Y染色体浓度'] >= 0.04).astype(int)
feat_all = ['检测孕周','孕妇BMI','年龄','身高','体重']
features = [c for c in feat_all if c in df_cls.columns]
X = df_cls[features]
y = df_cls['达标']
2. 模型训练:多项式特征 + 标准化 + 逻辑回归 + 概率校准
作用:
- 用多项式/交互项捕捉非线性。
- 标准化提升逻辑回归稳定性。
- 网格搜索选择最优多项式阶数与C。
- 概率校准(isotonic)提升预测概率的可信度。
pipe = Pipeline([('poly', PolynomialFeatures(degree=2, include_bias=False)),('scaler', StandardScaler()),('clf', LogisticRegression(max_iter=1000, solver='lbfgs', class_weight='balanced'))
])
param = {'poly__degree': [1, 2], 'clf__C': [0.5, 1.0, 5.0]}
cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)
grid = GridSearchCV(pipe, param, cv=cv, scoring='roc_auc', n_jobs=-1)
grid.fit(X, y)calib = CalibratedClassifierCV(grid.best_estimator_, cv=5, method='isotonic')
calib.fit(X, y)p_hat = calib.predict_proba(X)[:,1]
auc = roc_auc_score(y, p_hat)
brier = brier_score_loss(y, p_hat)
print(f'AUC={auc:.3f}, Brier={brier:.4f}, BestParams={grid.best_params_}')
3. BMI 分组
作用:
- 以医学常用区间对 BMI 分组,便于组内择时与解释。
bins = [0,28,32,36,40,np.inf]
labels = ['[0,28)','[28,32)','[32,36)','[36,40)','40+']
df_cls['BMI组'] = pd.cut(df_cls['孕妇BMI'], bins=bins, labels=labels, right=False)
4. 概率阈值法择时(核心)
作用:
- 在每个 BMI 组,用“组内中位特征”构造代表个体。
- 在孕周网格 t∈[10,25] 上预测达标概率曲线 P(t)。
- 取最早使 P(t)≥q(如0.9)的孕周为推荐时点。
weeks = np.arange(10, 25.01, 0.1)
q = 0.9rows_q = []
for g in df_cls['BMI组'].dropna().unique():sub = df_cls[df_cls['BMI组']==g]if sub.empty: continuemed = sub[features].median(numeric_only=True)gridX = pd.DataFrame({'检测孕周': weeks})for col in features:if col != '检测孕周':gridX[col] = med.get(col, sub[col].median())p_curve = calib.predict_proba(gridX[features])[:,1]idx = np.where(p_curve >= q)[0]best_week = weeks[idx[0]] if len(idx)>0 else np.nanrows_q.append({'BMI组': str(g), f'推荐孕周_阈值q={q}': best_week})res_q = pd.DataFrame(rows_q).sort_values('BMI组')
print('\n按BMI组的推荐孕周(基于概率阈值法)')
print(res_q)
5. 风险最小化择时(可与阈值法交叉验证)
作用:
- 定义简单风险 R(t)=w1·(1−P(t)) + w2·早检惩罚(max(12−t,0))。
- 找每组使风险最小的孕周,作为备选/对照的推荐时点。
w1, w2 = 1.0, 0.1
rows_risk = []
for g in df_cls['BMI组'].dropna().unique():sub = df_cls[df_cls['BMI组']==g]if sub.empty:continuemed = sub[features].median(numeric_only=True)gridX = pd.DataFrame({'检测孕周': weeks})for col in features:if col != '检测孕周':gridX[col] = med.get(col, sub[col].median())p_curve = calib.predict_proba(gridX[features])[:,1]early_pen = np.maximum(12 - weeks, 0)risk = w1*(1 - p_curve) + w2*early_pent_star = weeks[np.argmin(risk)]rows_risk.append({'BMI组': str(g), '最小风险孕周': t_star})res_risk = pd.DataFrame(rows_risk).sort_values('BMI组')
print('\n按BMI组的最小风险孕周(简单风险函数)')
print(res_risk)
6. 概率曲线可视化
作用:
- 展示各 BMI 组的 P(t) 曲线与阈值线,便于说明推荐时点的依据。
plt.figure(figsize=(10,6))
for g in df_cls['BMI组'].dropna().unique():sub = df_cls[df_cls['BMI组']==g]if sub.empty:continuemed = sub[features].median(numeric_only=True)gridX = pd.DataFrame({'检测孕周': weeks})for col in features:if col != '检测孕周':gridX[col] = med.get(col, sub[col].median())p_curve = calib.predict_proba(gridX[features])[:,1]plt.plot(weeks, p_curve, label=str(g))
plt.axhline(q, color='gray', linestyle='--', label=f'阈值q={q}')
plt.xlabel('检测孕周'); plt.ylabel('达标概率P'); plt.title('各BMI组的达标概率曲线P(t)')
plt.legend(); plt.tight_layout(); plt.show()
7. 简易稳健性(误差敏感性)分析
作用:
- 对阈值0.04做随机抖动,重复计算推荐孕周。
- 给出区间估计,体现方案对误差的稳健性。
np.random.seed(42)
B = 200
sigma_thr = 0.004 # 阈值抖动(可用重复检测估计更精确的σ)
rows_ci = []for g in df_cls['BMI组'].dropna().unique():sub = df_cls[df_cls['BMI组']==g]if sub.empty:continuemed = sub[features].median(numeric_only=True)t_list = []for _ in range(B):gridX = pd.DataFrame({'检测孕周': weeks})for col in features:if col != '检测孕周':gridX[col] = med.get(col, sub[col].median())p_curve = calib.predict_proba(gridX[features])[:,1]thr = np.random.normal(0.04, sigma_thr)adj = np.clip(p_curve - 0.5*((max(thr-0.04,0))/0.04), 0, 1) # 保守调整idx = np.where(adj >= q)[0]t_list.append(weeks[idx[0]] if len(idx)>0 else np.nan)arr = np.array([t for t in t_list if np.isfinite(t)])if len(arr)>0:rows_ci.append({'BMI组': str(g),'推荐孕周_中位数': np.median(arr),'2.5%': np.percentile(arr, 2.5),'97.5%': np.percentile(arr, 97.5),'有效样本': len(arr)})res_ci = pd.DataFrame(rows_ci).sort_values('BMI组')
print('\n稳健性分析:推荐孕周的模拟区间')
print(res_ci)
小结
- 第1–2段得到“可信的达标概率模型”。
- 第3–4段按 BMI 分组并给出“最早 P≥q 的推荐孕周”。
- 第5段提供“最小风险孕周”作为对照,增强决策稳健性。
- 第6–7段用图与区间说明“为什么选这个时点”和“对误差是否稳健”。
需要的话,我可以把以上分段自动合并进你当前的 game.py
并一次性运行,输出结果表与图。