Hi,

Can i use the package in Rstudio of windows? My R version for R studio is 3.4

Hi,

Just a suggestion for a future update of the package. It would be nice if there would exist a function that could generate read counts directly from CRAM files. For us, it would be a great help as now we have to convert our CRAM files back to BAM files to run ExomeDepth.

Thank you for considering it.

All the best,
Toon

Hi

I was wondering if there is currently a way to limit the resources of ExomeDepth? When I run ExomeDepth on our server, it takes all the available resources, even if I indicate that my PBS jobscript may use only just 1 core... Inserting "ulimit -c 1" does not work either.

All the best
Toon

Hi,
do you have an idea how I could calculate the BF for all exons (that they are not a cnv/ are a cnv)?
You wrote (p. 7):

It is in fact the log10 of the likelihood ratio of data for the CNV call divided by the null (normal copy number)

I see that the all.exons object has a member "likelihood", "expected", "reference", "test". Is it possible to get the BF scores from these values?

Thank you,
Vinzenz

Hi @vplagnol,
i switched to the new version of ExomeDepth and saw that when getBamCounts finish to parse each bam file, i get a warning that was partially discussed in issue #22 but i think this specific problem was lost as some point in the thread.
What i wanted to understand better is what 'x' in this warning is referring to, the bed file? I'm using GRCh38 reference with decoy, my bed file contains contigs from chr1 to chrX-Y, but the x list in the warnings stops at 22 (when chrX and chrY are parsed with getBamCounts).
However, read counts at the end are formatted properly and are not empty, so can i skip this warnings or is it better to do some further test to be sure?

Packages versions:

packageVersion("ExomeDepth")
[1] ‘1.1.15’
packageVersion("GenomicRanges")
[1] ‘1.40.0’

Command used:

my.bam = c(my.ctrl,my.cases)
my.counts <- getBamCounts(bed.file = bed_file,
bam.files = my.bam,
include.chr = FALSE,
referenceFasta = reference)

The warning i get:

50: In .Seqinfo.mergexy(x, y) :
Each of the 2 combined objects has sequence levels not in the other:

• in 'x': 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22
• in 'y': chrM, chr1_KI270706v1_random, chr1_KI270707v1_random, chr1_KI270708v1_random, chr1_KI2709v1_random, chr1_KI270710v1_random, chr1_KI270711v1_random, chr1_KI270712v1_random, chr1_KI270713vrandom, chr1_KI270714v1_random, chr2_KI270715v1_random, chr2_KI270716v1_random, chr3_GL000221v1_ranm, chr4_GL000008v2_random, chr5_GL000208v1_random, chr9_KI270717v1_random, chr9_KI270718v1_random, r9_KI270719v1_random, chr9_KI270720v1_random, chr11_KI270721v1_random, chr14_GL000009v2_random, chr_GL000225v1_random, chr14_KI270722v1_random, chr14_GL000194v1_random, chr14_KI270723v1_random, chr1KI270724v1_random, chr14_KI270725v1_random, chr14_KI270726v1_random, chr15_KI270727v1_random, chr16I270728v1_random, chr17_GL000205v2_random, chr17_KI270729v1_random, chr17_KI270730v1_random, chr22_270731v1_random, chr22_KI270732v1 [... truncated]

Many thanks,
Paola

Dear all,

is it possible to convert ExomeDepth output to a VCF?
I looked in the manual, but there is not reference to a VCF...

thanks,

Hi,

I'm using exomedepth (1.1.10) with great satisfaction, however I notice that the plotting of a locus is a bit suboptimal, and I'd like to change it.

In the documentation you show a nice plot:

What I get is:

The main difference is the gene bodies in the lower section; you have a clear genebody, where I get only the individual exon locations. With denser genes, this section becomes really cluttered. I'm trying to get to the plotting code but cannot find where this is done.

Can you point me on how to get the genebodies plotted as in your example?
thanks & regards, Ies Nijman

I had a question regarding how you are modeling the effect of duplications and deletions on expected copy number that I could not determine from your paper, I was hoping you could help.

In the code, the get_loglike_matrix function calculates the expected log likelihood for the default copy number of 2 in a straightforward way.

https://github.com/evolvedmicrobe/ExomeDepth/blob/master/src/CNV_estimate.cpp#L75-L77

However, I do not understand how you are changing this estimate to account for copy numbers of 1 and 3, and in particular the role of the bestSd variable. Based on the paper and my understanding of the code, it seems that if a target is deleted, the expectation changes from X / (X + Y) into .5X / (.5X + Y), this seems like a straightforward way to change the expected proportion of reads in the test set under the deletion hypothesis.

However, this change to the expectation seems to only occur in one place and there is a second variable, bestSd that is not updated.

https://github.com/evolvedmicrobe/ExomeDepth/blob/master/src/CNV_estimate.cpp#L73

Since the a1 and a2 variables used to calculate the likelihood are derived from the expected proportion and overdispersion parameters, shouldn't they also be updated to reflect the difference in the expected proportion? Any help understanding this appreciated.

Hello everyone. I'm trying to use ExomeDepth for my CNV analysis but I'm falling on this error:

"RangedData objects are deprecated. Please migrate your code to use
 GRanges or GRangesList objects instead."

I'm using the most recent version of ExomeDepth. The vignette says: "Version 1.1.12: •Updated for GRanges changes." but the error persists.

Does anyone has a fix for this issue?

Installed with R4.1.0.

Hello, I am trying to install ExomeDepth (install_github("vplagnol/ExomeDepth")) in order to install and use DECoN but unfortunately, it is not working.
I use the R version 4.2.1

Does anyone already has this issue? Thank you in advance.
Sincerely,
Lou

