How to refine this distribution comparison script

Question

My data look like this

They are compound similarity estimations(the whole file is around 10GB). What I am trying to achieve is to compare the similarity distributions of each compound using the Kolmogorov-Smirnoff test so as to use the data for AIC calculation.

This is my code:

import pandas as pd
from scipy.stats import ks_2samp
import re

with open('ecfp4_file.csv', 'r') as f, open('Metrics.tsv', 'a') as f_out:
    f_out.write('compound_1' + '\t' + 'compound_2' + '\t' + 'Similarity' + '\t' + 'KS Distance' + '\n')
    df = pd.read_csv(f, delimiter = ',', lineterminator = '\n', header = 0)

    for m in range(10, 99):
        dc1 = []
        dc2 = []
        m1 = m/100
        m2 = m/100+0.009
        df2 = df.loc[(df['ECFP4'] >= m1) & (df['ECFP4'] < m2)]          #Work on different thresholds each time
    
        for i in range(0, df2.shape[0]):
            if len(df2.index[df2['id1'] == df2.iloc[i, 0]].tolist()) >= 5:
                dc1.append(df2.iloc[i, 0])
                dc2.append(df2.iloc[i, 1])
        for i in range(0, len(dc1)):
            dc1_l = []
            dc2_l = []
            df3 = df.loc[df['id1'] == dc1[i]]
            dc1_l = df3.iloc[:, 2].tolist()
            df4 = df.loc[df['id1'] == dc2[i]]
            dc2_l = df4.iloc[:, 2].tolist()
            x1 = re.findall(r"statistic=(.*)\,.*$", str(ks_2samp(dc1_l, dc2_l)))
            f_out.write(str(dc1[i]) + '\t' + str(dc2[i]) + '\t' + str(m1) + '\t' + str(x1[0]) + '\n')

which so far has taken 10 hours to run and it has not yet passed the 0.1 threshold, so it is highly unpractical. Does anyone know of a better way to do this? Or can help me refine my code for faster and better results?

I also tried to refine it and run it in R, but I get a vector memory exhausted (limit reached) error when I run it

install.packages("ggplot2")
library(ggplot2)
PS22_ECFP4 <- read.csv("./PS22_ECFP4.csv")

df <- data.frame(matrix(ncol = 4, nrow = 0))
cols <- c("id1", "id2", "Similarity", "KS Distance")
colnames(df) <- cols

for (i in PS22_ECFP4) {
  sim_df <- PS22_ECFP4[PS22_ECFP4$id1 == i,]
  dist1 <- sim_df[ , 3]
  for (k in 1:nrow(sim_df)) {
    sim <- sim_df[k , 3]
    sim_r <- floor(sim * 100) / 100
    comp2 <- sim_df[k, 2]
    comp2_df <- PS22_ECFP4[PS22_ECFP4$id1==comp2,]
    dist2 <- comp2_df[ , 3]
    if (length(dist1) > 1 & length(dist2) > 1) {
      ks <- ks.test(dist1, dist2)
      df[nrow(df) + 1,] = c(i, comp2, sim_r, ks$statistic)
    }
  }
}

write.csv(df, "./Metrics.csv", row.names = FALSE)

Answer 1

I am not sure how useful this answer is. But may be you could try this.

For your python script:

You use f_out.write within a loop. I think usually this I/O operation is resource intensive. I guess the speed would improve significantly if you store the string in a variable and write the output only once outside the for loop like this:

out_str = ''
for i in range(0, len(dc1)):
            dc1_l = []
            dc2_l = []
            df3 = df.loc[df['id1'] == dc1[i]]
            dc1_l = df3.iloc[:, 2].tolist()
            df4 = df.loc[df['id1'] == dc2[i]]
            dc2_l = df4.iloc[:, 2].tolist()
            x1 = re.findall(r"statistic=(.*)\,.*$", str(ks_2samp(dc1_l, dc2_l)))
            out_str +=  (str(dc1[i]) + '\t' + str(dc2[i]) + '\t' + str(m1) + '\t' + str(x1[0]) + '\n')

f_out.write(out_str)

Regarding the memory error in your R script, are you using 32bit build of R?

Please see this:

http://web.mit.edu/~r/current/arch/i386_linux26/lib/R/library/base/html/Memory-limits.html#

They say "Under most 64-bit versions of Windows the limit for a 32-bit build of R is 4Gb: for the oldest ones it is 2Gb. The limit for a 64-bit build of R (imposed by the OS) is 8Tb."

Conclusion: Please try using your python script, but use the write operation outside the loop.

Answer 2

If I was doing something like this in R (10G file), I'd use data.table

I can't quite work out what you're R code is doing (does it do what you want it to?)

I think the below might do what you want

library(data.table)
dat <- fread("./PS22_ECFP4.csv")
dat[,Similarity:=floor(ECFP4*100)/100]
dat[,'KS Distance':=lapply(1:nrow(dat),function(i)
  ks.test(dat[id1==dat$id1[i],ECFP4],dat[id2==dat$id2[i],ECFP4])$statistic  
)]
dat[,ECFP4:=NULL]
fwrite(dat,"./Metrics.csv")

Edit: The script loops through the data, and for each row creates two subsets of the data matching id1 and id2. Then does the ks.test on the similarity columns of the two subsets.