# Computation of concordant and discordant pairs

I would like to compute the number of concordant and discordant pairs of two large vectors. This function is a nice try, and works quite efficiently.

pairs <- function(x,y){
n <- length(x)
ix <- order(x)
x <- x[ix]
y <- y[ix]
Nc <- sum(sapply(1:(n-1),function(i) sum(x[i]<x[(i+1):n] & y[i]<y[(i+1):n])))
Nd <- sum(sapply(1:(n-1),function(i) sum(x[i]<x[(i+1):n] & y[i]>y[(i+1):n])))
return(list(Nc=Nc,Nd=Nd))
}

x <- runif(10000)
y <- runif(10000)
system.time(pairs(x,y))


Do you have any idea on how it would be possible to boost the function?

• – m-dz
Jul 6 '17 at 15:49
• You've already sorted x right? So you don't need to calculate x[i]<x[(i+1):n]? Removing this calculation will half your computation time (through my testing).
– Chi Pak
Jul 6 '17 at 15:50
• It's a clever idea, and will work with the aforementioned example, but let's imagine we have some ties in the data, x <- sample(1:10, 1000, replace = T);y <- sample(1:10, 1000, replace = T). In this case the number of concordant and discordant pairs would be not correctly counted.
– And
Jul 6 '17 at 16:11
• @And Removing dupes for that matter?
– M--
Jul 7 '17 at 19:58

These sorts of operations are tough to vectorize because you need to compare every element to all later elements. If the counting operation is a performance bottleneck in your code, then one possibility would be to implement it in C++ using the Rcpp package:

library(Rcpp)
cppFunction(
"IntegerVector CDcount(NumericVector x, NumericVector y) {
IntegerVector counts(2, 0);
int n = x.size();
for (int i=0; i < n; ++i) {
for (int j=i+1; j < n; ++j) {
counts += (x[i] < x[j]) && (y[i] < y[j]);
counts += (x[i] < x[j]) && (y[i] > y[j]);
}
}
return counts;
}")

pairs2 <- function(x, y) {
n <- length(x)
ix <- order(x)
x <- x[ix]
y <- y[ix]
counts <- CDcount(x, y)
return(list(Nc=counts, Nd=counts))
}


For your test data, this results in a 50 times speedup on my system.

x <- runif(10000)
y <- runif(10000)
system.time(pairs(x,y))
#    user  system elapsed
#   4.892   1.447   6.584
system.time(pairs2(x, y))
#    user  system elapsed
#   0.119   0.001   0.120
identical(pairs(x, y), pairs2(x, y))
#  TRUE


The Rcpp package takes some time to compile the CDcount function, so this option is probably only helpful if you are operating on very large vectors or if you are repeating the operation many times.