预测模型的局部评价?

预测模型的局部评价

• 为什么要进行局部评价？
  首先是临床决策曲线分析通常会给预测模型的使用规定一个阈值范围，相应地预测模型的评价也应该局限在这个范围之内才是合理的；
  其次，全局性地评价往往不够敏感，即好的模型和坏的模型之间分值差别不大。

局部评价的指标：

目前了解到的并可以实现的指标有三种，

• 局部ROC曲线下面积：又有三种计算方式，pROC的partial.auc、双向AUC和整合pAUC。前两种有专门的R包可以实现，后一种的代码如下（笔者自行撰写，欢迎大家考察）。
  参考文献：Carrington AM, Fieguth PW, Qazi H, Holzinger A, Chen HH, Mayr F, Manuel DG. A new concordant partial AUC and partial c statistic for imbalanced data in the evaluation of machine learning algorithms. BMC Med Inform Decis Mak. 2020 Jan 6;20(1):4. doi: 10.1186/s12911-019-1014-6. PMID: 31906931; PMCID: PMC6945414.

#计算局部AUC的函数，input：概率、真实值、阈值范围
pAUCc=function(pred_prob,true_val,thres_range){
    #pROC，获取threhold，specificity和sensitivity的数据，并取两位小数，便于与输入的数据匹配。
    roc_obj=roc(true_val,pred_prob)
    df=coords(roc_obj)
    df<- df%>%lapply(function(x) round(x,2))%>%as.data.frame()%>%unique()
    # print(df%>%filter(threshold>=0.1 & threshold<=0.11))#threhold小 spec小
    #将输入的threshold匹配成specificity(x)或者sensitivity(y)
    spec_ls_thres_1=df%>%filter(threshold==thres_range[1])%>%select(specificity)#取小的，序号1
    # print(spec_ls_thres_1)
    spec_ls_thres_2=df%>%filter(threshold==thres_range[2])%>%select(specificity)#取大的序号length()
    # print(spec_ls_thres_2)
    #由于precrec这个R包识别的X是1-specificity,所以要转成1-
    spec_end=1-min(spec_ls_thres_1)#
    spec_start=1-max(spec_ls_thres_2)#
    
    sen_end=df%>%filter(specificity==min(spec_ls_thres_1))%>%select(sensitivity)#sen取大的，取1
    sen_start=df%>%filter(specificity==max(spec_ls_thres_2))%>%select(sensitivity)#sen

    # print(spec_start)
    # print(spec_end)
    sen_start=min(sen_start)#pass
    sen_end=max(sen_end)#pass
    # print(sen_start)
    # print(sen_end)
    #计算pAUCc
    curve=precrec::evalmod(scores=pred_prob,labels=true_val)  #计算
    curve_part_x=precrec::part(curve,xlim=c(spec_start,spec_end),curvetype='ROC')#计算纵向的AUC，传统的，限制specficity
    curve_part_y=precrec::part(curve,ylim=c(sen_start,sen_end),curvetype='ROC')#计算横向的AUC，文献知道的，限制sensitivity
    pAUCy=precrec::pauc(curve_part_x)%>% select(paucs)%>%na.omit()
    plot(curve_part_x)
    pAUCx=precrec::pauc(curve_part_y)%>% select(paucs)%>%na.omit()
    plot(curve_part_y)
    pAUCc=(pAUCy[1]+pAUCx[1])/2 
    return (pAUCc[1])
}

• 局部校准曲线，这个的实现是val.prob函数中限定lim即可实现。

• 局部决策曲线下面积：代码如下（笔者根据文献自行撰写，欢迎大家反馈）
  参考文献：Talluri R, Shete S. Using the weighted area under the net benefit curve for decision curve analysis. BMC Med Inform Decis Mak. 2016 Jul 18;16:94. doi: 10.1186/s12911-016-0336-x. PMID: 27431531; PMCID: PMC4949771.

##第三个评价指标，净阈值下面积
library(tidymodels)
tidymodels_prefer()
#使用以下函数，需要先将dataframe转成tibble，以tibble名称，预测概率列（需要加引号），实际分类列（0，1，不能加引号）
NBC<-function(pt,df=result_df,pred_prob='predict_logi_prob',real=status){
	result_df=as_tibble(result_df)#data.frame 转tibble
	y_pred_label=df[[pred_prob]]>pt
	df$y_pred_label<-ifelse(y_pred_label=='FALSE',0,1)%>%factor(levels=c(0,1))
	cm<-conf_mat(df,{{real}},y_pred_label)%>%tidy()
	tp=cm$value[4]
	fp=cm$value[3]
	n = length(df[[pred_prob]])
	# prevalence=1011/(90025+1011)
#这里计算的是演示的NB值，如果想计算标准化的NB值，需要除以发病率（prevalence）
	net_benefit = (tp / n) - (fp / n) * (pt / (1 - pt))
	# net_benefit_std=net_benefit/prevalence
	return (net_benefit)
}
#WA-NBC,使用的积分函数，
stats::integrate(function(x){NBC(x,df=result_df)},interest_range[1],interest_range[2])

代码一键运行和更详细信息参考

https://www.heywhale.com/mw/project/643ceb7c446c45f4595416f1