ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

基于Keras与BoxCox的糖尿病遗传风险预测实战

基于Keras与BoxCox的糖尿病遗传风险预测实战 简介这份资源面向医疗健康数据分析与机器学习入门者提供一套基于天池大数据竞赛的糖尿病遗传风险预测完整方案用于糖尿病早期筛查与风险评估的辅助研究。项目以Keras框架搭建神经网络模型并采用BoxCox变换改善数据分布特性提升预测的准确性与稳定性同时配套数据可视化分析工具以图表、曲线、热图等形式直观呈现分析结果。压缩包共159个文件约21.09MB包含34个Python脚本、84张可视化图片、15份PDF文档、8个CSV数据集及ipynb笔记本、模型文件等覆盖数据处理、模型训练与结果展示各环节。已有44人学习。读者可从中获取赛题数据预处理脚本、神经网络实现代码、训练好的模型以及设计原理与使用说明文档便于理解完整建模流程并迁移到同类医疗预测任务中。1. 糖尿病遗传风险预测从竞赛数据到可落地的筛查辅助工具糖尿病早期筛查的痛点在于空腹血糖和糖化血红蛋白往往在胰岛β细胞功能已明显受损时才亮红灯而遗传风险是更前置的信号。天池大数据竞赛里有一道经典赛题要求参赛者基于人口统计学、生活方式和家族史等字段预测个体未来患糖尿病的概率。这个标题指向的正是把竞赛方案工程化用 Keras 搭建前馈神经网络配合 BoxCox 变换处理偏态特征再通过可视化把模型输出翻译成医生和患者都能理解的风险分层。它适合三类人想入门结构化数据建模的算法新手、需要给体检机构做风险评分原型的工程师、以及准备打天池类竞赛但不知道如何把 baseline 推到可用水平的选手。整条链路不依赖图像或文本一台普通笔记本就能跑通但特征工程和阈值选择里的细节决定了模型是玩具还是工具。2. 数据到手先别急着喂网络BoxCox 与特征分箱的取舍2.1 为什么连续变量直接标准化会翻车竞赛数据里像年龄、BMI、腰围、收缩压这类连续字段分布往往右偏。如果直接做 Z-score 标准化均值和方差会被长尾拉偏神经网络在反向传播时对极端值的梯度更新会异常放大表现为损失曲线震荡、验证集 AUC 卡在 0.7 上不去。BoxCox 变换的核心是引入一个幂参数 λ把非正态分布往正态拉。当 λ0 时退化为对数变换λ1 时近似恒等变换λ0.5 时是平方根变换。我一般先用scipy.stats.boxcox_normmax对每个连续特征单独估计 λ再统一变换而不是拍脑袋选 0.25 或 0.3。import numpy as np from scipy.stats import boxcox_normmax, boxcox import pandas as pd def apply_boxcox(df, cols): 对指定连续列做 BoxCox 变换返回变换后的 DataFrame 和 lambda 字典 df_out df.copy() lambdas {} for col in cols: # 只对正值列做 BoxCox负值或零值需要先平移 min_val df_out[col].min() shift 0 if min_val 0: shift abs(min_val) 1e-6 df_out[col] df_out[col] shift # 估计最优 lambda lam boxcox_normmax(df_out[col].dropna(), methodmle) lambdas[col] {lambda: lam, shift: shift} # 应用变换 df_out[col], _ boxcox(df_out[col].dropna(), lmbdalam) return df_out, lambdas # 示例对年龄、BMI、腰围做变换 continuous_cols [age, bmi, waist, sbp, dbp] df_transformed, lambda_dict apply_boxcox(df_raw, continuous_cols) print(lambda_dict)这段代码里boxcox_normmax用最大似然估计找 λshift是为了处理零值或负值——BoxCox 要求输入严格为正。变换后每个特征的偏度会显著下降用df_transformed[col].skew()验证一般能从 1.5 以上降到 0.3 以内。注意 λ 是在训练集上估计的验证集和测试集必须复用同一个 λ 和 shift否则会引入数据泄露。2.2 分箱还是连续遗传风险字段的特殊处理家族史字段在竞赛数据里通常是「父母是否糖尿病」「兄弟姐妹是否糖尿病」这类二值或三值变量。有人喜欢做独热编码有人直接映射成 0/1/2 的序数。我的经验是如果树模型为主序数编码够用但前馈神经网络对序数关系不敏感独热编码能让输入层权重更稳定。对于「糖尿病家族史数量」这种计数特征可以保留连续形式但要做截断——比如超过 3 个亲属患病统一归为 3避免长尾样本把嵌入层带偏。# 家族史特征处理 family_cols [fam_history_parents, fam_history_siblings] for col in family_cols: # 独热编码drop_first 避免共线性 dummies pd.get_dummies(df_transformed[col], prefixcol, drop_firstTrue) df_transformed pd.concat([df_transformed, dummies], axis1) df_transformed.drop(col, axis1, inplaceTrue) # 家族史总数截断 df_transformed[fam_total] df_transformed[[fam_history_parents_1, fam_history_siblings_1]].sum(axis1) df_transformed[fam_total] df_transformed[fam_total].clip(upper3)独热编码后维度会增加但家族史字段本身基数小不会造成维度爆炸。clip(upper3)是防止极端值影响嵌入层——如果你的网络第一层是 Embedding截断能减少需要学习的嵌入向量数量。2.3 缺失值不是填个均值就完事竞赛数据里血糖、胰岛素等字段常有缺失。直接填均值会低估方差填中位数会扭曲分布形状。我一般分两步先看缺失比例超过 40% 的字段直接丢弃低于 40% 的用 KNNImputer 基于其他特征做插补但插补前必须把数据分成训练/验证集只在训练集上 fit避免验证集信息泄露。from sklearn.impute import KNNImputer from sklearn.model_selection import train_test_split # 先划分数据集 X_train, X_val, y_train, y_val train_test_split( df_transformed.drop(label, axis1), df_transformed[label], test_size0.2, stratifydf_transformed[label], random_state42 ) # 只在训练集上 fit imputer imputer KNNImputer(n_neighbors5, weightsuniform) X_train_imputed imputer.fit_transform(X_train) X_val_imputed imputer.transform(X_val)n_neighbors5是经验值样本量小于 5000 时可以降到 3避免过度平滑。weightsuniform表示等权如果特征量纲差异大先用 BoxCox 变换再插补效果更好。插补后建议用missingno库画一张缺失矩阵图确认没有整列被错误填充。3. Keras 前馈网络搭建从输入层到风险概率输出3.1 网络结构不是越深越好结构化数据上前馈神经网络超过 5 层后性能通常不升反降。我用的 baseline 是三层全连接输入层维度等于特征数两个隐藏层分别 64 和 32 个神经元激活函数用 ReLU输出层 1 个神经元配 Sigmoid。Dropout 加在隐藏层之后rate 设 0.3 到 0.5 之间。BatchNormalization 放在激活函数之前能加速收敛但对小批量数据batch_size 32反而会引入噪声这时候可以去掉。import tensorflow as tf from tensorflow.keras import layers, models, callbacks def build_model(input_dim): model models.Sequential([ layers.Input(shape(input_dim,)), layers.Dense(64, kernel_initializerhe_normal), layers.BatchNormalization(), layers.Activation(relu), layers.Dropout(0.4), layers.Dense(32, kernel_initializerhe_normal), layers.BatchNormalization(), layers.Activation(relu), layers.Dropout(0.3), layers.Dense(1, activationsigmoid) ]) return model model build_model(X_train_imputed.shape[1]) model.compile( optimizertf.keras.optimizers.Adam(learning_rate1e-3), lossbinary_crossentropy, metrics[AUC, accuracy] ) model.summary()he_normal初始化适合 ReLU 激活能缓解梯度消失。BatchNormalization 放在 Dense 之后、激活之前是标准做法。Dropout rate 第一个隐藏层设 0.4第二个设 0.3因为越靠近输出层丢弃太多信息对最终概率影响越大。Adam 学习率 1e-3 是起点如果损失震荡就降到 5e-4。3.2 早停与学习率衰减两个必须加的后悔药没有早停的模型训练就像没有刹车的车。我一般监控验证集 AUCpatience 设 10 到 15 个 epoch如果 15 轮内 AUC 没有提升就停止并恢复最佳权重。学习率衰减用 ReduceLROnPlateaufactor 设 0.5patience 设 5min_lr 设 1e-6。early_stop callbacks.EarlyStopping( monitorval_auc, patience15, modemax, restore_best_weightsTrue, verbose1 ) reduce_lr callbacks.ReduceLROnPlateau( monitorval_loss, factor0.5, patience5, min_lr1e-6, verbose1 ) history model.fit( X_train_imputed, y_train, validation_data(X_val_imputed, y_val), epochs200, batch_size64, callbacks[early_stop, reduce_lr], class_weight{0: 1, 1: 3}, # 正样本加权 verbose1 )class_weight是处理类别不平衡的关键。糖尿病阳性样本通常只占 10% 到 20%如果不加权模型会倾向于预测多数类AUC 看起来还行但召回率极低。权重设为 1:3 是经验值具体可以根据正负样本比例调整公式是n_negative / n_positive。batch_size64在几千条样本的数据集上比较稳太小会导致 BatchNorm 统计量不准太大则收敛慢。3.3 阈值移动0.5 不是金标准Sigmoid 输出的是概率默认 0.5 作为分类阈值在筛查场景下往往不合适。糖尿病筛查追求高召回率宁可误报也不能漏报。我一般用验证集画 ROC 曲线找到 Youden 指数最大的点作为阈值或者直接指定召回率不低于 0.85 时的最大阈值。from sklearn.metrics import roc_curve, recall_score y_pred_proba model.predict(X_val_imputed).ravel() fpr, tpr, thresholds roc_curve(y_val, y_pred_proba) youden_idx np.argmax(tpr - fpr) best_threshold thresholds[youden_idx] # 或者指定召回率 target_recall 0.85 recalls [recall_score(y_val, (y_pred_proba t).astype(int)) for t in thresholds] valid_idx [i for i, r in enumerate(recalls) if r target_recall] if valid_idx: best_threshold thresholds[max(valid_idx)] y_pred (y_pred_proba best_threshold).astype(int) print(f最佳阈值: {best_threshold:.4f}, 召回率: {recall_score(y_val, y_pred):.4f})阈值移动后准确率可能会下降但召回率提升对筛查场景更有价值。这个阈值要保存下来推理时直接用不能每次重新算。4. 可视化分析让风险评分能被非技术人员看懂4.1 特征重要性排序的三种画法神经网络不像随机森林有现成的 feature_importances_但可以用 permutation importance 或 SHAP 值来近似。Permutation importance 更直观打乱某一列的值看验证集 AUC 下降多少下降越多说明该特征越重要。from sklearn.inspection import permutation_importance import matplotlib.pyplot as plt def perm_importance(model, X, y, feature_names, n_repeats10): 基于验证集 AUC 的排列重要性 baseline_auc model.evaluate(X, y, verbose0)[1] importances [] for i in range(X.shape[1]): scores [] for _ in range(n_repeats): X_perm X.copy() np.random.shuffle(X_perm[:, i]) auc model.evaluate(X_perm, y, verbose0)[1] scores.append(baseline_auc - auc) importances.append(np.mean(scores)) # 排序并画图 idx np.argsort(importances)[::-1] plt.figure(figsize(10, 6)) plt.barh(range(len(idx)), [importances[i] for i in idx]) plt.yticks(range(len(idx)), [feature_names[i] for i in idx]) plt.xlabel(AUC 下降幅度) plt.title(特征重要性排序Permutation Importance) plt.gca().invert_yaxis() plt.tight_layout() plt.savefig(feature_importance.png, dpi150) return importances feature_names X_train.columns.tolist() importances perm_importance(model, X_val_imputed, y_val.values, feature_names)n_repeats10是平衡计算量和稳定性的经验值样本量小的时候可以加到 20。画图时用横向条形图特征名放 y 轴避免文字重叠。保存 dpi 设 150 以上方便放进报告。4.2 风险分层把概率切成可行动的区间医生不需要知道具体概率是 0.73 还是 0.78他们需要知道「高风险、中风险、低风险」。我一般按验证集概率分布的三分位数切低于 33% 分位为低风险33% 到 66% 为中风险高于 66% 为高风险。但更严谨的做法是按临床可接受的召回率反推阈值再结合业务定分层。# 基于验证集概率分布分层 low_thresh np.percentile(y_pred_proba, 33) high_thresh np.percentile(y_pred_proba, 66) def risk_stratify(proba): if proba low_thresh: return 低风险 elif proba high_thresh: return 中风险 else: return 高风险 df_val_result pd.DataFrame({ 真实标签: y_val.values, 预测概率: y_pred_proba, 风险分层: [risk_stratify(p) for p in y_pred_proba] }) # 交叉表看分层效果 print(pd.crosstab(df_val_result[风险分层], df_val_result[真实标签], normalizeindex))交叉表能看出每个风险层的阳性率。理想情况下高风险层的阳性率应该显著高于低风险层。如果中风险层阳性率和低风险层差不多说明分层阈值需要调整可以考虑用决策树或分位数回归来切。4.3 训练过程可视化损失和 AUC 双轴图训练历史里藏着过拟合的信号。我习惯把训练集和验证集的 loss 画在同一张图AUC 画在另一张双轴对比。如果训练 loss 持续下降但验证 loss 在某个 epoch 后抬头就是过拟合的典型表现需要加 Dropout 或减层。fig, (ax1, ax2) plt.subplots(1, 2, figsize(14, 5)) # Loss 曲线 ax1.plot(history.history[loss], label训练 Loss) ax1.plot(history.history[val_loss], label验证 Loss) ax1.set_xlabel(Epoch) ax1.set_ylabel(Loss) ax1.legend() ax1.set_title(损失曲线) # AUC 曲线 ax2.plot(history.history[auc], label训练 AUC) ax2.plot(history.history[val_auc], label验证 AUC) ax2.set_xlabel(Epoch) ax2.set_ylabel(AUC) ax2.legend() ax2.set_title(AUC 曲线) plt.tight_layout() plt.savefig(training_history.png, dpi150)如果验证 AUC 在 50 个 epoch 后还在缓慢上升但验证 loss 已经开始波动说明模型在过拟合边缘可以适当增加 Dropout 或加 L2 正则。如果训练 AUC 和验证 AUC 差距超过 0.1说明模型容量过大需要减层或减神经元。5. 避坑与排查竞赛方案落地时的五个血泪教训5.1 现象验证集 AUC 0.85测试集只有 0.65原因BoxCox 的 λ 在全体数据上估计或者 KNNImputer 在划分数据集之前 fit导致验证集信息泄露到训练过程。解决所有预处理步骤必须在训练集上 fit验证集和测试集只做 transform。用sklearn.pipeline.Pipeline把变换和模型串起来能从根本上避免这类错误。5.2 现象模型训练到 30 个 epoch 后 loss 变成 NaN原因学习率太大或者 BoxCox 变换后某些特征方差极小梯度爆炸。解决先检查变换后特征的方差如果小于 1e-6说明该特征几乎是常数直接丢弃。然后把学习率降到 1e-4加梯度裁剪clipnorm1.0。5.3 现象高风险层阳性率只有 30%低风险层也有 20%原因风险分层阈值用的是验证集概率分位数但验证集和测试集的概率分布可能不一致。解决分层阈值应该在测试集上重新按分位数切或者用等距分箱代替等频分箱。更稳的做法是训练一个浅层决策树用概率和几个关键特征一起做分层规则。5.4 现象Keras 模型保存后加载预测结果和保存前不一致原因自定义层或自定义损失函数没有注册或者 BatchNormalization 的移动均值和方差没有正确保存。解决保存时用model.save(model.h5)而不是只保存权重加载时用tf.keras.models.load_model。如果用了自定义组件在 load_model 时传custom_objects参数。5.5 现象 permutation importance 跑一次要半小时原因n_repeats设太大或者验证集样本量太大。解决把验证集采样到 1000 条以内再算重要性n_repeats降到 5。如果还是慢改用 SHAP 的DeepExplainer它对神经网络有优化但需要额外安装 shap 库。6. 把模型变成筛查工具阈值校准与增量更新模型训练完只是半成品真正落地还要解决两个问题概率校准和增量更新。神经网络输出的概率往往不是真实概率Sigmoid 输出 0.8 不代表 80% 的阳性率。我一般用 Platt Scaling 或 Isotonic Regression 做校准在验证集上拟合一个映射函数把模型输出映射到校准后的概率。from sklearn.isotonic import IsotonicRegression from sklearn.calibration import calibration_curve # 在验证集上拟合校准器 calibrator IsotonicRegression(out_of_boundsclip) calibrator.fit(y_pred_proba, y_val.values) # 校准后的概率 y_pred_calibrated calibrator.predict(y_pred_proba) # 画校准曲线 fraction_pos, mean_pred calibration_curve(y_val.values, y_pred_calibrated, n_bins10) plt.plot(mean_pred, fraction_pos, s-, label校准后) plt.plot([0, 1], [0, 1], --, label理想校准) plt.xlabel(预测概率) plt.ylabel(实际阳性率) plt.legend() plt.savefig(calibration_curve.png, dpi150)校准后概率的绝对值才有意义。比如校准后输出 0.3意味着这个人群里大约 30% 是阳性。这对医生解释风险至关重要。增量更新方面我习惯每季度用新数据重新 fit BoxCox 的 λ 和 KNNImputer但模型权重用 warm_start 继续训练而不是从头训。Keras 里可以用model.fit时传initial_epoch或者直接加载旧权重再训几个 epoch。学习率要调低到 1e-4避免新数据把旧知识冲掉。最后说一个我踩过的坑有一次我把阈值定在 0.5模型在验证集上召回率 0.9上线后医生反馈漏了好几个高风险。后来发现验证集是随机划分的但实际筛查人群的年龄分布更偏大特征分布漂移了。从那以后我每次上线前都会用最近一个月的数据做一次分布对比用 KS 检验看关键特征有没有显著偏移。如果偏移超过 0.1就重新校准阈值。这个习惯帮我省了很多后悔药。希望帮到你。本文还有配套的精品资源点击获取
返回列表