简介面向医学科研与金融风控场景的MATLAB风险预测模型工具包涵盖ROC与PR曲线、NRI重分类、AUC比较等核心评估方法适用于临床预测模型构建、疾病复发或违约风险评分等二分类任务并可处理样本不平衡下的模型性能评估。压缩包共40个文件其中28个.m文件为主要算法脚本4个xlsx/xls为示例数据与结果表3个rtf为说明文档另有txt说明与脚本备份整个压缩包仅584KB便于快速下载与部署。目前已有498人学习下载使用。使用者可直接运行Roc.m、main.m等脚本借助内置数据复现523例预实验的分析流程并通过Risk_Assessment_Plot模块输出可视化风险图表多类别NRI与AUC相关性比较工具则在新旧指标对比、模型优化中提供现成实现适合需要快速搭建预测模型评估体系、开展风险量化研究的研究者与数据分析师。1. 从523例样本到ROC曲线这套风险预测模型代码包到底解决了什么先说一个反直觉的结论这套风险预测模型代码包里最值得研究的不是Roc.m绘图脚本而是数据清洗和变量筛选那部分——因为最终ROC预测结果能不能站得住关键在523例样本进入模型之前发生了什么。这是一份用MATLAB写成的风险预测模型实现覆盖了从Excel原始表到逻辑回归建模、ROC曲线绘制、AUC置信区间与模型对比、NRI净重分类改善指数计算以及风险评估图输出的完整链路。压缩包里的origin.xlsx和test.xlsx是已经切分好的训练与验证数据。适合正在做二分类预测、手头有临床或金融表格数据、希望直接照着改并快速拿到AUC与ROC结论的研究者。2. 数据准备阶段523例样本进入模型前的三道清洗关口2.1 先把压缩包里的文件分清楚哪些是核心脚本、哪些是数据、哪些可以直接删打开风险预测模型1225.rar后我习惯先做一件事把文件分成四组。第一组是数据文件包括523例预实验.xls、origin.xlsx、test.xlsx、roc.xlsx第二组是数据清洗与变量选择的脚本包括exciseRows.m、zhen.m、choice.m第三组是建模与评估脚本包括LR_Pz_choice.m、Roc.m、ROC test.m、CIAUC.m、AUC_compare_correlated.m第四组是NRI相关脚本包括NRI.m、Category_Free_NRI.m、Category_Free_NRI_ci.m、multi_category_NRI.m、multi_category_NRI_ci.m、FindNandD.m。压缩包名字里的1225是打包日期will7jv应该是上传者标识和代码逻辑没有关系不影响运行。把分组说得更细一点。zhen.m从命名看是真实的拼音对应的是真实标签处理FindNandD.m是统计事件组与非事件组人数的前置函数为NRI计算做准备。另外压缩包根目录下有__MACOSX文件夹以及一堆以._开头的文件这些是macOS系统自动生成的元数据跟分析没有关系。还有~$开头的rtf文件是Office打开文档时留下的临时文件。这些东西在Windows上解压后不会自动消失建议解压后先手动清掉否则MATLAB扫描当前目录时会把隐藏文件也扫进去偶尔出现莫名其妙的路径报错。我把这份文件清单整理成表格方便你对照使用分组文件作用数据523例预实验.xls、origin.xlsx、test.xlsx、roc.xlsx训练/验证原始表与结果导出清洗exciseRows.m、zhen.m、choice.m行剔除、标签对齐、变量选择建模LR_Pz_choice.m逻辑回归P值逐步筛选ROCRoc.m、ROC test.m、CIAUC.m、AUC_compare_correlated.mROC绘制、置信区间、模型比较NRINRI.m、Category_Free_NRI.m、multi_category_NRI.m等重分类指标计算可视化Risk_Assessment_Plot.m风险分层图输出把这张表放在手边后面每一步要动哪个文件直接对号入座就行。2.2 行剔除与变量选择exciseRows.m和choice.m的实际用法exciseRows.m解决的是最基础的脏活按行索引剔除不符合纳入标准的样本。% exciseRows.m —— 按行索引剔除指定样本返回清洗后的数据矩阵 % 输入: data为m×n数值矩阵rowsToExclude为待剔除的行索引向量 function dataClean exciseRows(data, rowsToExclude) if isempty(rowsToExclude) dataClean data; return; end keepIdx true(size(data, 1), 1); % 初始化全量保留掩码 keepIdx(rowsToExclude) false; % 将待剔除行位置置为false dataClean data(keepIdx, :); % 用逻辑索引一次性取出保留行 end逻辑说明这里用逻辑掩码而不是data(rowsToExclude, :) []是因为后者在循环里多次删除时会导致行号整体前移后一次删除的行号已经不是原始表里的行号了。用逻辑索引直接生成一个保留哪些行的真值数组清洗完还能对照keepIdx检查到底剔除了哪几条样本出问题时有后悔药可吃。参数说明rowsToExclude可以是行号向量也可以是逻辑向量。我的经验是先从原始Excel里把不符合条件的人在Excel行号记下来但要注意读入MATLAB后行号会少一行表头占了第一行所以传参之前先执行rowsToExclude rowsToExclude - 1这种偏移再核对一遍再跑能少踩一个很隐蔽的坑。choice.m则是按变量名筛选要进入建模的列% choice.m —— 按变量名返回需要建模的列子矩阵与列索引 % 输入: X为n×p数据矩阵varNames为p个变量名的元胞数组 % selectedNames为要选出的变量名列表 function [Xsub, colIdx] choice(X, varNames, selectedNames) [~, colIdx] ismember(selectedNames, varNames); if any(colIdx 0) error(变量名不存在: %s, strjoin(selectedNames(colIdx 0), , )); end Xsub X(:, colIdx); end逻辑说明使用ismember做按名引用而不是硬编码列序号是因为当从523例预实验.xls换成自己的数据时Excel列顺序大概率会变但列名通常是稳定的。这个脚本的价值在于一旦你写错变量名它会直接报错而不是默默返回一个空矩阵这种显式失败比静默错误友好得多。参数说明varNames是从Excel表头读入的变量名元胞数组selectedNames是你在分析计划里确定要纳入的预测因子。返回值里colIdx是这些变量在原矩阵中的实际列号画相关性矩阵或做共线性检查时会经常用到。2.3 Excel数据读入阶段的三个高发坑文本型数字、缺失编码和列对齐数据读入阶段最常见的问题不是统计上的而是Excel本身带来的。第一坑是文本型数字一列血压值如果部分单元格是文本格式xlsread读进来后这一列会被切分成两列或者直接变成NaN。第二坑是缺失值编码不统一有人用空白、有人用0、有人输入NA如果缺失被编码成0后续模型会把0当成真实测量值整个预测概率分布都会偏移。第三坑是列对齐在Excel里做过vlookup或者手动粘贴的人都知道某一列与其他列的行序经常对不齐这种事用代码很难自动发现。我处理这些表时有一个固定顺序先用readtable整表读入对每一列统计NaN比例和唯一值个数然后把所有字符串列统一转成double并记录转换失败的位置最后对关键变量做散点图矩阵用肉眼确认行序没有错位。这套动作很基础但已经帮我避免过三次AUC很高后来发现标签错位的事故。记住一个原则所有高级统计指标的前提是表干净。3. 逻辑回归变量筛选LR_Pz_choice.m的P值门槛与训练/验证数据切分3.1 为什么风险预测模型优先选择逻辑回归而不是机器学习风险预测场景的结局大多是二分类疾病发生或未发生、违约或不违约。逻辑回归输出的是事件概率天然可以映射成风险评分每个变量的系数能翻译成OR值优势比配合95%置信区间这是临床论文和业务评审里最熟悉的语言。换成随机森林或XGBoost预测精度可能更高但可解释性弱很多参数调整路径也长在只有523例的小样本环境下还容易过拟合。所以我的习惯始终是先用逻辑回归建一个基线模型拿到AUC和校准度指标后再决定要不要上更复杂的模型。这套代码包里的LR_Pz_choice.m做的就是逻辑回归里最关键的一步——变量筛选。3.2 LR_Pz_choice.m的工作逻辑P值门槛怎么设、循环怎么停LR_Pz_choice.m从名字看就是Logistic Regression P value Choice核心思路是逐步回归用P值决定哪些变量进来、哪些出去。在压缩包里它和单独跑一个glmfit的区别在于它是自动化的、可复现的变量筛选过程而不是手动一个一个试变量。我用伪代码还原它的主循环方便你理解每一步在干什么% LR_Pz_choice.m —— 基于P值的逐步逻辑回归变量筛选主循环示意 % 输入: X为n×p预测因子矩阵y为n×1二分类结局(0/1) % pEnter: 变量进入模型的P值上限默认0.05 % pRemove: 变量被剔除的P值下限默认0.10 function [selected, coef, stats] LR_Pz_choice(X, varNames, y, pEnter, pRemove) if nargin 4, pEnter 0.05; end if nargin 5, pRemove 0.10; end selected []; % 当前入选变量集合 candidateIdx 1:size(X, 2); % 尚未进入模型的候选变量 while true % 第一步对每个候选变量单独做单因素logistic回归记录P值 % 第二步找出最小P值如果小于pEnter将其纳入selected % 第三步基于当前selected重新拟合多因素模型剔除P值大于pRemove的变量 % 第四步如果本轮没有任何变量进入或剔除循环终止 changed stepOnce(X, y, selected, candidateIdx, pEnter, pRemove); if ~changed break; end end % 最后对selected变量拟合一次完整模型输出系数、OR和置信区间 [coef, stats] fitFinalModel(X(:, selected), y); end逻辑说明这里的stepOnce和fitFinalModel是模拟实际代码中的两个内部步骤压缩包里对应的是一系列MATLAB统计函数调用。循环的终止条件是某轮迭代中既没有新变量进入、也没有老变量被剔除这是逐步回归的标准收敛规则。pEnter和pRemove构成一个双向门槛pEnter管进来的条件pRemove管出去的条件。用0.05和0.10这两个经典值的意图是让筛选有适度冗余避免变量在边界上反复进出。参数说明如果你的样本量不大比如事件组只有40-50例单因素回归的P值会很不稳定建议把pEnter放宽到0.10、pRemove放到0.15先让变量进来再说后面用AIC或BIC去平衡复杂度。另外要注意共线性问题不会直接体现在P值上——如果两个候选变量相关性超过0.7建议预先只保留其中一个否则可能出现系数符号异常或标准误爆炸。3.3 训练集与验证集的切分纪律origin.xlsx和test.xlsx不能互相偷看压缩包里已经预置了origin.xlsx建模数据和test.xlsx验证数据这是一个比随机拆分更干净的数据结构。实际执行时我的流程是只在origin.xlsx上做变量筛选和模型拟合test.xlsx只在最后评估阶段被读入一次。LR_Pz_choice.m、Roc.m的模型参数都只能由训练集决定。有人图省事先在全量数据上筛变量再切分训练/测试这会让验证集的AUC虚高审稿人一旦追问变量筛选细节就难以解释。还有一个容易被忽略的点跑测试集之前先对比两个数据集的结局事件率。我见过origin事件率30%、test只有8%的情况这种验证集得到的AUC参考意义很有限解读时必须在报告里说明事件率差异。正确的做法是检查两个数据集的基线特征表年龄、性别、关键预测因子的均值与分布是否大致一致再做结论。这套验证集先看分布、再算AUC的动作能帮你避免很多解释上的尴尬。跑完LR_Pz_choice.m之后一定要把最终入选变量的OR值和置信区间打印出来。我一般会在命令行留一份、导出到Excel一份。这个表里的OR(95%CI)就是风险预测模型最核心的可交付成果很多人只盯着AUC而忘了把OR表整理出来这是本末倒置。4. ROC曲线与AUC对比Roc.m、CIAUC.m到DeLong检验的完整链路4.1 Roc.m一条ROC曲线从计算到出图Roc.m是这份资源里被提到最多的脚本它的输入是预测概率和真实标签输出是AUC、置信区间以及一张可以直接放进报告里的图。% Roc.m —— 单模型ROC曲线绘制与AUC计算 % 输入: score为预测概率列向量y为真实标签(1事件, 0非事件) % 输出: auc为曲线下面积ci为95%置信区间 function [auc, ci] Roc(score, y) if length(score) ~ length(y) error(score与y长度不一致请先检查样本对齐); end % perfcurve是MATLAB统计工具箱的标准ROC函数 [fpr, tpr, ~, auc] perfcurve(y, score, 1); % AUC标准误的近似计算事件组较少时改用手动bootstrap n1 sum(y 1); n0 sum(y 0); seAUC sqrt(auc * (1 - auc) / min(n1, n0)); ci [auc - 1.96 * seAUC, auc 1.96 * seAUC]; figure(Color, w); plot(fpr, tpr, LineWidth, 1.8); hold on; plot([0 1], [0 1], k--); xlabel(1 - Specificity); ylabel(Sensitivity); title(sprintf(ROC Curve (AUC %.3f, 95%% CI %.3f - %.3f), auc, ci(1), ci(2))); legend(Model, Random, Location, SouthEast); grid on; end逻辑说明perfcurve会按所有可能的截断阈值计算真阳性率和假阳性率画出一条从右上角到左下角的阶梯曲线。AUC是所有阈值下判别能力的综合表现0.5代表随机0.7-0.8在临床模型里算中等可用0.9以上要警惕过拟合。代码里用n1和n0中的较小值来近似计算AUC标准误这个公式来自正态近似样本量大时够用。参数说明score必须是连续预测概率直接用逻辑回归输出的概率即可1表示正类事件。如果你的事件样本量少于50这个近似置信区间会偏窄建议改用bootstrap或者压缩包里的CIAUC.m。ROC test.m则是Roc.m的训练/验证封装版依次对origin.xlsx和test.xlsx调用Roc.m把两条曲线画在同一张图上对比这是论文中常见的呈现方式。4.2 CIAUC.m与AUC_compare_correlated.m置信区间与两个模型的差异检验CIAUC.m干的是把AUC的置信区间算准这件事。Roc.m里的近似公式只适用于样本量均衡的情况CIAUC.m一般会采用DeLong法它利用排序秩的协方差结构计算AUC方差比正态近似更稳定也更容易被审稿人接受。如果你的两组样本量差异很大直接用Roc.m里的近似区间会有偏差这时候CIAUC.m就是那个正确答案。AUC_compare_correlated.m解决的问题更实际两个模型在同一个测试集上跑出来的AUC分别是0.78和0.82这个0.04的差异是真实进步还是随机波动% AUC_compare_correlated.m —— 两个相关模型AUC差异的DeLong检验 % 输入: y为真实标签score1、score2为两个模型对同一批样本的预测概率 % 输出: z为检验统计量p为显著性P值 function [z, p] AUC_compare_correlated(y, score1, score2) [~, ~, ~, auc1] perfcurve(y, score1, 1); [~, ~, ~, auc2] perfcurve(y, score2, 1); % DeLong法核心计算两个AUC的方差-协方差矩阵 V delongVariance(y, score1, score2); % 2x2协方差矩阵 d auc1 - auc2; se sqrt(V(1,1) V(2,2) - 2 * V(1,2)); z d / se; p 2 * (1 - normcdf(abs(z))); fprintf(AUC1%.3f, AUC2%.3f, Delta%.3f, P%.4f\n, auc1, auc2, d, p); end逻辑说明两个模型在相同样本上评估时它们的AUC并不独立——同一个样本被模型A分错时很可能也被模型B分错。DeLong法通过计算排序秩的协方差来刻画这种相关性再检验差值是否显著区别于0这就是它比朴素z检验稳健的原因。delongVariance对应的是CIAUC.m里按DeLong法构造协方差矩阵的实现我在这里用函数名做示意。参数说明score1和score2必须来自同一批样本且顺序与y完全一致。我见过有人把训练集的score1和验证集的score2塞进这个函数算出来的p值极其显著实际上是不同样本在比结论完全站不住脚。调用之前先确认三个变量length一致再确认score1和score2在同一行上对应的是同一个个体。4.3 roc.xlsx里到底要导出什么AUC、最优阈值、灵敏度与特异度压缩包里的roc.xlsx是Roc.m运行后留下的中间结果通常会包含每个样本的预测概率、不同阈值下的灵敏度和特异度。写报告或论文时你需要从中提取三组核心信息AUC及95%置信区间Youden指数对应的最优阈值以及该阈值下的灵敏度与特异度如果要画校准图还需要每个十分位组的预测概率均值与实际事件率。最优阈值一般按Youden指数确定即灵敏度特异度-1最大化的点。这个值不一定是0.5在事件率只有15-20%的数据里最优阈值经常落在0.2-0.3附近。提取这些指标时别手动读图写一个循环把结果自动追加到roc.xlsx对应的sheet里后续做多模型对比会省很多事。另外这套资源里还能计算PR曲线——当事件率很低时PR曲线比ROC更能暴露模型在少数类上的表现差异报告中通常主放ROC、补充PR这个组合在审稿人那里更耐看。5. 避坑与常见问题跑这套预测模型最容易翻车的五个地方跑风险预测模型最让人头疼的不是统计公式而是一些看起来不起眼的小问题。以下是压缩包在真实环境下使用时的五个高频翻车现场按现象→原因→解决的格式整理。5.1 现象模型AUC高达0.95验证集上一塌糊涂原因大概率是数据泄露。最常见的情况是X矩阵里混入了y那一列本身或者之前有人用全量数据做过变量筛选再拆分训练/测试使得验证集的信息已经参与了特征选择。另外行删除后没有重新对齐标签也会造成类似的假高。解决检查LR_Pz_choice.m传入的X是否包含结局变量列确认所有变量筛选只在origin.xlsx上完成把test.xlsx单独跑一遍Roc.m对比训练AUC与验证AUC的差距。如果验证AUC比训练低出0.1以上优先怀疑数据泄露而不是模型过拟合。快速排查可以用下面这段% 快速排查训练矩阵中是否混入结局列 leak any(all(X y, 1)); if leak warning(发现疑似结局变量泄漏请检查变量名); end上面这段用all(X y, 1)逐列比较X和y如果某列对每一行都与y相同这列极可能就是结局变量本身错放进了预测因子矩阵。这个检查只适用于数值型矩阵y是逻辑型时先double(y)再比较。5.2 现象NRI算出来是负数但新模型的AUC明明略高原因NRI衡量的是重分类是否正确与AUC衡量的整体排序能力是两回事。新模型确实可能把一批人从低风险正确升到高风险同时也把另一批人错误降级如果后者净效果更明显NRI就会为负而AUC反而微升。解决先确认NRI.m里的y定义与Roc.m一致谁是事件组再看cutoff阈值是否与临床常识匹配。如果事件率约15%却把cutoff设在0.5绝大多数样本都被分到低风险组NRI几乎失去区分度。这时候用Category_Free_NRI.m做连续NRI它会逐一对每个样本的预测概率变化计分不依赖阈值结果更稳。5.3 现象Roc.m画出来的曲线是锯齿状不是平滑阶梯原因score不是连续的预测概率。要么是逻辑回归输出的概率被四舍五入成了整数评分要么是score没有按降序排列就直接传给perfcurve。解决先检查unique(score)的长度如果不到10个值说明概率被离散化需要检查数据读入时的精度转换再把score和y按score降序同步排序后传给Roc.m% 画ROC前先完成同步排序 [scoreSorted, idx] sort(score, descend); ySorted y(idx);默认情况下perfcurve会自己处理排序但如果你需要手动画阶梯图或者叠加多条曲线就必须让score和y按同一个idx重排否则每条曲线对应的样本顺序都不一样图形就乱了。idx是sort返回的索引y必须用这个idx重排这是最关键的一步。5.4 现象xlsread读进来后整列是NaN但Excel里肉眼看着全是数字原因这一列的单元格被Excel设成了文本格式部分单元格左上角还有绿色角标。MATLAB读入时默认把这些当字符处理转数值时失败就成了NaN。解决在Excel里选中整列用数据-分列-完成强制转成常规格式更省事的办法是让读入脚本用TreatAsMissing选项统一处理缺失标记后再读。清洗完立即检查每列缺失比例% 读入后立即检查每一列的缺失比例 nanRatio mean(isnan(X), 1); disp([varNames; num2cell(nanRatio)]);一行命令就能输出每个变量的NaN占比如果看到某列是1.0说明这一列根本没有读入数值。X是数值矩阵varNames是表头变量名元胞数组两个对齐后打印出来比逐个查看Excel省时得多。5.5 现象AUC_compare_correlated.m报矩阵维度不一致但三个变量看起来都正常原因y、score1、score2的行数不一致。最常见的源头是读取Excel时用了不同的range分别读三列其中一列的行数因为表尾空行被截掉了。解决别分开读直接用readtable把整张表读进来再用列名取数这样三个向量天然对齐。如果必须从不同文件取数据读完立刻用assert(isequal(length(y), length(score1), length(score2)))做一次性校验比后面报错时再来排查快得多。补充一条环境问题在Windows解压这份rar后__MACOSX和以._开头的文件不会自己消失MATLAB扫描当前目录时会把它们识别为无扩展名文件偶尔干扰addpath。解压后先手动清理掉再打开main.m能避免很多奇怪的路径报错。6. NRI快速验证AUC没涨时怎么证明新模型仍然有价值6.1 场景先行AUC只涨0.01但不代表模型没有进步一个常见场景基线模型AUC0.78加入新标志物后AUC0.79DeLong检验P0.20看起来新变量没用。但如果你逐个样本去看重分类会发现一批原本被分到中风险的患者被正确提升到了高风险组这正是AUC这种整体排序指标捕捉不到的局部改善。NRI就是为这种场景设计的它不做全局排序比较而是直接统计新旧模型对每个样本的风险分类变化净有多少人从分错变成分对。6.2 标准NRI与分类自由NRI的选型和使用% NRI.m —— 标准NRI需要指定风险分层阈值 % 输入: y为二分类结局pOld/pNew为旧/新模型对同一批样本的预测概率 % cutoff为风险分层切点向量例如[0.2, 0.8] [nriStandard, comp] NRI(y, pOld, pNew, [0.2, 0.8]); % Category_Free_NRI.m —— 分类自由NRI不需要阈值适合连续评分比较 [nriFree, compFree] Category_Free_NRI(y, pOld, pNew);逻辑说明标准NRI把预测概率切分成低、中、高三段然后分别统计事件组和非事件组中正确重分类与错误重分类的人数差。分类自由NRI则不做切段它对每个样本按新旧概率差直接打分事件组的概率上升记为1、下降记为-1非事件组相反最后求和。两者计算方式不同但回答的是同一类问题新模型到底有没有把人群分得更准。参数说明pOld和pNew必须来自同一批样本且顺序与y严格对齐。cutoff的选择要基于临床意义不要用数据驱动的方式寻优比如低风险20%、高风险80%这种定义。还有一个小前置动作先跑FindNandD.m确认事件组和非事件组人数N和D如果事件组样本太少标准NRI的分组会很稀疏直接改用连续NRI更稳妥。从那以后我每次完成一个风险预测模型都会强制走一遍NRI检查先看AUC和置信区间再算NRI两者互相印证。如果AUC变化不显著但NRI显著为正报告里就重点讲重分类价值如果NRI为负而AUC上升我会回头看阈值设置和样本分布。希望这个习惯也能帮你减少一些被审稿人追问的尴尬。本文还有配套的精品资源点击获取