I have two datasets (DataA and DataB) and I want to find the Spearman correlation between genes and also pull out the gene names (stored in first column of dataset) in R. I am using fread from read.table to read the file and cor.test to find Rho and p-value. The code works fine for a smaller data set of 10 obs, 10 variables. But, my actual data set is larger (DataA - 795 obj, 542 variables and DataB - 925 obs, 542 variables). When I try to run the code on my actual data, it is consuming a large amount of time (it has been 23 hrs since my code is running). Is there a way to optimize this code? Or is there a bug in the code?

A.data <- fread("DataA.csv", header = TRUE, data.table = FALSE, stringsAsFactors = FALSE)
B.data <- fread("DataB.csv", header = TRUE, data.table = FALSE, stringsAsFactors = FALSE)

correlation.result <- data.frame()

## Computing Correlation

for (i in 1:nrow(A.data)) {
  for (j in 1:nrow(B.data)) {
       correln <- cor.test(as.numeric(A.data[i,2:ncol(A.data)]), as.numeric(B.data[j,2:ncol(B.data)]), method = "spearman", exact = FALSE)
       rho <- correln$estimate
       p.val <- correln$p.value
       Gene <- paste(A.data[i,1], B.data[j,1], sep="_", collapse=" ") # To pullout gene name pairwise
       comb.data <- data.frame(rho, p.val, Gene)
       correlation.result <- rbind(correlation.result,comb.data)

Your code looks fine, but rbind to progressively build and resize a data frame can be gluttonous in terms of memory and time, especially given the A-row x B-row times you are doing this operation, which is fine for very small datasets, but not often realistic for real-world data.

Instead, perhaps cat results to a text file, where results are written in tab-delimited columns. You can then read in that text file later on.

In addition to serializing or saving your results to text, which make them easier to work with in non-R contexts, you won't have to deal with R's memory allocation problems when running rbind repeatedly. You just write out the answer to a file and read it in later on, in one pass.

Another approach is to preallocate memory to a presized data frame, and then write to that allocated space.

You know the size of your inputs ahead of time, which makes this approach tenable. In this case, it looks like you are storing two floats and a string (rho, p.val, and Gene) and that you are going to store A-row x B-row number of them.

Here's one way you might go about this:

# create a presized data frame
Arows <- nrow(A.data)
Brows <- nrow(B.data)
correlation.result <- data.frame(t(matrix(c(0.0, 0.0, "foo"), nrow=3, ncol=Arows * Brows)), stringsAsFactors=FALSE)
colnames(correlation.result) <- c("rho", "p.val", "Gene")

# assign values to that presized data frame 
Acols <- ncol(A.data)
Bcols <- ncol(B.data)
for (i in 1:Arows) {
  # get A data
  Ai <- as.numeric(A.data[i, 2:Acols])
  Aname <- A.data[i, 1]
  for (j in 1:Brows) {
    # calculate row offset in presized data frame
    ij <- (i - 1) * Brows + j
    # get B data
    Bj <- as.numeric(B.data[j, 2:Bcols])
    Bname <- B.data[j, 1]
    # run test
    correln <- cor.test(Ai, Bj, method = "spearman", exact = FALSE)
    rho <- correln$estimate
    p.val <- correln$p.value
    Gene <- paste(Aname, Bname, sep="_", collapse=" ")
    # assign value to row of result data frame
    correlation.result[ij,] <- c(rho, p.val, Gene)

This avoids R needing to resize and recopy the data frame contents as they grow, which I suspect is the major reason it is taking so long to run.

Chapter or "Circle" 2 of Patrick Burns' book The R Inferno discusses the problems with growing objects with rbind and ways to handle it — you might find this a useful read. This book is a PDF download via your search engine of choice.

  • $\begingroup$ Is there a reason for the downvote(s)? $\endgroup$ – Alex Reynolds Jan 9 '19 at 15:24
  • $\begingroup$ I think it is bad for a community when downvoting without any explanation. Please explain yourself. $\endgroup$ – Alex Reynolds Jan 11 '19 at 2:41

You can use corr.test() from psych package, e.g.

ct <- corr.test(t(A.data), t(B.data), method='spearman', adjust='none')

then ct$p and ct$r would give you (uncorrected) p-value and correlation matrices of size 795 x 925. Note that t() is required to compare matrices row-wise.

That would be orders of magnitude faster than loop-based approach.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.