2024-04-24
nonp.auc <- function(x,y) {
Ker=outer(x,y,function(x,y) (sign(y-x)+1)/2)
auc=mean(Ker)
return(auc)
}
nonp.auc.check <- function(data.1,data.2) {
auc.vec=numeric(ncol(data.1))
for (i in 1:ncol(data.1)) {
new.1=data.1[,i]
new.2=data.2[,i]
auc.vec[i]=nonp.auc(new.1,new.2) }
return(auc.vec)
}
nonpar.combine2.auc <- function(alpha,rate,data.1,data.2) {
new.1 = data.1%*%c(alpha,rate)
new.2 = data.2%*%c(alpha,rate)
nonp.auc(new.1,new.2)
}
nonpar.combine2.coef <- function(dat_1,dat_2,evalnum=201) {
rate=seq(-1,1,length=evalnum)
alpha=rev(rate)[-1]
auc.rate_x = sapply(rate, nonpar.combine2.auc, alpha=1,data.1=dat_1,data.2=dat_2)
auc.alpha_x = sapply(alpha, nonpar.combine2.auc, rate=1,data.1=dat_1,data.2=dat_2)
auc.0 = c(auc.rate_x,auc.alpha_x)
vmax.idx = which.max(auc.0)
if(vmax.idx<=evalnum) return(c(alpha=1,rate=rate[vmax.idx], auc.max=auc.0[vmax.idx]))
if(vmax.idx> evalnum) return(c(alpha=alpha[vmax.idx-evalnum],rate=1,auc.max=auc.0[vmax.idx]))
}
step.coef <- function(new.1,new.2,design='step-down') {
n1 = nrow(new.1)
n2 = nrow(new.2)
VARnum = ncol(new.1)
combcoef = matrix(0,nrow=VARnum-1,ncol=2)
if (design=='step-down') {
auc.order = sort(nonp.auc.check(data.1=new.1,data.2=new.2),index.return=T,decreasing=T)$ix } else {
auc.order = sort(nonp.auc.check(data.1=new.1,data.2=new.2),index.return=T,decreasing=F)$ix }
combmarker.1=new.1[,auc.order[1]]
combmarker.2=new.2[,auc.order[1]]
final.coef=1
for (i in 2:VARnum) {
combmarker.1 = cbind(combmarker.1,new.1[,auc.order[i]])
combmarker.2 = cbind(combmarker.2,new.2[,auc.order[i]])
temp.info = nonpar.combine2.coef(combmarker.1,combmarker.2)
combcoef[i-1,] = temp.info[1:2]
combmarker.1 = combmarker.1%*%temp.info[1:2]
combmarker.2 = combmarker.2%*%temp.info[1:2]
final.coef = c(final.coef*combcoef[i-1,1],combcoef[i-1,2])
}
final.coef = final.coef[sort(auc.order,index.return=T)$ix]
return(list(coef=as.numeric(final.coef),
auc.combined=as.numeric(temp.info[3])))
}
library(mvtnorm)
df_neg = rmvnorm(100,mean=rep(0,6),sigma=diag(6))
df_pos = rmvnorm(100,mean=c(4,3,2,8,6,1)/3, sigma=diag(6)+0.2)
nonp.auc.check(df_neg,df_pos)
sort(nonp.auc.check(df_neg,df_pos),index.return=T,decreasing=T)
step.coef(df_neg,df_pos)
nonp.auc(df_neg%*%step.coef(df_neg,df_pos)$coef,df_pos%*%step.coef(df_neg,df_pos)$coef)
What do we do when outcome is right censored survival time?