CPvar=function(time,cens,group=rep(1,length(time)))
{
####################################################
# this function calculates                         #
#     the CIF for the event of interest            #
#     the CIF for competing risks (fcr)            #  
#     the conditional probability                  #
#     the variance for the conditional probability #
#           based on Pepe and Mori-s article in    #
#           Statistics in Medicine, 1993.          #
####################################################
lg=labels(table(group))$group
ng=length(lg)
dd=table(time,cens,group)
tt=sort(unique(time))

for (i in 1:ng)
{

dd1=dd[,,i] ## table of censor, ev, cr
dd2=apply(dd1,1,sum) ## sum for each time point
nrisk=sum(dd2)-cumsum(dd2)
nrisk=c(sum(dd2),nrisk[1:(length(nrisk)-1)])
dev=dd1[,2]
dcr=dd1[,3]
dall=dev+dcr
si=(nrisk-dall)/nrisk
s=cumprod(si)
sminus=c(1,s[1:(length(s)-1)])
fi=dev/nrisk*sminus
f=cumsum(fi)
fcri=dcr/nrisk*sminus
fcr=cumsum(fcri)
fcrminus=c(0,fcr[1:(length(fcr)-1)])
cp=f/(1-fcr)

t1i=dev*(1-fcr)^2
t2i=f^2*(dall-dev)
v1i=(t1i+t2i)/(nrisk*(nrisk-1))
v1i[nrisk==1]=0
v1=cumsum(v1i)
v=sminus^2*v1/((1-fcr)^4)

res=data.frame(time=tt,cif=f,fcr,cp,varCP=v,
group=rep(lg[i],length(v)))
res1=res[dall>0,]
if (i==1) cpvar=res1
else cpvar=list(cpvar,res1)
 
}

return(cpvar)
}