Score a c-index

c-index.R

get_cindex <- function(
	predictions = NULL,
	test_truth = NULL
){
	
	start_time <- Sys.time()
	
	# Prep data
	predictions_names <- c("patientid","prediction" )
	if(!identical(predictions_names,names(predictions))) stop(paste(
		"predictions not tibble with",predictions_names
	))
	
	test_truth_names <- c("patientid", "timeto",    "had_event")
	if(!identical(test_truth_names,names(test_truth))) stop(paste(
		"test_truth not tibble with",test_truth_names
	))
	
	if (nrow(predictions) != nrow(test_truth)) 
		stop("predictions must have same length as test_truth")
	
	if(length(setdiff(test_truth$patientid,predictions$patientid)) != 0){
		stop("patients present in one but not the other")
	}
	
	inputdata <- dplyr::left_join(
		predictions  %>%
			# if they didn't rank (e.g. all 1), force it (not required, but we asked for a ranking)!
			arrange(prediction) %>%
			mutate(
				prediction = row_number()
			) %>%
			arrange(patientid), 
		test_truth %>%
			arrange(patientid),
		by = "patientid"
	)
	
	# sub function for scoring
	sub_score_cindex <- function(
		prediction = NULL,
		test_timeto = NULL,
		test_status = NULL
	){
		
		miss <- is.na(prediction) | is.na(test_timeto) | is.na(test_status)
		nmiss <- sum(miss)
		if (nmiss > 0) {
			miss <- !miss
			prediction <- prediction[miss]
			y <- y[miss]
			test_status <- test_status[miss]
		}
		n <- length(prediction)
		ne <- sum(test_status)
		storage.mode(prediction) <- "double"
		storage.mode(test_timeto) <- "double"
		storage.mode(test_status) <- "logical"
		z <- .Fortran(Hmisc:::F_cidxcn, prediction, test_timeto, test_status, length(prediction), nrel = double(1), 
									nconc = double(1), nuncert = double(1), c.index = double(1), 
									gamma = double(1), sd = double(1), as.logical(FALSE))
		r <- c(z$c.index, z$gamma, z$sd, n, nmiss, ne, z$nrel, z$nconc, 
					 z$nuncert)
		names(r) <- c("C Index", "Dxy", "S.D.", "n", "missing", "uncensored", 
									"Relevant Pairs", "Concordant", "Uncertain")
		
		return(r["C Index"])
	}
	
	# score!
	score <- sub_score_cindex(
			prediction = inputdata$prediction,
			test_timeto = inputdata$timeto,
			test_status = inputdata$had_event
		)
	
	list(
		c_index = score[["C Index"]],
		calculation_time = Sys.time() - start_time
	)
}