mskilab-org / fragCounter

GC and mappability corrected fragment coverage for paired end whole genome sequencing
MIT License
7 stars 11 forks source link

run test data and get an DNACopy error #13

Closed dagezh closed 1 year ago

dagezh commented 2 years ago

Hi, i run the test data and get an error. The command line is "frag -b /usr/local/lib/R/site-library/fragCounter/extdata/chr21.bam -d /usr/local/lib/R/site-library/fragCounter/extdata/gcMAP21 -w 200 -o ./"

i got an output:

Correcting coverage at individual scales Correcting coverage at 200 bp scale, with 234722 intervals Quantile filtering response and covariates Response min quantile: 0 max quantile: 0 Nominated 0 of 234722 data points for loess fitting Converting to GRanges Made GRanges Warning messages: 1: replacing previous import ‘GenomicRanges::shift’ by ‘data.table::shift’ when loading ‘bamUtils’ 2: replacing previous import ‘GenomicAlignments::second’ by ‘data.table::second’ when loading ‘bamUtils’ 3: replacing previous import ‘GenomicAlignments::last’ by ‘data.table::last’ when loading ‘bamUtils’ 4: replacing previous import ‘GenomicAlignments::first’ by ‘data.table::first’ when loading ‘bamUtils’ 5: In for (what in notDone) installClassMethod(get(what, envir = thisClass@refMethods), : closing unused connection 3 (/usr/local/lib/R/site-library/fragCounter/extdata/chr21.bam) 6: In value[3L] : DNACopy error moving on without segmenting 7: In FUN(as.data.frame(grs[[i]]), fields, seg = seg) : Not enough samples for loess fitting - check to see if missing or truncated data?


cov.rds is a truncated file whose size is 866K and can't be load in R.

when test my own data I also got an error.

bam.cov.tile.st test.bam chunk 21 num fragments processed 2.1e+07 mean cov: 7 per bin, estimated tot fragments: 21.85 million fragments, processing 371543.9 fragments/second time elapsed: 0.02 hours, estimated time remaining: 0 hours , estimated total time 0.02 hours Finished computing coverage, and making GRanges Finished acquiring coverage Loaded GC and mappability length cov is 3137207, length gc is 3074148, length map is 3074148 Synced coverage, GC, and mappability Modified gc / mappability correction Converting to data.table Grouping intervals Defining 100 fold collapsed ranges Presegmenting at NA bp scale Error in [.data.table(cov.dt, , list(chr = seqnames[1], start = min(start), : 'by' appears to evaluate to column names but isn't c() or key(). Use by=list(...) if you can. Otherwise, by=eval(get(paste("lev", numlevs, sep = ""))) should work. This is for efficiency so data.table can detect which columns are needed. Calls: fragCounter -> multicoco -> seg2gr -> is -> [ -> [.data.table In addition: Warning messages: 1: replacing previous import ‘GenomicRanges::shift’ by ‘data.table::shift’ when loading ‘bamUtils’ 2: replacing previous import ‘GenomicAlignments::second’ by ‘data.table::second’ when loading ‘bamUtils’ 3: replacing previous import ‘GenomicAlignments::last’ by ‘data.table::last’ when loading ‘bamUtils’ 4: replacing previous import ‘GenomicAlignments::first’ by ‘data.table::first’ when loading ‘bamUtils’ 5: In for (i in seq_along(defined)) { : closing unused connection 3 (test.bam) 6: In max(width(cov)) : no non-missing arguments to max; returning -Inf 7: In cat("Presegmenting at ", as.integer(WID * base^(numlevs)), " bp scale \n") : NAs introduced by coercion to integer range Execution halted

williamphu commented 1 year ago

Any updates on this issue? I've also been experiencing the same problem with the test data.

zining01 commented 1 year ago

Hi, apologies for the delay and thanks for bringing this to our attention. This issue seems to have happened because the example chr21.bam was empty. On our dev branch I have uploaded a new chr21.bam which has nonzero coverage for a small 100 Kbp region of chr 21 ("21:30000000-30100000"). I was able to run the new test example from my end... maybe you could clone and reinstall from the dev branch and retry?

williamphu commented 1 year ago

Thanks - I think that works; I was able to run the test command with the following warnings at the end:

Warning messages:
1: replacing previous import ‘GenomicRanges::shift’ by ‘data.table::shift’ when loading ‘bamUtils’
2: replacing previous import ‘GenomicAlignments::second’ by ‘data.table::second’ when loading ‘bamUtils’
3: replacing previous import ‘GenomicAlignments::last’ by ‘data.table::last’ when loading ‘bamUtils’
4: replacing previous import ‘GenomicAlignments::first’ by ‘data.table::first’ when loading ‘bamUtils’
5: closing unused connection 3 (/home/wip4002/github/fragCounter/inst/extdata/chr21.bam)
6: In FUN(as.data.frame(grs[[i]]), fields, seg = seg) :
  Not enough samples for loess fitting - check to see if missing or truncated data?

Does this match your output? If so, then I think the issue has been resolved.

zining01 commented 1 year ago

Yes, this matches the output that we get. Thank you for helping us to figure out this problem!