mquinodo / OFF-PEAK

CNV detection tool for WES data
4 stars 0 forks source link

Error in seq.default(min(correlation), 1, length.out = 71) : #3

Open jsl-netizen opened 1 week ago

jsl-netizen commented 1 week ago

Hi i have got that error message while trying to run the Rscript 03_OFF-PEAK.R on 58 sample. Rscript 03_OFF-PEAK.R \ --output output_directory \ --data output_directory/ALL.target.tsv \ --databasefile data/data-hg38.RData

Attaching package: ‘gplots’

The following object is masked from ‘package:stats’:

lowess

Type 'citation("pROC")' for a citation.

Attaching package: ‘pROC’

The following objects are masked from ‘package:stats’:

cov, smooth, var

Error in seq.default(min(correlation), 1, length.out = 71) : 'from' must be a finite number Calls: heatmap.2 -> duplicated -> seq -> seq.default Execution halted

the head of my excel ALL.target.tsv look like this.

Chr     begin   end     name    GC      P58     P59     P60     P61     P2      P4      P5      P6      P7      P8      P9      P10     P11     P12     P13     P14     P15     P16     P17     P18     P19     P20      P21     P22     P23     P24     P25     P26     P27     P28     P29     P30     P31     P32     P33     P34     P35     P36     P37     P38     P39     P40     P41     P42     P43     P44     P45     P46      P47     P48     P49     P50     P51     P52     P53     P54     P55     P56     P57
chr1    10000   59416   Off-target_chr1:10000-59416     0.455750        0       0       0       150     0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       151     0       0       0       198     0       248     0       0
chr1    59416   108832  Off-target_chr1:59416-108832,OR4F5_NM_001005484.2_exon1,OR4F5_NM_001005484.2_exon3,OR4F5_NM_001005484.2_exon2   0.395950        0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       151      0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0
chr1    108832  158248  Off-target_chr1:108832-158248   0.492850        0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0
chr1    158248  207666  Off-target_chr1:158248-207666   0.518700        0       133     0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       179     0       00       0       0       0       0       0       0       0       0       0       0       0       0       0
chr1    257666  297968  Off-target_chr1:257666-297968   0.409250        0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       109     0       0       89      123     0       0       0       0       0       0       135     0       0       00       0       0       0       0       163     149     0       0       0       0       0       0       198
chr1    347968  394973  Off-target_chr1:347968-394973   0.504200        0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0
chr1    394973  441978  Off-target_chr1:394973-441978   0.355900        181     0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0
chr1    441978  488983  Off-target_chr1:441978-488983,OR4F29_NM_001005221.2_exon1       0.365800        0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0
chr1    488983  535988  Off-target_chr1:488983-535988   0.441800        0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       0       00       0       0       0       0       0       0       0       0       0       0       0       0       0

I dont find the solution wihtout modifying the script is it normal or did i make a mistake

mquinodo commented 1 week ago

Hi, Thank you for using OFF-PEAK and for reporting the issue. I modified one line in 03_OFF-PEAK.R that should fix the bug. Can you try with the new script? Best, Mathieu

jsl-netizen commented 1 week ago

Hi thank you for your quick answer with the new script i get : Rscript 03_OFF-PEAK.R --output output_directory --data output_directory/ALL.target.tsv --databasefile data/data-hg38.RData

Attaching package: ‘gplots’

The following object is masked from ‘package:stats’:

lowess

Type 'citation("pROC")' for a citation.

Attaching package: ‘pROC’

The following objects are masked from ‘package:stats’:

cov, smooth, var

Error in seq.default(min(max(min(correlation), 0), 0.95), 1, length.out = 71) : 'from' must be a finite number Calls: heatmap.2 -> duplicated -> seq -> seq.default Execution halted

jsl-netizen commented 1 week ago

when id o this modification :

# Correlation between samples (heatmap and text output)
{
  # res is matrix for correlations
  num = dim(dataALL)[2]
  if (num - 4 <= 5) {
    stop(paste("You need to analyze at least 6 samples to run OFF-PEAK. You provided ", num - 4, " samples. Exit.", sep = ""))
  }
  correlation = matrix(nrow = (num - 4), ncol = (num - 4))
  colnames(correlation) = colnames(dataALL)[5:dim(dataALL)[2]]
  rownames(correlation) = colnames(dataALL)[5:dim(dataALL)[2]]

  # selecting targets with signal and noise limits
  meanALL = apply(dataALL[, 5:num], 1, mean)
  sdALL = apply(dataALL[, 5:num], 1, sd)
  ok = which(meanALL > minsignal & dataALL[, 1] != "chrX" & log10(sdALL / meanALL) < maxvar)

  # downsampling targets to 10'000 for quick computation
  if (length(ok) > 10000) {
    ok = ok[sample(1:length(ok), 10000, replace = FALSE)]
  }
  selection = dataALL[ok, 5:num]

  # computing correlations
  for (i in 1:(num - 4)) {
    for (j in 1:(num - 4)) {
      correlation[i, j] = cor(selection[, i], selection[, j])
    }
  }

  # Check and handle non-finite values in the correlation matrix
  correlation[!is.finite(correlation)] = 0

  # Calculate the minimum value of the correlation matrix after replacing non-finite values
  min_correlation = min(correlation, na.rm = TRUE)

  # plotting heatmap of correlation and writing output to text file
  pdf(file = paste(folder, "/01_general-stats/Heatmap-correlations-all.pdf", sep = ""), width = 7, height = 7)
  par(mfrow = c(1, 1))
  heatmap.2(
    correlation,
    cexRow = 0.1,
    cexCol = 0.1,
    breaks = seq(min_correlation, 1, length.out = 71),
    trace = "none",
    col = colorRampPalette(c("red", "orange", "yellow", "green"))(n = 70)
  )
  dev.off()
  correlation2 = cbind(rownames(correlation), correlation)
  correlation2 = rbind(colnames(correlation2), correlation2)
  write.table(correlation2, file = paste(folder, "/01_general-stats/Pairwise-correlations-all.tsv", sep = ""), quote = FALSE, sep = "\t", row.names = FALSE, col.names = FALSE)
}

The correlation[!is.finite(correlation)] = 0 line replaces any non-finite values (like NA, NaN, Inf, -Inf) in the correlation matrix with zero. The min_correlation = min(correlation, na.rm = TRUE) line calculates the minimum value of the correlation matrix after replacing non-finite values. The seq(min_correlation, 1, length.out = 71) ensures that the breaks parameter in the heatmap.2 function uses a finite minimum value.

the process go through but i dont know if i should change the line like that

jsl-netizen commented 1 week ago

I am in the process for testing with hg19 because on your publication it seems thats it works with it i will tell you

jsl-netizen commented 1 week ago

the tools works well with hg19 so i think maybe it is a problem with the hg38 build

mquinodo commented 1 week ago

Thank you for checking! In fact the script 03_OFF-PEAK.R is not different for hg19 and hg38 so I don't understand the issue. Could you add the following line at line 973: write.table(correlation,file=paste(folder,"/01_general-stats/correlation-test.tsv",sep=""),quote=F,sep="\t") Then have a look at the numbers and also at the plot 01_general-stats/Heatmap-correlations-all.pdf

jsl-netizen commented 1 week ago

thank you for your time,

I get this file i think something sketchy is happening in my bam. I will investigate this correlation-test.csv