> install_github("vplagnol/ExomeDepth")
Downloading GitHub repo vplagnol/ExomeDepth@HEAD
Skipping 5 packages not available: GenomicAlignments, GenomicRanges, Rsamtools, IRanges, Biostrings
v  checking for file '/tmp/RtmpZeBznm/remotes1179f51e5699a/vplagnol-ExomeDepth-80da0cb/DESCRIPTION' ...
-  preparing 'ExomeDepth':
v  checking DESCRIPTION meta-information ...
-  cleaning src
-  checking for LF line-endings in source and make files and shell scripts
-  checking for empty or unneeded directories
-  looking to see if a 'data/datalist' file should be added
-  building 'ExomeDepth_1.1.15.tar.gz'
   
Installing package into '/DECoN/Linux/renv/library/R-4.2/x86_64-pc-linux-gnu'
(as 'lib' is unspecified)
* installing *source* package 'ExomeDepth' ...
** using staged installation
** libs
** arch - 
R CMD SHLIB hmm.cpp fdiv.c stream.c hypot.c error.c ldfrexp.c VP_infnan.c elementary.c invhyp.c log1p.c VP_log.c exp.c math.c VP_zeta.c VP_trig.c VP_psi.c VP_gamma.c beta.c CNV_estimate.cpp -o ExomeDepth.so 
make[1]: Entering directory '/tmp/RtmpSq3pBX/R.INSTALL133fd341b3e96/ExomeDepth/src'
g++ -std=gnu++14 -I"/usr/share/R/include" -DNDEBUG      -fpic  -g -O2 -ffile-prefix-map=/build/r-base-wSz8CH/r-base-4.2.1=. -flto=auto -ffat-lto-objects -flto=auto -ffat-lto-objects -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2  -c hmm.cpp -o hmm.o
gcc -I"/usr/share/R/include" -DNDEBUG      -fpic  -g -O2 -ffile-prefix-map=/build/r-base-wSz8CH/r-base-4.2.1=. -flto=auto -ffat-lto-objects -flto=auto -ffat-lto-objects -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2  -c fdiv.c -o fdiv.o
In file included from fdiv.c:21:
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:85:1: error: expected ';' before 'extern'
   85 | __MATHCALL_VEC (acosh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:87:1: error: expected ';' before 'extern'
   87 | __MATHCALL_VEC (asinh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:89:1: error: expected ';' before 'extern'
   89 | __MATHCALL_VEC (atanh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:119:1: error: expected ';' before 'extern'
  119 | __MATHCALL_VEC (expm1,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:122:1: error: expected ';' before 'extern'
  122 | __MATHCALL_VEC (log1p,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:147:1: error: expected ';' before 'extern'
  147 | __MATHCALL_VEC (hypot,, (_Mdouble_ __x, _Mdouble_ __y));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:85:1: error: expected ';' before 'extern'
   85 | __MATHCALL_VEC (acosh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:87:1: error: expected ';' before 'extern'
   87 | __MATHCALL_VEC (asinh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:89:1: error: expected ';' before 'extern'
   89 | __MATHCALL_VEC (atanh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:119:1: error: expected ';' before 'extern'
  119 | __MATHCALL_VEC (expm1,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:122:1: error: expected ';' before 'extern'
  122 | __MATHCALL_VEC (log1p,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:147:1: error: expected ';' before 'extern'
  147 | __MATHCALL_VEC (hypot,, (_Mdouble_ __x, _Mdouble_ __y));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:85:1: error: expected ';' before 'extern'
   85 | __MATHCALL_VEC (acosh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:87:1: error: expected ';' before 'extern'
   87 | __MATHCALL_VEC (asinh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:89:1: error: expected ';' before 'extern'
   89 | __MATHCALL_VEC (atanh,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:119:1: error: expected ';' before 'extern'
  119 | __MATHCALL_VEC (expm1,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:122:1: error: expected ';' before 'extern'
  122 | __MATHCALL_VEC (log1p,, (_Mdouble_ __x));
      | ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:147:1: error: expected ';' before 'extern'
  147 | __MATHCALL_VEC (hypot,, (_Mdouble_ __x, _Mdouble_ __y));
      | ^~~~~~~~~~~~~~
make[1]: *** [/usr/lib/R/etc/Makeconf:168: fdiv.o] Error 1
make[1]: Leaving directory '/tmp/RtmpSq3pBX/R.INSTALL133fd341b3e96/ExomeDepth/src'
make: *** [Makefile:15: CNV_estimate.so] Error 1
ERROR: compilation failed for package 'ExomeDepth'
* removing '/DECoN/Linux/renv/library/R-4.2/x86_64-pc-linux-gnu/ExomeDepth'
Warning message:
In i.p(...) :
  installation of package '/tmp/RtmpZeBznm/file1179f56ce92d5/ExomeDepth_1.1.15.tar.gz' had non-zero exit status

Hi,

I was running Exomedepth package for cnv detection from germline CNV's, its running fine till RData file creation but it is throwing an error when its running CallCNVs function, the error is:

Now computing the likelihood for the different copy number states
Error in CallCNVs(x = all.exons, transition.probability = 10^-4, chromosome = refcounts.dafr$space, :
Chromosome, start and end vector must have the same lengths.

Can you please explain what this error means? I am not able to debug this error, this will be very helpful.

The documentation for select.reference set defines bin.length as "Length (in bp) of each of the regions" (my emphasis).

However, the vignette gives as an example for bin.length:

bin.length = (ExomeCount.dafr$end - ExomeCount.dafr$start)/1000

... which looks to me to be the length in kb rather than bp. Which should it be, or does it not matter as long as it is proportional for all of the regions?

Also-- since ExomeCount.dafr is derived from a GRanges object, I think these are 1-based inclusive coordinates, so the bin.length would be more accurately ExomeCount.dafr$end - ExomeCount.dafr$start + 1) (adding 1), right?

hey,
Does this tool support to calculate the absoulte copy number ?
an if yes, how to get that value ?

thank you !

Dear Vincent,
I am using ExomeDepth for a big batch of BAM file.
I am now having this error:

It looks like the test samples has only 2 bins with more than 5 reads. The coverage is too small to perform any meaningful inference so no likelihood will be computed.
Error in if ((i > 2) && (res.data.frame$mean.p[i] < 0.05)) break :
missing value where TRUE/FALSE needed**

How can I deal with it? The optimization worked for all other files but not for this BAM.

Best Regards,
Francesco

Dear Vincent,

I am facing a problem with fasta files

bis1 <- getBamCounts(bed.frame = regions.df,
bam.files = file.path("./BAMs_all","qG16018101bis.final.bam"),
include.chr = FALSE,
referenceFasta = "Z:/RNA_seq/Ref_Genomes_fa/hg19.fa")
Reference fasta file provided so ExomeDepth will compute the GC content in each window
Error in value[3L] : record 1 (1:1138878-1139350) failed
file: Z:/RNA_seq/Ref_Genomes_fa/hg19.fa

I have tried also with hg38.fa (both fasta files are being used by other methods with no problems.
If I ommit the fasta files it works but results are not what I expect so I am trying to make the method adjust by GC content as a first trial.
Best,

Lara

sessionInfo()

R version 3.4.3 (2017-11-30)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 7 x64 (build 7601) Service Pack 1

Matrix products: default

locale:
[1] LC_COLLATE=Spanish_Spain.1252 LC_CTYPE=Spanish_Spain.1252 LC_MONETARY=Spanish_Spain.1252 LC_NUMERIC=C
[5] LC_TIME=Spanish_Spain.1252

attached base packages:
[1] parallel stats4 stats graphics grDevices utils datasets methods base

other attached packages:
[1] GenomicRanges_1.30.0 GenomeInfoDb_1.14.0 IRanges_2.12.0 S4Vectors_0.16.0 BiocGenerics_0.24.0 XLConnect_0.2-14
[7] XLConnectJars_0.2-14 ExomeDepth_1.1.10

loaded via a namespace (and not attached):
[1] XVector_0.18.0 splines_3.4.3 zlibbioc_1.24.0 GenomicAlignments_1.14.1 BiocParallel_1.12.0
[6] lattice_0.20-35 tools_3.4.3 SummarizedExperiment_1.8.1 grid_3.4.3 Biobase_2.38.0
[11] matrixStats_0.52.2 yaml_2.1.16 Matrix_1.2-12 rJava_0.9-9 GenomeInfoDbData_1.0.0
[16] bitops_1.0-6 RCurl_1.95-4.10 aod_1.3 VGAM_1.0-5 DelayedArray_0.4.1
[21] compiler_3.4.3 Rsamtools_1.30.0 Biostrings_2.46.0

Hi,

I am trying to get the counts in my data but when I run getBamCounts I have the followed error:
Error in countBamInGRanges.exomeDepth(bam.file = bam, index = index, granges = target, :
Some sequences in the target data frame cannot be found in the index of the BAM file

Can someone help understand why this happen?

Thanks

Dear Vincent,

ExomeDepth can end up with only one control after the 'correlation' search.
Having a single sample as a control can be problematic.
Would it be possible to add a minimal number of controls as an option ?

I downloaded it from the GitHub Exomedepth and got all dependencies sorted out. now I am having an issue with the processing of the ExomeDepth vignette.

It gives me the following errors. Can it be sorted out.

These are the last lines when I try make build.
###############################################
R CMD build --resave-data working/ExomeDepth

• checking for file ‘working/ExomeDepth/DESCRIPTION’ ... OK
• preparing ‘ExomeDepth’:
• checking DESCRIPTION meta-information ... OK
• cleaning src
• installing the package to build vignettes
• creating vignettes ... ERROR
  --- re-building ‘ExomeDepth-vignette.Rnw’ using knitr
  Error: processing vignette 'ExomeDepth-vignette.Rnw' failed with diagnostics:
  Running 'texi2dvi' on 'ExomeDepth-vignette.tex' failed.
  Messages:
  /usr/bin/texi2dvi: TeX neither supports -recorder nor outputs \openout lines in its log file
  --- failed re-building ‘ExomeDepth-vignette.Rnw’

SUMMARY: processing the following file failed:
‘ExomeDepth-vignette.Rnw’

Error: Vignette re-building failed.
Execution halted
make: *** [Makefile:7: build] Error 1

I really appreciate any help you can provide.

Hi,

A couple of questions regarding the ExomeDepth correlation threshold 0.97:

• How did you go about getting the threshold of 0.97 for correlation of read counts between the test and reference samples? As the vignette says, if the correlation is lower than 0.97, “consider the output of ExomeDepth as less reliable (i.e. most likely a high false positive rate)”.

• Do you have any data on how well ExomeDepth performs with correlations less than 0.97, in terms of sensitivity and specificity? In other words, at what point does the performance really start to taper off?

Thanks very much,
Rajat

Dear creators of ExomeDepth,

I am facing very curious Error form my Agilent WES BAM files.

Can anybody tell me why I have ths kind of error?

Error in CallCNVs(x = all_exons, transition.probability = 10^-4, chromosome = countdf$chromosome,  : 
  Chromosome, start and end vector must have the same lengths.
Calls: CallCNVs -> CallCNVs
In addition: Warning messages:
1: In aod::betabin(data = data.for.fit, formula = as.formula(formula),  :
  The data set contains at least one line with weight = 0.

2: In model.matrix.default(mt, mfb, contrasts) :
  non-list contrasts argument ignored
3: In aod::betabin(data = data.for.fit, formula = as.formula(formula),  :
  The data set contains at least one line with weight = 0.

4: In model.matrix.default(mt, mfb, contrasts) :
  non-list contrasts argument ignored
5: In aod::betabin(data = data.for.fit, formula = as.formula(formula),  :
  The data set contains at least one line with weight = 0.

6: In model.matrix.default(mt, mfb, contrasts) :
  non-list contrasts argument ignored
7: In model.matrix.default(mt, mfb, contrasts) :
  non-list contrasts argument ignored
8: In model.matrix.default(mt, mfb, contrasts) :
  non-list contrasts argument ignored
9: In aod::betabin(data = data.for.fit, formula = as.formula(formula),  :
  The data set contains at least one line with weight = 0.

10: In model.matrix.default(mt, mfb, contrasts) :
  non-list contrasts argument ignored
Execution halted

Thank you so much for any help,

Paul.

Hello,
at first thank you for very good tools for CNVs. I just analyzed a few WES and I found negative BF value at some positions.
How can I this interpreting?
Example:

BF= -75
READS_EXPECTED= 500
READ_OBSERVED= 750
RATIO= 1.44

And status is "deletion". 

Thank you very much for help with interpretation of this result.

Paul.

Whenever I try to add a reference fasta file for computing the GC content, i get an error like this
Reference fasta file provided so ExomeDepth will compute the GC content in each window Error in value[[3L]](cond) : record 38397 (chr11:134201957-134202041) failed file: data/hg19.fa
I have tried both with the UCSC and Ensembl build but I get the same error on a different chromosome.

R version 3.6.3 (2020-02-29)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 18363)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252    LC_MONETARY=English_United States.1252
[4] LC_NUMERIC=C                           LC_TIME=English_United States.1252    

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] GenomicRanges_1.38.0 GenomeInfoDb_1.22.1  IRanges_2.20.2       S4Vectors_0.24.4     BiocGenerics_0.32.0 
[6] ExomeDepth_1.1.15    readr_1.3.1         

Hi,

While generating 6 different reference pools (for gender, sequencers and target kits) we have run into the same error every time. It looks like this:

*** caught segfault ***
address 0xffffffff8f200760, cause 'memory not mapped'

Traceback:
 1: .Call(.bamfile_open, path(con), index, "rb")
 2: doTryCatch(return(expr), name, parentenv, handler)
 3: tryCatchOne(expr, names, parentenv, handlers[[1L]])
 4: tryCatchList(expr, classes, parentenv, handlers)
 5: tryCatch({    .io_check_exists(path(con))    index <- sub("\\.bai$", "", index(con, asNA = FALSE))    con$.extptr <- .Call(.bamfile_open, path(con), index, "rb")}, error = function(e
 6: open.BamFile(BamFile(file, index = index, asMates = asMates),     "rb")
 7: open(BamFile(file, index = index, asMates = asMates), "rb")
 8: .open_BamFile(file, index = index, param = param)
 9: GenomicAlignments::readGAlignments(file = bam.file, index = index,     param = my.param.single)
10: GenomicAlignments::readGAlignments(file = bam.file, index = index,     param = my.param.single)
11: countBamInGRanges.exomeDepth(bam.file = bam, index = index, granges = target,     min.mapq = min.mapq, read.width = read.width)
12: getBamCounts(bed.file = opt$target, bam.files = bams, include.chr = opt$chr,     referenceFasta = opt$reference, read.width = opt$width, )
An irrecoverable exception occurred. R is aborting now ...

This seems to happen at random times, as it occurred for example after having parsed 104 input samples and a runtime of 14h19m, but also after having parsed 214 input samples and a runtime of 33h15mins.

We can confirm that this error is not a RAM limitation.

Furthermore we would like your input regarding the size of reference pools. Since ExomeDepth selects the best matching samples from the reference pool, larger pools would always be better as they increase the chance of finding good matches. Is that interpretation correct?

Thanks!

Hi,

What is the default min MQ value is when running getBamCounts? I've tried setting it to 20 and to 1, but the results are not the same as when I run it without specifying the min.mapq argument.

Thank you,
Josianne

Hello,

I was hoping to discuss this specific option in the getBamCounts function, as it is generally not mentioned in the paper or in the vignette.

The description of the read.width option in the reference manual is:

numeric, maximum distance between the side of the target region and the middle
of the paired read to include the paired read into that region.

As this is essentially a wrapper around the countBamInGRanges function from exomeCopy, here is its description from exomeCopy:

The width of a read, used in counting overlaps of mapped reads with the genomic
ranges. The default is 1, resulting in the counting of only read starts in
genomic ranges. If the length of fixed width reads is used, e.g. 100 for 100bp
reads, then the function will return the count of all overlapping reads with the genomic
ranges. However, counting all overlapping reads introduces dependency
between the counts in adjacent windows.

It also says:

With the default setting (read.width=1), only the read starts are used for counting purposes (the
leftmost position regardless of the strandedness of the read).

My questions/remarks are:

1. There seems to be discordance as to whether the leftmost position or the middle of the read is considered when measuring the distance from a range.
2. By using 300 as a default, if two adjacent targets are less that 600 bp apart, then reads might be counted twice, once for each adjacent target. Is this assumption correct? It seems quite important but is never mentioned anywhere else except the vignette.

Thanks a lot in advance!

Hi @vplagnol , I am a bioinformatician and I am facing some troubles while running ExomeDepth.

I am actually trying to use ExomeDepth to call CNVs on Panel data of 24 mitochondrial patients (our gene panels are very big, almost resuming a real Exome).

These are the commands I use, following ExomeDepth official R vignettes:

-library(ExomeDepth)

-my.counts <- getBamCounts(bed.file = "panel.bed", bam.files = c("MT1.bam", "MT2_S12.bam", "MT3.bam", "MT4.bam", "MT5.bam", "MT6.bam", "MT7.bam", "MT8.bam", "MT9.bam", "MT10.bam", "MT11.bam", "MT12.bam", "MT13.bam", "MT14.bam", "MT15.bam", "MT16.bam", "MT17.bam", "MT18.bam", "MT19.bam", "MT20.bam", "MT21.bam", "MT22.bam", "MT23.bam", "MT24.bam"), include.chr = TRUE, referenceFasta = "hg19.fa")

-ExomeCount.dafr <- as(my.counts, 'data.frame')

-ExomeCount.dafr$chromosome <- gsub(as.character(ExomeCount.dafr$chromosome), pattern = 'chr', replacement = '')

-my.test <- my.counts$MT1.bam

-my.ref.samples <- c("MT2_S12.bam", "MT3.bam", "MT4.bam", "MT5.bam", "MT6.bam", "MT7.bam", "MT8.bam", "MT9.bam", "MT10.bam", "MT11.bam", "MT12.bam", "MT13.bam", "MT14.bam", "MT15.bam", "MT16.bam", "MT17.bam", "MT18.bam", "MT19.bam", "MT20.bam", "MT21.bam", "MT22.bam", "MT23.bam", "MT24.bam")

-my.reference.set <- as.matrix(ExomeCount.dafr[, my.ref.samples])

-my.choice <- select.reference.set (test.counts = my.test, reference.counts = my.reference.set, bin.length = (ExomeCount.dafr$end -

-ExomeCount.dafr$start)/1000, n.bins.reduced = 10000)

-my.matrix <- as.matrix( ExomeCount.dafr[, my.choice$reference.choice])

-my.reference.selected <- apply(X = my.matrix, MAR = 1, FUN = sum)

-all.exons <- new('ExomeDepth', test = my.test, reference = my.reference.selected, formula = 'cbind(test, reference) ~ 1')

-all.exons <- CallCNVs(x = all.exons, transition.probability = 10^-4, chromosome = ExomeCount.dafr$chromosome, start = ExomeCount.dafr$start, end = ExomeCount.dafr$end, name = ExomeCount.dafr$exon)

Regarding this I have some questions:

1. How is it possible that CallCNV() function call just 1 CNV always for the same chromosome? In my case, whatever test sample I use as input it always result with: "Number of calls for chromosome 16 : 1". Am I doing something wrong?
2. I examinated the ExomeCount data.frame that you provide as data example in the vignettes. Why when I run the getBamCounts() function the resulting dataframe is different from the one thet you provide? For instance, seqnames, ranges and strand as well as names are not present in my.counts once I run the command with my input bam files. What am I missing to do?
3. When I loop over multiple samples as suggested in the vignettes I use these commands:
   -ExomeCount.mat <- as.matrix(ExomeCount.dafr[, grep(names(ExomeCount.dafr), pattern = 'Exome.*')])
   -nsamples <- ncol(ExomeCount.mat)
   My question is: why the ExomeCount.mat object is empty when I use the same commands?

Thank you in advance for the attention,

I hope you can help me soon

Does ExomeDepth require multiple control samples for each run? Would having only 1 control sample work?

Thanks

Error in CallCNVs(x = all.exons, transition.probability = trans_prob, :
Chromosome, start and end vector must have the same lengths.

my bed file seems to be well tabulated as:

chr1 69090 70008 OR4F5
chr1 450739 451678 OR4F16
chr1 450739 451678 OR4F29
chr1 450739 451678 OR4F3
chr1 685715 686654 OR4F16
chr1 685715 686654 OR4F29
chr1 685715 686654 OR4F3
chr1 925941 926013 SAMD11
chr1 930154 930336 SAMD11
chr1 931038 931089 SAMD11

what's wrong ?

Dear Vincent,

Thank you for developing and maintaining ExomeDepth!

Are there any plans to support GRCh38, more specifically have you crated by any chance
exons.GRCh38 Positions of exons on build GRCh38 of the human genome
exons.GRCh38 Positions of exons on build GRCh38 of the human genome and on chromosome X
genes.GRCh38 Positions of genes on build GRCh38 of the human genome

Best,
Mike

Hello, I'm a new to ExomeDepth and is currently try it on my custom targeted capture panel. At the end of getBamCounts, the following error shows up.

Now parsing ~/OPERATION/RedCellNGS_hg19/BAM_testing/Normal9.bam
Parsing chromosome chr1
Parsing chromosome chr10
Parsing chromosome chr11
Parsing chromosome chr12
Parsing chromosome chr13
Parsing chromosome chr14
Parsing chromosome chr15
Parsing chromosome chr16
Parsing chromosome chr17
Parsing chromosome chr18
Parsing chromosome chr19
Parsing chromosome chr2
Parsing chromosome chr20
Parsing chromosome chr21
Parsing chromosome chr22
Parsing chromosome chr3
Parsing chromosome chr4
Parsing chromosome chr5
Parsing chromosome chr6
Parsing chromosome chr7
Parsing chromosome chr8
Parsing chromosome chr9
Number of counted fragments : 836974
There were 50 or more warnings (use warnings() to see the first 50)

warnings()
Warning messages:
1: In .Seqinfo.mergexy(x, y) :
Each of the 2 combined objects has sequence levels not in the other:

• in 'x': 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22
• in 'y': chrX, chrY, chr6_ssto_hap7, chr6_mcf_hap5, chr6_cox_hap2, chr6_mann_hap4, chr6_apd_hap1, chr6_qbl_hap6, chr6_dbb_hap3, chr17_ctg5_hap1, chr4_ctg9_hap1, chr1_gl000192_random, chrUn_gl000225, chr4_gl000194_random, chr4_gl000193_random, chr9_gl000200_random, chrUn_gl000222, chrUn_gl000212, chr7_gl000195_random, chrUn_gl000223, chrUn_gl000224, chrUn_gl000219, chr17_gl000205_random, chrUn_gl000215, chrUn_gl000216, chrUn_gl000217, chr9_gl000199_random, chrUn_gl000211, chrUn_gl000213, chrUn_gl000220, chrUn_gl000218, chr19_gl000209_random, chrUn_gl000221, chrUn_gl000214, chrUn_gl000228, chrUn_gl000227, chr1_gl000191_random, chr19_gl000208_random, chr9_gl000198_random, chr17_gl000204_random, chrUn_gl000233, chrUn_gl000237, chrUn_gl000230, chrUn_gl000242, chrUn_gl000243, chrUn_gl000241, chrUn_gl000236, chrUn_gl000240, chr17_gl0 [... truncated]
  2: In .Seqinfo.mergexy(x, y) :
  Each of the 2 combined objects has sequence levels not in the other:
• in 'x': 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22
• in 'y': chrX, chrY, chr6_ssto_hap7, chr6_mcf_hap5, chr6_cox_hap2, chr6_mann_hap4, chr6_apd_hap1, chr6_qbl_hap6, chr6_dbb_hap3, chr17_ctg5_hap1, chr4_ctg9_hap1, chr1_gl000192_random, chrUn_gl000225, chr4_gl000194_random, chr4_gl000193_random, chr9_gl000200_random, chrUn_gl000222, chrUn_gl000212, chr7_gl000195_random, chrUn_gl000223, chrUn_gl000224, chrUn_gl000219, chr17_gl000205_random, chrUn_gl000215, chrUn_gl000216, chrUn_gl000217, chr9_gl000199_random, chrUn_gl000211, chrUn_gl000213, chrUn_gl000220, chrUn_gl000218, chr19_gl000209_random, chrUn_gl000221, chrUn_gl000214, chrUn_gl000228, chrUn_gl000227, chr1_gl000191_random, chr19_gl000208_random, chr9_gl000198_random, chr17_gl000204_random, chrUn_gl000233, chrUn_gl000237, chrUn_gl000230, chrUn_gl000242, chrUn_gl000243, chrUn_gl000241, chrUn_gl000236, chrUn_gl000240, chr17_gl0 [... truncated]
  3: In .Seqinfo.mergexy(x, y) :
  Each of the 2 combined objects has sequence levels not in the other:
• in 'x': 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22
• in 'y': chrX, chrY, chr6_ssto_hap7, chr6_mcf_hap5, chr6_cox_hap2, chr6_mann_hap4, chr6_apd_hap1, chr6_qbl_hap6, chr6_dbb_hap3, chr17_ctg5_hap1, chr4_ctg9_hap1, chr1_gl000192_random, chrUn_gl000225, chr4_gl000194_random, chr4_gl000193_random, chr9_gl000200_random, chrUn_gl000222, chrUn_gl000212, chr7_gl000195_random, chrUn_gl000223, chrUn_gl000224, chrUn_gl000219, chr17_gl000205_random, chrUn_gl000215, chrUn_gl000216, chrUn_gl000217, chr9_gl000199_random, chrUn_gl000211, chrUn_gl000213, chrUn_gl000220, chrUn_gl000218, chr19_gl000209_random, chrUn_gl000221, chrUn_gl000214, chrUn_gl000228, chrUn_gl000227, chr1_gl000191_random, chr19_gl000208_random, chr9_gl000198_random, chr17_gl000204_random, chrUn_gl000233, chrUn_gl000237, chrUn_gl000230, chrUn_gl000242, chrUn_gl000243, chrUn_gl000241, chrUn_gl000236, chrUn_gl000240, chr17_gl0 [... truncated]
  4: In .Seqinfo.mergexy(x, y) :

Furthermore, the output of getBamCounts in my case is a list and when I tried to convert it to a data frame, it failed with the following error.

ExomeCount.dafr <- as(ExomeCount, 'data.frame')
Error in as(ExomeCount, "data.frame") :
internal problem in as(): “tbl_df” is(object, "data.frame") is TRUE, but the metadata asserts that the 'is' relation is FALSE
typeof(ExomeCount)
[1] "list"

My script is as follows.
data(exons.hg19)
hg19<-file.path("/OPERATION/RedCellNGS_hg19/ref","hg19.fa")
bamName <- list.files("/OPERATION/RedCellNGS_hg19/BAM", pattern = '*.bam$')
bamFile <- file.path("~/OPERATION/RedCellNGS_hg19/BAM",bamName)
ExomeCount <- getBamCounts(bed.frame = exons.hg19,bam.files = bamFile,include.chr = TRUE,referenceFasta = hg19)

My sessionInfo() is as follows.
R version 3.6.2 (2019-12-12)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 18.04.3 LTS

Matrix products: default
BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.7.1
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.7.1

locale:
[1] LC_CTYPE=C.UTF-8 LC_NUMERIC=C LC_TIME=C.UTF-8
[4] LC_COLLATE=C.UTF-8 LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
[7] LC_PAPER=C.UTF-8 LC_NAME=C LC_ADDRESS=C
[10] LC_TELEPHONE=C LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] stats graphics grDevices utils datasets methods base

other attached packages:
[1] ExomeDepth_1.1.13

loaded via a namespace (and not attached):
[1] matrixStats_0.55.0 lattice_0.20-38
[3] IRanges_2.20.1 Rsamtools_2.2.1
[5] Biostrings_2.54.0 GenomicAlignments_1.22.1
[7] bitops_1.0-6 grid_3.6.2
[9] GenomeInfoDb_1.22.0 stats4_3.6.2
[11] magrittr_1.5 zlibbioc_1.32.0
[13] XVector_0.26.0 S4Vectors_0.24.1
[15] Matrix_1.2-18 aod_1.3.1
[17] BiocParallel_1.20.1 tools_3.6.2
[19] Biobase_2.46.0 RCurl_1.95-4.12
[21] DelayedArray_0.12.1 parallel_3.6.2
[23] compiler_3.6.2 BiocGenerics_0.32.0
[25] GenomicRanges_1.38.0 SummarizedExperiment_1.16.1
[27] GenomeInfoDbData_1.2.2

Thank you very much for your help!

Thanks a lot for making such a powerful tool.
when I use it with 20 normal samples,
one MET 10 copy standard sample(4 different input hybird library)
, but it continue throw the error

Error in [<-.data.frame(*tmp*, , column.name, value = c(NA_character_, :
replacement has 2 rows, data has 0

after I continuously remove the error column, it finally went through, and the error column starts from the MET 10 copy sample, and it just kept one MET 10 copy sample, and I can not know the copy number from the result

aftrer remove 20th

Hi,

I have used ExomeDepth for detecting CNV on WES, and I found there were several CNVs are different type of CNV but with same start coordinate, like the example below:
chr10 46109621 46115329 gain N8037M
chr10 46109621 46123270 loss N8037M
chr10 46109621 46131358 gain N8037M

I changed the format to bed format. How can I know this area is a gain or loss?

Many thanks,
Min

Dear Vincent,

Thanks for v1.1.6 update. I notice the vignette pdf is messed up in section 9. The content in p10 should go under section 9 in p9, right?

Best wishes,
Fengyuan

i try to install DEcON,but it cant install ExomeDepth

Retrieving 'https://cran.ma.imperial.ac.uk/src/contrib/Archive/ExomeDepth/ExomeDepth_1.1.15.tar.gz' ...
Warning: failed to find source for 'ExomeDepth 1.1.15' in package repositories
Warning: failed to retrieve 'https://cran.ma.imperial.ac.uk/src/contrib/Archive/ExomeDepth/ExomeDepth_1.1.15.tar.gz' [error code 22]
Error: failed to retrieve package '[email protected]'
In addition: There were 11 warnings (use warnings() to see them)
Traceback (most recent calls last):

try to rearch ExomeDepth on cran ,I get that

https://cran.r-project.org/web/packages/ExomeDepth/index.html

When looping through a large number of samples, I found that ExomeDepth would induce a crash on very specific files. This manifests as a termination of R with no error message. I have dug down into the function code and have attached data and code to reproduce the error.

The error only happens with the specific order of inputs (i.e. column 2 in "fail" as the test and column 1 as the reference). When flipping the two, the code executes without any issue. I have attempted to trace the error deeper but when running the ExomeDepth constructor piecemeal, a crash cannot be induced.

It is perplexing that specific data would cause such a crash. There is nothing unique about the sample associated with the test column even on manual inspection of the data.

library(ExomeDepth)

fail = read.delim("failure_file.txt", sep = "\t", stringsAsFactors = F)

formula = "cbind(test, reference) ~ 1"
phi.bins = 1
data = NULL


my.mod <- new("ExomeDepth", test = fail[,2], 
            reference = fail[,1], formula = formula, data = data, 
            phi.bins = phi.bins, verbose = TRUE)

failure_file.txt

Hello,
I have been able to produce very promising CNV calls with ExomeDepth, but now I would like to use the results to check the consistency in given trios (parents and child). Some CNV calls might show small differences due to differences in the hidden states and consequently lack a proper comparability.
Q: (How) can I use ExomeDepth to give it a list of possible CNVs (from previous CNV calls) and then receive this list aligned with the probabilities (or any other score) of each CNV to be present in a given sample?
That way, I could make the CNV calls of two different but related samples comparable and thereby check mendelian consistency within a trio.

Thank you very much and best regards,
Vinzenz

Hey,
I use bwa-mem2 to get some bams , for one bam, getBamCounts returns 0 , but for other bams , the getBamCounts result is correct.
so whats wrong with this function ?

Hi,

Is it possible to get reads ratios for regions that are not called as CNVs by the function CallCNVs() ? Or can we get them from the getBamCounts() or some other way?

Thanks very much,
Rajat

Hi
I just noticed an issue where regions with more reads than expected are being called as deletions (reads.ratio >> 1). In the example there's also an issue with overlapping regions being called as both deletions and duplications. Anyone have any ideas on what might be happening here?

124085,124093,"deletion",9,55142239,55156999,"chr7","chrchr7:55142239-55156999",26.5,4474,17460,3.9
124085,124113,"duplication",20,55142239,55201432,"chr7","chrchr7:55142239-55201432",57.8,12671,45280,3.57
124085,124126,"deletion",13,55142239,55366179,"chr7","chrchr7:55142239-55366179",115,20719,81052,3.91
124085,124134,"duplication",8,55142239,55431346,"chr7","chrchr7:55142239-55431346",85.4,23742,90328,3.8
124085,124136,"deletion",2,55142239,55433799,"chr7","chrchr7:55142239-55433799",132,25303,97612,3.86
124085,124187,"duplication",51,55142239,55795890,"chr7","chrchr7:55142239-55795890",178,45707,146585,3.21

Best regards
-Harald

Hi @vplagnol,

Is exomedepth compatible with R verson >4 in ubuntu 22?
I have tried to run it, but it requested that I reinstall it, then failed:

g++ -std=gnu++14 -I"/usr/share/R/include" -DNDEBUG -fpic -g -O2 -ffile-prefix-map=/build/r-base-4A2Reg/r-base-4.1.2=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g -c CNV_estimate.cpp -o CNV_estimate.o
gcc -I"/usr/share/R/include" -DNDEBUG -fpic -g -O2 -ffile-prefix-map=/build/r-base-4A2Reg/r-base-4.1.2=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g -c ExomeDepth_init.c -o ExomeDepth_init.o
gcc -I"/usr/share/R/include" -DNDEBUG -fpic -g -O2 -ffile-prefix-map=/build/r-base-4A2Reg/r-base-4.1.2=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g -c VP_gamma.c -o VP_gamma.o
In file included from gsl_math.h:22,
from VP_gamma.c:23:
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:85:1: error: expected ‘;’ before ‘extern’
85 | __MATHCALL_VEC (acosh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:87:1: error: expected ‘;’ before ‘extern’
87 | __MATHCALL_VEC (asinh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:89:1: error: expected ‘;’ before ‘extern’
89 | __MATHCALL_VEC (atanh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:119:1: error: expected ‘;’ before ‘extern’
119 | __MATHCALL_VEC (expm1,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:122:1: error: expected ‘;’ before ‘extern’
122 | __MATHCALL_VEC (log1p,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:147:1: error: expected ‘;’ before ‘extern’
147 | __MATHCALL_VEC (hypot,, (Mdouble __x, Mdouble __y));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:85:1: error: expected ‘;’ before ‘extern’
85 | __MATHCALL_VEC (acosh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:87:1: error: expected ‘;’ before ‘extern’
87 | __MATHCALL_VEC (asinh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:89:1: error: expected ‘;’ before ‘extern’
89 | __MATHCALL_VEC (atanh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:119:1: error: expected ‘;’ before ‘extern’
119 | __MATHCALL_VEC (expm1,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:122:1: error: expected ‘;’ before ‘extern’
122 | __MATHCALL_VEC (log1p,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:147:1: error: expected ‘;’ before ‘extern’
147 | __MATHCALL_VEC (hypot,, (Mdouble __x, Mdouble __y));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:85:1: error: expected ‘;’ before ‘extern’
85 | __MATHCALL_VEC (acosh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:87:1: error: expected ‘;’ before ‘extern’
87 | __MATHCALL_VEC (asinh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:89:1: error: expected ‘;’ before ‘extern’
89 | __MATHCALL_VEC (atanh,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:119:1: error: expected ‘;’ before ‘extern’
119 | __MATHCALL_VEC (expm1,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:122:1: error: expected ‘;’ before ‘extern’
122 | __MATHCALL_VEC (log1p,, (Mdouble __x));
| ^~~~~~~~~~~~~~
/usr/include/x86_64-linux-gnu/bits/mathcalls.h:147:1: error: expected ‘;’ before ‘extern’
147 | __MATHCALL_VEC (hypot,, (Mdouble __x, Mdouble __y));
| ^~~~~~~~~~~~~~
make: *** [/usr/lib/R/etc/Makeconf:168: VP_gamma.o] Error 1
ERROR: compilation failed for package ‘ExomeDepth’

`

After successful generation of several plots, an error occurs
"Error in if (exon.array$short.name[nrow(exon.array)] != exon.array$short.name[nrow(exon.array) - :
argument is of length zero"
What could be the reason and how can I fix it?
I have version 1.1.10 from CRAN.

Hi,
I'm trying to use ExomeDepth but I have the following error message:

my.counts <- getBamCounts(bed.frame = exons.hg19,bam.files= sample1.bam, include.chr=TRUE)
Parse 1 BAM files
[1] "sample1.bam"
Now parsing sample1.bam
[bam_header_read] EOF marker is absent. The input is probably truncated.
Parsing chromosome chr1
[bam_header_read] EOF marker is absent. The input is probably truncated.
[bam_index_load] fail to load BAM index.
Error in value[3L] :
failed to open BamFile: failed to load BAM index
file:sample1.bam
Inoltre: Warning messages:
1: 'isNotPrimaryRead' is deprecated.
Use 'isSecondaryAlignment' instead.
See help("Deprecated")
2: 'isNotPrimaryRead' is deprecated.
Use 'isSecondaryAlignment' instead.
See help("Deprecated")

help("Deprecated")
Can you help me?
Thanks

where it says
name = ExomeCount.dafr$names
it should say
name = ExomeCount.dafr$exon
in here

also, is the transition.probability = 10^-4 value still recommended? other users are working with higher values (such as 10^-3, or even DeCon's 10^-2). any thoughts on this?

Hi vplgnol
Thanks for this wonderful tool. we are trying to setup CNV pipeline for germline samples, but we are facing some issue mentioned below. Can you please look into the issue and let us know the solution. This would really help us to accelerate in setting up our pipeline.

1. R version is 3.6.3
2. samtools = Version: 1.9 (using htslib 1.9)
3. bwa = Version: 0.7.12-r1039

thanks alot

my.counts <- getBamCounts(bed.frame = exons.hg19, bam.files = "bams.list", include.chr = TRUE)
Parse 1 BAM files
[1] "bams.list"
Now parsing bams.list
[samopen] no @sq lines in the header.
Error in value[3L] :
failed to open BamFile: SAM/BAM header missing or empty
file: 'bams.list'