tzhu-bio / cisDynet

An integrated platform for modeling gene-regulatory dynamics and networks
MIT License
25 stars 3 forks source link

Peak2Gene bug when corr is not numeric #8

Open IStevant opened 5 months ago

IStevant commented 5 months ago

Hi,

The getPeak2Gene function fails to escape the z.test when the correlation is "NA":

      if (is.numeric(corr)){
      p <- BSDA::z.test(corr_random, mu = corr, sigma.x = 15)$p.value
      } else{
        corr = 0
        p <- BSDA::z.test(corr_random, mu = corr, sigma.x = 15)$p.value
      }

https://github.com/tzhu-bio/cisDynet/blob/c5b99d37f12e7d2d8ef1fde766e37efd639f6904/R/Peak2Gene.R#L50C1-L54C8

Currently, when the correlation is "NA", it returns an error and stop the calculation:

2024-04-29 11:59:33 Remove the gene with all expression value is 0.
2024-04-29 11:59:35 Make the pseudo data for permutation...
[===============>-------------------------------------------------------]  22%Error in BSDA::z.test(corr_random, mu = corr, sigma.x = 15) :
  mu must be a single number
Calls: getPeak2Gene ... getGenePeaks -> calculate_pvalue -> apply -> FUN -> <Anonymous>

I modified the function as following to make it run:

      if (is.na(corr)){
        corr <- 0
        p <- 1
      } else {
        p <- BSDA::z.test(corr_random, mu = corr, sigma.x = 15)$p.value
      }
tzhu-bio commented 5 months ago

Thank you for the help you provided. But please note that the calculated corelation is NA is usually a problem with your matrix, so please be careful to check the contents of the matrix rather than simply changing it to corr <- 0.

IStevant commented 5 months ago

I analysed my bulk ATAC and RNA-seq outside of cisDynet, and I don't have the exact same number of replicates between the RNA and ATACseq, so I built matrices with the mean of expression and mean signal of ATAC per condition. By chance it happened that there is no change in expression, nor change in accessibility for a particular comparison, thus returning "NA" because no correlation can be computed, and the original function just stop, while I would like it to continue, that's why I modified the fuction as mentioned above.

tzhu-bio commented 5 months ago

Okay, I understand.