Hi everyone,

I am attempting to copy code taken from the article "Decision curve analysis: a technical note" (https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6123195/)

When I try to make the simple model taken from the article, I get this error:

"Error in ntbft(data = dt, outcome = outcome, from = model1, xstart = xstart, :

unused argument (from = model1)"

I have no explanation for why this is happening, I am 99% sure I am copying the code correctly. I would appreciate the help!

This code requires nothing to download from except one package, so if anyone would like to copy the code themselves to see if they get this error, here it is:

set.seed(123)

n<-500

rr<-round(abs(rnorm(n,30,10)))

hr<-round(abs(rnorm(n,90,20)))

crp<-round(abs(rnorm(n,150,80)))

install.packages("dummies")

beta0=-7;betarr=0.05

betahr=0.02;betacrp=0.02

linpred<-cbind(1,rr,hr,crp)%*%

c(beta0,betarr,betahr,betacrp)

pi<-exp(linpred)/(1+exp(linpred))

sepsis.tag<-rbinom(n=n,size=1,prob=pi)

dt<-data.frame(rr,hr,crp,sepsis.tag)

ntbft<-function(data,outcome, frm = NULL,

exterdt=NULL, pred=NULL, xstart=0.01,

xstop=0.99, step=0.01, type="treated"){

pt<-seq(from=xstart,to=xstop,by=step)

lpt<-length(pt)

if(type=="treated") coef<-cbind(rep(1,lpt),rep(0,lpt))

if(type=="untreated") coef<-cbind(rep(0,lpt),rep(1,lpt))

if(type=="overall") coef<-cbind(rep(1,lpt),rep(1,lpt))

if(type=="adapt") coef<-cbind(1-pt,pt)

reponse<-as.vector(t(data[outcome]))

if(is.data.frame(exterdt))response<-as.vector(t(exterdt[outcome]))

event.rate<-mean(response)

nball<-event.rate-(1-event.rate)*pt/(1-pt)

nbnone<-1-event.rate-event.rate*(1-pt)/pt

if(is.null(pred)){

model<-glm(frm,data=data,family=binomial("logit"))

pred<-model$fitted.values

if(is.data.frame(exterdt))

pred<-predict(model,newdata=exterdt,type="response")

}

N<-length(pred)

nbt<-rep(NA,lpt)

nbu<-rep(NA,lpt)

for(t in 1:lpt){

tp<-sum(pred>=pt[t]&response==1)

fp<-sum(pred>=pt[t]&response==0)

fn<-sum(pred<pt[t]&response==1)

tn<-sum(pred<pt[t]&response==0)

nbt[t]<-tp/N-fp/N*(pt[t]/(1-pt[t]))

nbu[t]<-tn/N-fn/N*((1-pt[t])/pt[t])

}

nb<-data.frame(pt)

names(nb)<-"threshold"

nb["all"]<-coef[,1]*nball

nb["none"]<-coef[,2]*nbnone

nb["pred"]<-coef[,1]*nbt+coef[,2]*nbu

return(nb)

}

outcome<-"sepsis.tag"

model1<-sepsis.tag~rr+hr

model2<-sepsis.tag~rr+hr+crp

xstart<-0.01;xstop<-0.99;step<-0.01

type<-"treated"

nb.simple<-ntbft(data=dt, outcome=outcome,

from=model1,xstart=xstart,xstop=xstop,

step=step, type=type)