[Diagnosismed-commits] r18 - in pkg/DiagnosisMed: R man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Mar 1 23:02:24 CET 2010
Author: pedrobrasil
Date: 2010-03-01 23:02:24 +0100 (Mon, 01 Mar 2010)
New Revision: 18
Added:
pkg/DiagnosisMed/R/summary.diag.R
Removed:
pkg/DiagnosisMed/R/diagnosisI.r
pkg/DiagnosisMed/man/plot.diag.Rd
Modified:
pkg/DiagnosisMed/R/LRgraph.r
pkg/DiagnosisMed/R/ROC.r
pkg/DiagnosisMed/R/diagnosis.r
pkg/DiagnosisMed/R/plot.ROC.r
pkg/DiagnosisMed/R/plot.diag.r
pkg/DiagnosisMed/R/print.ROC.r
pkg/DiagnosisMed/R/print.diag.r
pkg/DiagnosisMed/R/zzz.r
pkg/DiagnosisMed/man/LRgrgaph.Rd
pkg/DiagnosisMed/man/ROC.Rd
pkg/DiagnosisMed/man/TGROC.Rd
pkg/DiagnosisMed/man/diagnosis.Rd
Log:
Modified: pkg/DiagnosisMed/R/LRgraph.r
===================================================================
--- pkg/DiagnosisMed/R/LRgraph.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/LRgraph.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -1,16 +1,17 @@
-LRgraph<-function(a=cbind(t1,t2),lwd=2,lty=1,cex=1,leg.cex=1.5,pt.cex=2,...){
- plot(1-a[[6,1]],a[[4,1]],xlim=c(0,1),ylim=c(0,1),xlab="False positive rate",ylab="True positive rate",col=1,cex=cex,lwd=lwd,lty=lty)
- abline(coef=c(0,((a[[4,1]])/(1-a[[6,1]]))),lwd=lwd)
- abline(coef=c(1-1*((1-a[[4,1]])/(1-(1-a[[6,1]]))),(1-a[[4,1]])/(1-(1-a[[6,1]]))),lwd=lwd)
- abline(v=1-a[[6,1]],lty=6,col="lightgray",lwd=lwd)
- abline(h=a[[4,1]],lty=6,col="lightgray",lwd=lwd)
- fill.col<-c(1)
- symbol<-c(1)
- for(i in 2:ncol(a))
- {
- points(1-a[[4,i]],a[[6,i]],col=i,pch=i,cex=cex,lwd=lwd,lty=lty)
- fill.col<-c(fill.col,i)
- symbol<-c(symbol,i)
- }
- legend("bottomright",legend=colnames(a),col=fill.col, pch=symbol, bty="n",cex=leg.cex,pt.cex=pt.cex,pt.lwd=lwd)
-}
+LRgraph <- function (tests, lwd = 2, lty = 1, cex = 1, leg.cex = 1.5, pt.cex = 2, ...){
+ plot(1 - tests[[6, 1]], tests[[4, 1]], xlim = c(0, 1), ylim = c(0,1), xlab = "False positive rate", ylab = "True positive rate",
+ col = 1, cex = cex, lwd = lwd, lty = lty)
+ abline(coef = c(0, ((tests[[4, 1]])/(1 - tests[[6, 1]]))), lwd = lwd)
+ abline(coef = c(1 - 1 * ((1 - tests[[4, 1]])/(1 - (1 - tests[[6,1]]))), (1 - tests[[4, 1]])/(1 - (1 - tests[[6, 1]]))), lwd = lwd)
+ abline(v = 1 - tests[[6, 1]], lty = 6, col = "lightgray", lwd = lwd)
+ abline(h = tests[[4, 1]], lty = 6, col = "lightgray", lwd = lwd)
+ fill.col <- c(1)
+ symbol <- c(1)
+ for (i in 2:ncol(tests)) {
+ points(1 - tests[[6, i]], tests[[4, i]], col = i, pch = i, cex = cex, lwd = lwd, lty = lty)
+ fill.col <- c(fill.col, i)
+ symbol <- c(symbol, i)
+ }
+ legend("bottomright", legend = colnames(tests), col = fill.col, pch = symbol, bty = "n", cex = leg.cex, pt.cex = pt.cex,
+ pt.lwd = lwd)
+}
Modified: pkg/DiagnosisMed/R/ROC.r
===================================================================
--- pkg/DiagnosisMed/R/ROC.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/ROC.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -5,14 +5,14 @@
Prevalence=0,
Plot=TRUE,
Plot.point="Min.ROC.Dist",
- cex.sub=.85,
- cex=1,
- lwd=1,
p.cex=1,
- Print.full=FALSE,
+ Full=FALSE,
Print=TRUE
){
# A simple warning ...
+ if(any(is.na(test) | is.na(gold))){
+ stop('It seems there are NAs either in the index test or in the reference test. Consider imputing or removing NAs!')
+ }
test.table<-table(test,gold)
if (dim(test.table)[2] != 2){
stop("It seems that your gold standard has more than 2 categories")
@@ -147,12 +147,12 @@
class(reteval)<-"ROC"
if(Print==TRUE){
- if(Print.full==TRUE){ print(reteval,Full=TRUE) }
+ if(Full==TRUE){ print(reteval,Full=TRUE) }
else{ print(reteval) }
}
# the plot commands
if(Plot==TRUE){
- plot(reteval,Plot.point=Plot.point,cex.sub=cex.sub,cex=cex,lwd=lwd,p.cex=p.cex)
+ plot(reteval,Plot.point=Plot.point,p.cex=p.cex)
}
invisible(reteval)
}
\ No newline at end of file
Modified: pkg/DiagnosisMed/R/diagnosis.r
===================================================================
--- pkg/DiagnosisMed/R/diagnosis.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/diagnosis.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -1,23 +1,70 @@
-diagnosis <- function(gold,test,CL=0.95,print=TRUE,plot=FALSE){
+diagnosis <- function(a,b=NULL,c=NULL,d=NULL,CL=0.95,print=TRUE,plot=FALSE){
#require(epitools)
- #require(epicalc)
- # to do ...
- # by option (multicenter)
- # test with 3 categories (indeterminate results)
- # testef <- as.factor(teste)
- # if(nlevels(teste) 2)
- # {
- # analysis this way
- # }
- tab<-table(test,gold,dnn=c(deparse(substitute(test)),deparse(substitute(gold))))
- dimnames(tab) <- list(test = c("negative" , "positive"), gold.standard = c("negative","positive"))
- tabmarg<-addmargins(tab)
+ if(is.numeric(a)){
+ if(all(length(a)==1 & length(b)==1 & length(c)==1 & length(d)==1 & !is.matrix(a))){
+ reference.name <- 'Not informed'
+ index.name <- 'Not informed'
+ tab<-as.table(cbind(rbind(d,c),rbind(b,a)))
+ dimnames(tab)<-list(index.test=c("negative","positive"),reference.standard=c("negative","positive"))
+ TN<-d
+ FN<-b
+ FP<-c
+ TP<-a
+ }
+ if(all(length(a) > 1 & length(b) > 1 & is.null(c) & is.null(d) & !is.matrix(a))){
+ if(any(is.na(a),is.na(b))){stop('There are NAs either in index test or reference standard. Consider removing or inputing!')}
+ if(nlevels(as.factor(a))!=2 | nlevels(as.factor(b))!=2){
+ stop('It seems there are more levels then 0 and 1.')
+ }
+ if(!all(levels(as.factor(a))==c(0,1) & levels(as.factor(b))==c(0,1))){
+ stop('Either the index test or the reference test is not correctly coded. 0 and 1 were expected!')
+ }
+ else{reference.name <- deparse(substitute(a))
+ index.name <- deparse(substitute(b))
+ tab<-table(b,a,dnn=c(deparse(substitute(b)),deparse(substitute(a))))
+ TN<-tab[1,1]
+ FN<-tab[1,2]
+ FP<-tab[2,1]
+ TP<-tab[2,2]
+ }
+ }
+ if(any(is.table(a) | is.matrix(a))){
+ if(!all(dim(a)==c(2,2))){
+ stop('It seems the inputed table is not 2x2. Check your table output.')
+ }
+ else{tab<-a
+ reference.name <- names(dimnames(tab)[2])
+ index.name <- names(dimnames(tab)[1])
+ TN<-tab[1,1]
+ FN<-tab[1,2]
+ FP<-tab[2,1]
+ TP<-tab[2,2]
+ }
+ }
+ }
+ if(all(any(is.factor(a) | is.character(a)) & any(is.factor(b) | is.character(b)) & !is.matrix(a))){
+ if(any(is.na(a) | is.na(b))){
+ stop('There seem to be NAs either in the reference standard or index test. Consider removing or inputing!')
+ }
+ if(nlevels(as.factor(a))!=2 | nlevels(as.factor(b))!=2){
+ stop('It seems there are more levels then negative/absence and positive/presence.')
+ }
+ if(!all(levels(as.factor(a))==c("negative","positive") & levels(as.factor(b))==c("negative","positive")) &
+ !all(levels(as.factor(a))==c("absence","presence") & levels(as.factor(b))==c("absence","presence"))){
+ stop('It seems categories are not correctly coded in either the reference or index test.')
+ }
+ else{reference.name <- deparse(substitute(a))
+ index.name <- deparse(substitute(b))
+ tab<-table(b,a,dnn=c(deparse(substitute(b)),deparse(substitute(a))))
+ TN<-tab[1,1]
+ FN<-tab[1,2]
+ FP<-tab[2,1]
+ TP<-tab[2,2]
+ }
+ }
+
+ tabmarg<-addmargins(tab)
Conf.limit<-CL
- #dnn is option of table command - specifies the names of row and column
- TN<-tab[1,1]
- FN<-tab[1,2]
- FP<-tab[2,1]
- TP<-tab[2,2]
# sample size
n<-sum(tab)
# prevalence
@@ -73,11 +120,9 @@
# {ROC<-roc.from.table(tab, graph = FALSE)}
# gives same results as AUC<-(Se+Sp)/2
Youden<-Se+Sp-1
- Youden.inf.cl<-Youden-qnorm(CL/2)*sqrt(((Se * (1 - Se))/(TP+FN) +
- ((Sp * (1 - Sp))/(TN+FP))))
- Youden.sup.cl<-Youden+qnorm(CL/2)*sqrt(((Se * (1 - Se))/(TP+FN) +
- ((Sp * (1 - Sp))/(FP+TN))))
- #rm(ROC)
+ Youden.inf.cl<-Youden-qnorm(CL/2)*sqrt(((Se * (1 - Se))/(TP+FN) + ((Sp * (1 - Sp))/(TN+FP))))
+ Youden.sup.cl<-Youden+qnorm(CL/2)*sqrt(((Se * (1 - Se))/(TP+FN) + ((Sp * (1 - Sp))/(FP+TN))))
+# rm(ROC)
rm(tab)
# results evaluations
reteval <- list(tabmarg=tabmarg,n=n,p=p,Se=Se,Se.cl=Se.cl,Sp=Sp,Sp.cl=Sp.cl,PLR=PLR,
@@ -85,7 +130,7 @@
NLR.sup.cl=NLR.sup.cl,accu=accu,accu.cl=accu.cl,PPV=PPV,PPV.cl=PPV.cl,NPV=NPV,NPV.cl=NPV.cl,
DOR=DOR,DOR.inf.cl=DOR.inf.cl,DOR.sup.cl=DOR.sup.cl,ET=ET,ER=ER,ER.cl=ER.cl,
Youden=Youden,Youden.inf.cl=Youden.inf.cl,Youden.sup.cl=Youden.sup.cl,AUC=AUC,
- Conf.limit=Conf.limit)
+ Conf.limit=Conf.limit,reference.name=reference.name,index.name=index.name)
class(reteval) <- "diag"
if(print==TRUE) {print(reteval)}
invisible(reteval)
Deleted: pkg/DiagnosisMed/R/diagnosisI.r
===================================================================
--- pkg/DiagnosisMed/R/diagnosisI.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/diagnosisI.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -1,83 +0,0 @@
-diagnosisI <- function(TP,FN,FP,TN,CL=0.95,print=TRUE,plot=FALSE){
- #require(epitools)
- #equire(epicalc)
- tab<-as.table(cbind(rbind(TN,FP),rbind(FN,TP)))
- dimnames(tab)<-list(test=c("negative","positive"),gold.standard=c("negative","positive"))
- tabmarg<-addmargins(tab)
- Conf.limit<-CL
- #TN<-tab[1,1]
- #FN<-tab[1,2]
- #FP<-tab[2,1]
- #TP<-tab[2,2]
- # sample size
- n<-sum(tab)
- # prevalence
- p<-(TP+FN)/n
- # sensitivity and confidence limits
- Se<-TP/(TP+FN)
- Se.cl<-as.numeric(binom.wilson(TP, TP+FN, conf.level = CL)[4:5])
- # especificity and confidence limits
- Sp<-TN/(FP+TN)
- Sp.cl<-as.numeric(binom.wilson(TN, FP+TN, conf.level = CL)[4:5])
- # positive and negative likelyhood ratios and confidence limits
- PLR<-Se/(1-Sp)
- # LR confidence limists inspired in epi.tests{epiR}
- PLR.inf.cl<-exp(log(PLR)-(qnorm(1-((1-CL)/2),mean=0,sd=1))*sqrt((1-Se)/(
- (TP+FN)*Sp)+(Sp)/((FP+TN)*(1-Sp))))
- PLR.sup.cl<-exp(log(PLR)+(qnorm(1-((1-CL)/2),mean=0,sd=1))*sqrt((1-Se)/(
- (TP+FN)*Sp)+(Sp)/((FP+TN)*(1-Sp))))
- NLR<-(1-Se)/Sp
- NLR.inf.cl<-exp(log(NLR)-(qnorm(1-((1-CL)/2),mean=0,sd=1))*sqrt((Se)/((TP+
- FN)*(1-Se))+(1-Sp)/((FP+TN)*(Sp))))
- NLR.sup.cl<-exp(log(NLR)+(qnorm(1-((1-CL)/2),mean=0,sd=1))*sqrt((Se)/((TP+
- FN)*(1-Se))+(1-Sp)/((FP+TN)*(Sp))))
- #accuracy and confidence limits
- accu<-(TP+TN)/n
- accu.cl<-as.numeric(binom.wilson(TP+TN, n, conf.level = CL)[4:5])
- # positive and negative predictive values and confidence limits
- PPV<-TP/(TP+FP)
- PPV.cl<-as.numeric(binom.wilson(TP, TP+FP, conf.level = CL)[4:5])
- NPV<-TN/(TN+FN)
- NPV.cl<-as.numeric(binom.wilson(TN, TN+FN, conf.level = CL)[4:5])
- # diagnostic odds ratio and confidence limits
- OR<-oddsratio(tab,conf.level = CL)
- DOR<-OR$measure[2,1]
- #DOR<-(TP*TN)/(FP*FN)
- DOR.inf.cl<-OR$measure[2,2]
- DOR.sup.cl<-OR$measure[2,3]
- rm(OR)
- # error rate and error trade
- #ER<-((FN/(FN+TN))*p)+(((FP/(FP+TP))*(TN+FP))
- ER<-(FN+FP)/n
- ER.cl<-as.numeric(binom.wilson(FN+FP, n, conf.level = CL)[4:5])
- ET<-(FN/FP)
- # pre-test and pos-test odds (to do)
- # area under ROC curve
- AUC<-(Se+Sp)/2
- if(plot==TRUE){
- plot(1-Sp,Se,xlim=c(0,1),ylim=c(0,1))
- segments(0,0,1-Sp,Se,col="red")
- segments(1-Sp,Se,1,1,col="red")
- grid()
- }
- #if(plot==FALSE)
- # {ROC<-roc.from.table(tab, graph = FALSE)}
- # gives same results as AUC<-(Se+Sp)/2
- Youden<-Se+Sp-1
- Youden.inf.cl<-Youden-qnorm(CL/2)*sqrt(((Se * (1 - Se))/(TP+FN) +
- ((Sp * (1 - Sp))/(TN+FP))))
- Youden.sup.cl<-Youden+qnorm(CL/2)*sqrt(((Se * (1 - Se))/(TP+FN) +
- ((Sp * (1 - Sp))/(FP+TN))))
-# rm(ROC)
- rm(tab)
- # results evaluations
- reteval <- list(tabmarg=tabmarg,n=n,p=p,Se=Se,Se.cl=Se.cl,Sp=Sp,Sp.cl=Sp.cl,PLR=PLR,
- PLR.inf.cl=PLR.inf.cl,PLR.sup.cl=PLR.sup.cl,NLR=NLR,NLR.inf.cl=NLR.inf.cl,
- NLR.sup.cl=NLR.sup.cl,accu=accu,accu.cl=accu.cl,PPV=PPV,PPV.cl=PPV.cl,NPV=NPV,NPV.cl=NPV.cl,
- DOR=DOR,DOR.inf.cl=DOR.inf.cl,DOR.sup.cl=DOR.sup.cl,ET=ET,ER=ER,ER.cl=ER.cl,
- Youden=Youden,Youden.inf.cl=Youden.inf.cl,Youden.sup.cl=Youden.sup.cl,AUC=AUC,
- Conf.limit=Conf.limit)
- class(reteval) <- "diag"
- if(print==TRUE) {print(reteval)}
- invisible(reteval)
-}
\ No newline at end of file
Modified: pkg/DiagnosisMed/R/plot.ROC.r
===================================================================
--- pkg/DiagnosisMed/R/plot.ROC.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/plot.ROC.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -1,21 +1,21 @@
# the plot commands
-plot.ROC<-function(x,...,Plot.point="Min.ROC.Dist",cex.sub=.85,cex=1,lwd=1,p.cex=1){
+plot.ROC<-function(x,Plot.point="Min.ROC.Dist",cex.sub=.85,p.cex=1,...){
if(Plot.point!="None" & Plot.point!="Min.ROC.Dist" & Plot.point!="Max.Accuracy" &
Plot.point!="Max.DOR" & Plot.point!="Error.rate" & Plot.point!="Max.Accuracy.area" &
Plot.point!="Max.Sens+Spec" & Plot.point!="Max.Youden" & Plot.point!="Se=Sp" &
- Plot.point!="Min.ROC.Dist" & Plot.point!="Max.Efficiency" & Plot.point!="Min.MCT")
- stop("The Plot.point option is not correctly set! Type '?ROC' to check possible options.")
+ Plot.point!="Min.ROC.Dist" & Plot.point!="Max.Efficiency" & Plot.point!="Min.MCT"){
+ stop("The Plot.point option is not correctly set! Type '?ROC' to check possible options.")}
plot(1-x$test.diag.table$Specificity,x$test.diag.table$Sensitivity,type="l",
- col=1,lwd=lwd,xlab="1-Specificity",ylab="Sensitivity",xlim=c(0,1),ylim=c(0,1))
+ col=1,xlab="1-Specificity",ylab="Sensitivity",xlim=c(0,1),ylim=c(0,1))
grid()
segments(0,0,1,1,col="lightgray")
if(Plot.point=="None"){
legend("bottomright",legend=(
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- ),bty="n",cex=cex)}
+ ),bty="n")}
if(Plot.point=="Max.Accuracy")
{points(1-x$test.best.cutoff[1,5],x$test.best.cutoff[1,2],col=1,pch=19,cex=p.cex)
@@ -24,7 +24,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[1,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[1,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Max.DOR")
{points(1-x$test.best.cutoff[2,5],x$test.best.cutoff[2,2],col=1,pch=19,cex=p.cex)
@@ -34,7 +34,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[2,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[2,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Error.rate")
{points(1-x$test.best.cutoff[3,5],x$test.best.cutoff[3,2],col=1,pch=19,cex=p.cex)
@@ -44,7 +44,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[3,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[3,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Max.Accuracy.area")
{points(1-x$test.best.cutoff[4,5],x$test.best.cutoff[4,2],col=1,pch=19,cex=p.cex)
@@ -54,7 +54,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[4,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[4,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Max.Sens+Spec")
{points(1-x$test.best.cutoff[5,5],x$test.best.cutoff[4,2],col=1,pch=19,cex=p.cex)
@@ -64,7 +64,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[5,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[5,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Max.Youden")
{points(1-x$test.best.cutoff[6,5],x$test.best.cutoff[6,2],
@@ -75,7 +75,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[6,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[6,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Se=Sp")
{points(1-x$test.best.cutoff[7,5],x$test.best.cutoff[7,2],col=1,pch=19,cex=p.cex)
@@ -85,7 +85,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[7,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[7,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Min.ROC.Dist")
{points(1-x$test.best.cutoff[8,5],x$test.best.cutoff[8,2],col=1,pch=19,cex=p.cex)
@@ -95,7 +95,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[8,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[8,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Max.Efficiency")
{points(1-x$test.best.cutoff[9,5],x$test.best.cutoff[9,2],col=1,pch=19,cex=p.cex)
@@ -106,7 +106,7 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[9,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[9,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
if(Plot.point=="Min.MCT")
{points(1-x$test.best.cutoff[10,5],x$test.best.cutoff[10,2],col=1,pch=19,cex=p.cex)
@@ -118,5 +118,5 @@
paste("Sensitivity:",formatC(x$test.best.cutoff[10,2],digits=4,format="f")),
paste("Specificity:",formatC(x$test.best.cutoff[10,5],digits=4,format="f")),
paste("AUC:",formatC(x$AUC.summary[2],digits=4,format="f"))
- )),bty="n",cex=cex)}
+ )),bty="n")}
}
\ No newline at end of file
Modified: pkg/DiagnosisMed/R/plot.diag.r
===================================================================
--- pkg/DiagnosisMed/R/plot.diag.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/plot.diag.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -14,7 +14,7 @@
#error.rate<-numeric(100)
#ER<-((FN/(FN+TN))*p)+(((FP/(FP+TP))*(TN+FP))
#error.rate<-((1-x$NPV)*pre.test)+(((1-x$PPV)*(x$n-pre.test))
- Result<-cbind(pre.test,post.test)
+ Result<-as.data.frame(cbind(pre.test,post.test))
if(print==TRUE)
{print(Result)}
plot(pre.test,post.test,xlab="Pre-test probability (prevalence)"
@@ -23,5 +23,6 @@
# lines(error.rate,lty=2,col=2)
grid()
axis(3)
- legend("bottomright",legend=paste("Positive likelihood ratio: ",formatC(x$PLR, digits=4)))
+ legend("bottomright",legend=paste("Positive likelihood ratio: ",formatC(x$PLR, digits=4)),bty='n')
+ invisible(Result)
}
Modified: pkg/DiagnosisMed/R/print.ROC.r
===================================================================
--- pkg/DiagnosisMed/R/print.ROC.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/print.ROC.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -1,5 +1,5 @@
print.ROC<-function(x,Full=FALSE,...){
- if (Full==TRUE){ page(x$test.diag.table[,-c(2:5,24:34)],method="print")}
+ if (Full==TRUE){ View(x$test.diag.table[,-c(2:5,24:34)])}
cat(" Sample size:",x$sample.size,"\n")
cat(" Sample prevalence:",round(x$sample.prevalence,digits = 4),"\n")
cat("Population prevalence:",round(x$pop.prevalence,digits = 4)," - same as sample prevalence if not informed\n")
Modified: pkg/DiagnosisMed/R/print.diag.r
===================================================================
--- pkg/DiagnosisMed/R/print.diag.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/print.diag.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -1,22 +1,23 @@
print.diag <- function(x,...){
- cat("\n\n")
+ cat("Reference standard:",x$reference.name,"\n")
+ cat("Index test :",x$index.name,"\n")
+ cat("---------------------------------------------------------------\n")
print(x$tabmarg)
- cat("\n\n")
- cat("The test has the following parameters [",x$Conf.limit*100,"% confidence interval]\n")
+ cat("The test has the following parameters [",x$Conf.limit*100,"% confidence interval]\n",sep="")
cat("---------------------------------------------------------------\n")
- cat("Sample size: ",x$n,"\n")
- cat("Prevalence considered: ",formatC(x$p*100,digits=2,format="f"),"%\n")
- cat("Sensitivity: ",formatC(x$Se*100,digits=2,format="f"),"% [",formatC(x$Se.cl[1]*100,digits=2,format="f")," - ",formatC(x$Se.cl[2]*100,digits=2,format="f"),"]\n")
- cat("Specificity: ",formatC(x$Sp*100,digits=2,format="f"),"% [",formatC(x$Sp.cl[1]*100,digits=2,format="f")," - ",formatC(x$Sp.cl[2]*100,digits=2,format="f"),"]\n")
- cat("Positive predictive value: ",formatC(x$PPV*100,digits=2,format="f"),"% [",formatC(x$PPV.cl[1]*100,digits=2,format="f")," - ",formatC(x$PPV.cl[2]*100,digits=2,format="f"),"]\n")
- cat("Negative predictive value: ",formatC(x$NPV*100,digits=2,format="f"),"% [",formatC(x$NPV.cl[1]*100,digits=2,format="f")," - ",formatC(x$NPV.cl[2]*100,digits=2,format="f"),"]\n")
- cat("Positive likelihood ratio: ",formatC(x$PLR,digits=2,format="f")," [",formatC(x$PLR.inf.cl,digits=2,format="f")," - ",formatC(x$PLR.sup.cl,digits=2,format="f"),"]\n")
- cat("Negative likelihood ratio: ",formatC(x$NLR,digits=2,format="f")," [",formatC(x$NLR.inf.cl,digits=2,format="f")," - ",formatC(x$NLR.sup.cl,digits=2,format="f"),"]\n")
- cat("Diagnostic odds ratio: ",formatC(x$DOR,digits=2,format="f")," [",formatC(x$DOR.inf.cl,digits=2,format="f")," - ",formatC(x$DOR.sup.cl,digits=2,format="f"),"]\n")
- cat("Error trade off (FN : FP) ",round(x$ET,digits=2)," : 1 \n")
- cat("Error rate: ",formatC(x$ER*100,digits=2,format="f"),"% [",formatC(x$ER.cl[1]*100,digits=2,format="f")," - ",formatC(x$ER.cl[2]*100,digits=2,format="f"),"]\n")
- cat("Accuracy: ",paste(round(x$accu*100,digits=2)),"% [",paste(round(x$accu.cl[1]*100,digits=2))," - ",paste(round(x$accu.cl[2]*100,digits=2)),"]\n")
- cat("Youden index: ",formatC(x$Youden,digits=4,format="f")," [",formatC(x$Youden.inf.cl,digits=4,format="f")," - ",formatC(x$Youden.sup.cl,digits=4,format="f"),"]\n")
- cat("Area under ROC curve: ",round(x$AUC,digits=4),"\n")
+ cat("Sample size: ",x$n,"\n")
+ cat("Prevalence considered(%): ",formatC(x$p*100,digits=2,format="f"),"\n")
+ cat("Sensitivity(%): ",formatC(x$Se*100,digits=2,format="f")," [",formatC(x$Se.cl[1]*100,digits=2,format="f")," - ",formatC(x$Se.cl[2]*100,digits=2,format="f"),"]\n")
+ cat("Specificity(%): ",formatC(x$Sp*100,digits=2,format="f")," [",formatC(x$Sp.cl[1]*100,digits=2,format="f")," - ",formatC(x$Sp.cl[2]*100,digits=2,format="f"),"]\n")
+ cat("Positive predictive value(%): ",formatC(x$PPV*100,digits=2,format="f")," [",formatC(x$PPV.cl[1]*100,digits=2,format="f")," - ",formatC(x$PPV.cl[2]*100,digits=2,format="f"),"]\n")
+ cat("Negative predictive value:(%):",formatC(x$NPV*100,digits=2,format="f")," [",formatC(x$NPV.cl[1]*100,digits=2,format="f")," - ",formatC(x$NPV.cl[2]*100,digits=2,format="f"),"]\n")
+ cat("Positive likelihood ratio: ",formatC(x$PLR,digits=2,format="f")," [",formatC(x$PLR.inf.cl,digits=2,format="f")," - ",formatC(x$PLR.sup.cl,digits=2,format="f"),"]\n")
+ cat("Negative likelihood ratio: ",formatC(x$NLR,digits=2,format="f")," [",formatC(x$NLR.inf.cl,digits=2,format="f")," - ",formatC(x$NLR.sup.cl,digits=2,format="f"),"]\n")
+ cat("Diagnostic odds ratio: ",formatC(x$DOR,digits=2,format="f")," [",formatC(x$DOR.inf.cl,digits=2,format="f")," - ",formatC(x$DOR.sup.cl,digits=2,format="f"),"]\n")
+ cat("Error trade off (FN : FP) ",round(x$ET,digits=2)," : 1 \n",sep='')
+ cat("Error rate(%): ",formatC(x$ER*100,digits=2,format="f")," [",formatC(x$ER.cl[1]*100,digits=2,format="f")," - ",formatC(x$ER.cl[2]*100,digits=2,format="f"),"]\n")
+ cat("Accuracy(%): ",formatC(x$accu*100,digits=2,format="f")," [",formatC(x$accu.cl[1]*100,digits=2,format="f")," - ",formatC(x$accu.cl[2]*100,digits=2,format="f"),"]\n")
+ cat("Youden index: ",formatC(x$Youden,digits=4,format="f")," [",formatC(x$Youden.inf.cl,digits=4,format="f")," - ",formatC(x$Youden.sup.cl,digits=4,format="f"),"]\n")
+ cat("Area under ROC curve: ",round(x$AUC,digits=4),"\n")
cat("---------------------------------------------------------------\n")
}
Added: pkg/DiagnosisMed/R/summary.diag.R
===================================================================
--- pkg/DiagnosisMed/R/summary.diag.R (rev 0)
+++ pkg/DiagnosisMed/R/summary.diag.R 2010-03-01 22:02:24 UTC (rev 18)
@@ -0,0 +1,37 @@
+summary.diag <- function(object,...){
+ diag.tab <- matrix(
+ c(object$n,NA,paste(formatC(object$p*100,digits=2,format="f")),NA,
+ formatC(object$Se*100,digits=2,format="f"),
+ paste('[',formatC(object$Se.cl[1]*100,digits=2,format="f"),'-',formatC(object$Se.cl[2]*100,digits=2,format="f"),']'),
+ formatC(object$Sp*100,digits=2,format="f"),
+ paste('[',formatC(object$Sp.cl[1]*100,digits=2,format="f"),'-',formatC(object$Sp.cl[2]*100,digits=2,format="f"),']'),
+ formatC(object$PPV*100,digits=2,format="f"),
+ paste('[',formatC(object$PPV.cl[1]*100,digits=2,format="f"),'-',formatC(object$PPV.cl[2]*100,digits=2,format="f"),']'),
+ formatC(object$NPV*100,digits=2,format="f"),
+ paste('[',formatC(object$NPV.cl[1]*100,digits=2,format="f"),'-',formatC(object$NPV.cl[2]*100,digits=2,format="f"),']'),
+ formatC(object$PLR,digits=2,format="f"),paste('[',formatC(object$PLR.inf.cl,digits=2,format="f"),'-',formatC(object$PLR.sup.cl,digits=2,format="f"),']'),
+ formatC(object$NLR,digits=2,format="f"),paste('[',formatC(object$NLR.inf.cl,digits=2,format="f"),'-',formatC(object$NLR.sup.cl,digits=2,format="f"),']'),
+ formatC(object$DOR,digits=2,format="f"),paste('[',formatC(object$DOR.inf.cl,digits=2,format="f"),'-',formatC(object$DOR.sup.cl,digits=2,format="f"),']'),
+ paste(round(object$ET,digits=2),':1',sep=''),NA,
+ formatC(object$ER*100,digits=2,format="f"),
+ paste('[',formatC(object$ER.cl[1]*100,digits=2,format="f"),'-',formatC(object$ER.cl[2]*100,digits=2,format="f"),']'),
+ formatC(object$accu*100,digits=2,format="f"),
+ paste('[',formatC(object$accu.cl[1]*100,digits=2,format="f"),'-',formatC(object$accu.cl[2]*100,digits=2,format="f"),']'),
+ formatC(object$Youden,digits=4,format="f"),
+ paste('[',formatC(object$Youden.inf.cl,digits=4,format="f"),'-',formatC(object$Youden.sup.cl,digits=4,format="f"),']'),
+ round(object$AUC,digits=4),NA
+ ),nrow = 14, ncol=2, byrow=TRUE,
+ dimnames = list(c('Sample size:','Prevalence(%):','Sensitivity(%):','Specificity(%):',
+ 'Postive predictive value(%):','Negative predictive value(%):',
+ 'Positive likelihood ratio:','Negative likelihood ratio:',
+ 'Diagnostic Odds Ratio:','Error trade off (FN : FP):','Error rate(%):',
+ 'Accuracy(%):','Youden index:','Area under ROC curve:'),
+ c('Estimate', paste(object$Conf.limit,'Confidence limits')))
+ )
+ diag.tab <- format(diag.tab,justify = "centre",na.encode = FALSE,trim = TRUE)
+ colnames(diag.tab) <-format(colnames(diag.tab),justify = "centre",trim = TRUE)
+ cat("--------------------------------------------------------------------------\n")
+ print(diag.tab,quote=FALSE,na.print='')
+ cat("--------------------------------------------------------------------------\n")
+ invisible(diag.tab)
+}
Modified: pkg/DiagnosisMed/R/zzz.r
===================================================================
--- pkg/DiagnosisMed/R/zzz.r 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/R/zzz.r 2010-03-01 22:02:24 UTC (rev 18)
@@ -1,7 +1,7 @@
# first and last
.First.lib <- function(lib, pkg)
{
-#require("epitools","epicalc","TeachingDemos","tcltk",quietly=TRUE,warn.conflicts=FALSE)
+#require("epitools","TeachingDemos","tcltk",quietly=TRUE,warn.conflicts=FALSE)
# if (.Platform$OS.type=="windows")
see <- packageDescription(pkg,fields="Version")
cat("'DiagnosisMed' library",see," loaded\n",sep=" ")
Modified: pkg/DiagnosisMed/man/LRgrgaph.Rd
===================================================================
--- pkg/DiagnosisMed/man/LRgrgaph.Rd 2009-08-18 18:28:51 UTC (rev 17)
+++ pkg/DiagnosisMed/man/LRgrgaph.Rd 2010-03-01 22:02:24 UTC (rev 18)
@@ -2,14 +2,13 @@
\alias{LRgraph}
\title{Comparing diagnositic tests: a simple graphic using likelihood ratios.}
\description{
-LRgraph graphicaly compares two or more (all of them with the first test) diagnostic tests with binary results through their likelihood ratios, based on the rationale that the predictive hability of a test is a more interresting charachteristic than sensitivity and/or specificity. It is possible to see thruogh the graph that if the tests with smaller sensitivity or specifity may have superior predictive habilty, that is, increases the prediction hability with small sensitivity/specificity trade-off.
+LRgraph graphically compares two or more (all of them with the first test) diagnostic tests with binary results through their likelihood ratios, based on the rationale that the predictive ability of a test is a more interesting characteristic than sensitivity and/or specificity. It is possible to see through the graph that if the tests with smaller sensitivity or specificity may have superior predictive ability, that is, increases the prediction ability with small sensitivity/specificity trade-off.
}
\usage{
-LRgraph(a = cbind(t1, t2), lwd = 2, lty = 1, cex = 1, leg.cex = 1.5, pt.cex = 2,...)
+LRgraph(tests, lwd = 2, lty = 1, cex = 1, leg.cex = 1.5, pt.cex = 2,...)
}
\arguments{
- \item{a}{a is a matrix composed by tow or more tests. Therefore, the user should edit only what is between the parethesis in a=cbind().The user may insert as many tests as one wishes. See t1, t2 below.}
- \item{t1,t2}{ t1, t2... tn are objects created by the \code{\link{diagnosis}} function (see example below). Therefore, one should first analyze data for each single test and store it in a object. Later, stick these objects in the LRgraph command. The objects names will appear in the graph legend.}
+ \item{tests}{a is a object composed by two or more tests. This object should be created binding two objects created by diagnosis functions as 'cbind(mytest1,mytest2)'.The user may insert as many tests as one wishes. See below.}
\item{lwd}{Line width. See \link[graphics]{par},\link[graphics]{points},\link[graphics]{legend}}
\item{lty}{Line type. See \link[graphics]{par}}
\item{cex}{Symbols and text size. See \link[graphics]{par},\link[graphics]{points}}
@@ -17,9 +16,7 @@
\item{pt.cex}{Size of the symbols in the legend. See \link[graphics]{legend}}
\item{...}{Other graph parameters. See \link[graphics]{plot.default}}
}
-\details{When a diagnostic test has both sensitivity and specificity higher than a competing test is easy to see that the former is superior than the later. However, sometimes a test may have superior sensitivity and inferior specificity (or the other way around). In this case, a good decision may be toward the test that have a better prediction ability. The graph visually helps the user to see and compare these abilities. The graph is very similar to the ROC graph. The vertical and horizontal axis have the same length as the ROC graph. However, the diagnostic tests are represented as dots instead of curves. The solid line passing through (0,0) is the likelihood ratio positive-line and the solid line passing through (1,1) is the likelihood ratio negative-line. Both negative and positive likelihood are numerically equivalent to the slopes of the solid lines. The solid lines split the graph into four areas (run the example). Also, there are dashed lines representing the sensitivity and specificity of the first test plotted. One may see that there are areas that a test may have superior sensitivity (or specificity) and yet the dot may be below the likelihood solid line. That is because the sensitivity / specificity trade-off is not reasonable, making the test with less predictive ability.
-
-}
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/diagnosismed -r 18
More information about the Diagnosismed-commits
mailing list