https://github.com/cran/Hmisc
Tip revision: 256ea0f4576fef9b8641bd64e4d800ac668b805f authored by Frank E Harrell Jr on 29 November 2020, 05:00:16 UTC
version 4.4-2
version 4.4-2
Tip revision: 256ea0f
cidxcn.f
C Output from Public domain Ratfor, version 1.01
C------------------------------------------------------------------------
C Compute c-index (c) and Brown-Hollander-Krowar-Goodman-Kruskal-Somer
C rank correlation (gamma) between X and Y with censoring indicator E
C Also returns number of relevant, concordant, and uncertain pairs
C (nrel, nconc, nuncert) and estimated s.d. of gamma (sd) using
C Quade formula (see SAS PROC MATPAR). Pairs with tied x are
C excluded if outx=.TRUE.
C
C F. Harrell 27Nov90
C Modification of SAS Procedure KGKC (1980)
C-------------------------------------------------------------------------
subroutine cidxcn(x,y,e,n,nrel,nconc,nuncert,c,gamma,sd,outx)
implicit double precision (a-h,o-z)
double precision x(n),y(n),dx,dy
logical e(n),outx
double precision nrel,nuncert,nconc
nconc=0d0
nrel=0d0
nuncert=0d0
sumr=0d0
sumr2=0d0
sumw=0d0
sumw2=0d0
sumrw=0d0
do23000 i=1,n
wi=0d0
ri=0d0
do23002 j=1,n
if(j.ne.i)then
dx=x(i)-x(j)
dy=y(i)-y(j)
if(dx.ne.0. .or. .not.outx)then
if((e(i).and.dy.lt.0.).or.(e(i).and..not.e(j).and.dy.eq.0.))then
if(dx.lt.0.)then
nconc=nconc+1d0
wi=wi+1d0
else
if(dx.eq.0.)then
nconc=nconc+.5d0
else
wi=wi-1d0
endif
endif
nrel=nrel+1d0
ri=ri+1d0
else
if((e(j).and.dy.gt.0.).or.(e(j).and..not.e(i).and.dy.eq.0.))then
if(dx.gt.0.)then
nconc=nconc+1d0
wi=wi+1d0
else
if(dx.eq.0.)then
nconc=nconc+.5d0
else
wi=wi-1d0
endif
endif
nrel=nrel+1d0
ri=ri+1d0
else
if(.not.(e(i).and.e(j)))then
nuncert=nuncert+1d0
endif
endif
endif
endif
endif
23002 continue
23003 continue
sumr=sumr+ri
sumr2=sumr2+ri*ri
sumw=sumw+wi
sumw2=sumw2+wi*wi
sumrw=sumrw+ri*wi
23000 continue
23001 continue
c=nconc/nrel
gamma=2.*(c-.5)
Ccall dblepr('sumr',4,sumr,1)
Ccall dblepr('sumw',4,sumw,1)
Ccall dblepr('sumr2',5,sumr2,1)
Ccall dblepr('sumw2',5,sumw2,1)
Ccall dblepr('sumrw',5,sumrw,1)
sd=sumr2*sumw**2-2d0*sumr*sumw*sumrw+sumw2*sumr**2
sd=2.*dsqrt(sd)/sumr/sumr
return
end