abyzovlab / CNVnator

a tool for CNV discovery and genotyping from depth-of-coverage by mapped reads
Other
209 stars 66 forks source link

Complete usage example #155

Closed dfermin closed 3 years ago

dfermin commented 5 years ago

Hi

I need some help using CNVnator to identify DEL events. For testing I have 3 BAM files aligned to GRCh37 and I want to look for DEL events in chromosome 22. Ultimately I'm going to run it on 600 samples. But first I want to know if it's working and how to interpret the output.

These are the commands I'm running, cnvnator is in my path:

export HS37d5=/data/fastas/hs37d5

cnvnator -root chr22.root -genome GRCh37 -chrom 22 -tree \
    /nfs/wgs30x/bam.files/samp1.bam \
    /nfs/wgs30x/bam.files/samp2.bam \
    /nfs/wgs30x/bam.files/samp3.bam;

cnvnator -root chr22.root -genome GRCh37 -chrom 22 -his 1000 -d $HS37D5

cnvnator -root chr22.root -genome GRCh37 -chrom 22 -stat 1000

cnvnator -root chr22.root -genome GRCh37 -chrom 22 -partition 1000

cnvnator -root chr22.root -genome GRCh37 -chrom 22 -call 1000 > chr22.cnvnator.out

cnvnator -pe /nfs/wgs30x/bam.files/samp1.bam \
   /nfs/wgs30x/bam.files/samp2.bam \
   /nfs/wgs30x/bam.files/samp3.bam \
   -qual 20 -over 0.8 -root chr22.root -f new_chr22.support_reads;

cnvnator2VCF.pl -reference GRCh37 chr22.cnvnator.out $HS37D5 > final_out.vcf

Is this the correct pipeline? I'm sorry if this is a stupid question but I can't find a complete example. The github site provides an excellent information on executing the commands but doesn't tell me much about parameters or best practices.

For instance in the cnvnator -pe command what does the -qual argument do? What about the -over option? The -f command produces the final output? What is the optimal bin_size to use? Is there a way to compute it?

I'd be happy to share my completed pipeline example once I get it working.

Thanks in advance for any and all help.

abyzov commented 5 years ago

Hi, you have to do analysis per sample, e.g.,

cnvnator -root samp1.chr22.root -chrom 22 -tree /nfs/wgs30x/bam.files/samp1.bam

cnvnator -root chr22.root -chrom 22 -his 1000 -d $HS37D5

cnvnator -root chr22.root -chrom 22 -stat 1000

cnvnator -root chr22.root -chrom 22 -partition 1000

cnvnator -root chr22.root -chrom 22 -call 1000 > samp1.chr22.cnvnator.out

cnvnator -pe /nfs/wgs30x/bam.files/samp1.bam \ /nfs/wgs30x/bam.files/samp2.bam \ /nfs/wgs30x/bam.files/samp3.bam \ -qual 20 -over 0.8 -root chr22.root -f samp1.chr22.cnvnator.out

-qual option is the cut off for mapping quality for paired-reds supporting CNVs -over is the reciprocal fraction of overlapping between the span of paired read and CNV boundaries

We recommend to set up bin size such that average RD is about 4-5 time of its standard deviation. You can use option -eval to evaluate that.

Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org/ tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

dfermin commented 5 years ago

Thanks! This was very helpful.

I'm still not sure if this is working. The last command is giving me an error:

cnvnator -pe /nfs/wgs30x/bam.files/samp1.bam \
   /nfs/wgs30x/bam.files/samp2.bam \
   /nfs/wgs30x/bam.files/samp3.bam \
   -qual 20 -over 0.8 -root chr22.root -f samp1.chr22.cnvnator.out

Error in <TString::AssertElement>: out of bounds: i = -1, Length = 0
Invalid start .
/nfs/turbo/goldfinch/apps/CNVnator/CNVnator/cnvnator -root chr22.root -chrom 22 -call 150       0
Nowhere to scroll.

Not sure what it means. No file containing read counts is produced though.

abyzov commented 5 years ago

What is the content of the samp1.chr22.cnvnator.out file? Could you share a few head lines.

Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

dfermin commented 5 years ago

Here are the first few lines of samp1.chr22.cnvnator.out

deletion        22:1-16050150   1.60502e+07     3.76999e-06     9.92966e-15     0       9.93096e-15     0       -1
duplication     22:16050151-16064700    14550   1.70541 0       995056  0       3.14236e+06     0.126538
deletion        22:16082551-16090650    8100    0.668045        9.98931e-08     5.79435e+06     0.000157416     584567  0.719741
deletion        22:16091701-16096200    4500    0.502126        7.32228e-07     1.26348e-13     0.475064        0.0034244       0.958386
deletion        22:16106551-16113750    7200    0.496534        2.21351e-11     1.62003e-08     3.12495e-11     1.08537e-16     0.873777
duplication     22:16125001-16128450    3450    1.49688 0.0424363       1.19327e+09     21647.8 2.03625e+09     0.656954
duplication     22:16160251-16167900    7650    1.84143 4.66659e-09     2.8248e+08      0.000734881     5.33866e+08     0.445384
duplication     22:16173301-16186800    13500   1.50596 0       2.84163e+09     2.8016e-08      2.84618e+09     0.57962
duplication     22:16189651-16207200    17550   1.62141 0       2.86322e+09     0       2.86415e+09     0.352988
duplication     22:16241401-16250850    9450    1.42914 0.00624862      2.80371e+09     0.843094        2.80944e+09     0.448143
abyzov commented 5 years ago

Hi, not sure why you are getting the last error message but may be you have an extra header line in the file? I’ve fixed error with string being out of bound. The updated version is committed into the repository. I tried it with the output you sent to me and it works with no error on my system. BTW, note that -pe support will be show support of all discordant reads for all input bam files. If you want output per bam then you have to run the script for each bam file.

Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

hyBio commented 5 years ago

Hi, Can cnvnator software become more convenient to use? Although this software has a good reputation in call CNV, there are many people feel it difficult to use. I also encountered the same problem as dfermin, but I don't know where the problem is. Does this error affect the final result?

abyzov commented 5 years ago

Hi, we are working on making CNVnator better. Not sure which error you are referring to.

Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

hyBio commented 5 years ago

This is my error: "Error in : out of bounds: i = -1, Length = 0Error in : out of bounds: i = -1, Length = 0 Error in : out of bounds: i = -1, Length = 0 Error in : out of bounds: i = -1, Length = 0 Error in : out of bounds: i = -1, Length = 0" And this is my code: `~/speedseq/bin/annotate_rd.py \ ~/speedseq/bin/cnvnator \ -s sample1 \ -w 100 \ -r path_to_temp/cnvnator-temp/sample1.bam.hist.root \ -v sample1.sv.vcf \

sample1.sv.rd.vcf`

------------------ 原始邮件 ------------------ 发件人: "abyzov"notifications@github.com; 发送时间: 2019年7月3日(星期三) 下午5:16 收件人: "abyzovlab/CNVnator"CNVnator@noreply.github.com; 抄送: "HU YAN"568038247@qq.com;"Comment"comment@noreply.github.com; 主题: Re: [abyzovlab/CNVnator] Complete usage example (#155)

Hi, we are working on making CNVnator better. Not sure which error you are referring to.

Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

— You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.

hyBio commented 5 years ago

Hi, we are working on making CNVnator better. Not sure which error you are referring to. Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic ----------------------------- Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.org<http://www.abyzovlab.org> tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

This is my error: "Error in : out of bounds: i = -1, Length = 0 Error in : out of bounds: i = -1, Length = 0 Error in : out of bounds: i = -1, Length = 0 Error in : out of bounds: i = -1, Length = 0 Error in : out of bounds: i = -1, Length = 0" And this is my code: `~/speedseq/bin/annotate_rd.py \ ~/speedseq/bin/cnvnator \ -s sample1 \ -w 100 \ -r path_to_temp/cnvnator-temp/sample1.bam.hist.root \ -v sample1.sv.vcf \

sample1.sv.rd.vcf` Thanks very much!

abyzov commented 5 years ago

Hi, this error does not affect results, and it is also fixed in the github repository.

Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

hyBio commented 5 years ago

Hi, this error does not affect results, and it is also fixed in the github repository. Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic ----------------------------- Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.org<http://www.abyzovlab.org> tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

Ok,Thank you so much.

dfermin commented 5 years ago

Thanks for the code update to github. Based on your suggested changes this is my now working pipeline for running CNVnator. I'm only showing the commands for 1 BAM file sampN. But I would run the code on multiple samples 1-N. Is this correct? Ultimately I'd like to get a VCF file of the DEL/DUP calls for each sample. Thanks

cnvnator -root sampN.chr22.root -chrom 22 -tree bam_files/sampN.bam

export HS37D5=/data/fastas/hs37d5
cnvnator -root sampN.chr22.root -chrom 22 -his 1000 -d $HS37D5

cnvnator -root sampN.chr22.root -chrom 22 -stat 1000

cnvnator -root sampN.chr22.root -chrom 22 -partition 1000

cnvnator -root sampN.chr22.root -chrom 22 -call 1000 > sampN.chr22.cnvnator.out

## produce VCF file
cnvnator2VCF.pl -reference GRCh37 sampN.chr22.cnvnator.out $HS37D5 > sampN.chr22.cnvnator.vcf
abyzov commented 5 years ago

Hi, yes you have to run it for each bam. We currently don’t have a standard way to merge multiple callsets, as there is no consensus about this in the field.

Alexej Abyzov, Ph.D. Senior Associate Consultant, Assistant Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

Solyris83 commented 4 years ago

Hi Alexej, I have this error instead, any help is appreciated.

cnvnator -pe ../KO.recal.bam -qual 20 -over 0.8 -root KO.root -f KO.CNVnator.txt

Break segmentation violation Generating stack trace... 0x000055a7c7839a00 in from ../../../software/CNVnator/cnvnator 0x000055a7c782d953 in from ../../../software/CNVnator/cnvnator 0x000055a7c77f9f5e in from ../../../software/CNVnator/cnvnator 0x000055a7c77fb06b in from ../../../software/CNVnator/cnvnator 0x000055a7c77c5412 in from ../../../software/CNVnator/cnvnator 0x00007fc7d859ab97 in __libc_start_main + 0xe7 from /lib/x86_64-linux-gnu/libc.so.6 0x000055a7c77c6eea in from ../../../software/CNVnator/cnvnator

abyzov commented 4 years ago

Hi, can’t really say what is the reason for crash. Does it work when you run this function interactively?

Alexej Abyzov, Ph.D. Senior Associate Consultant, Associate Professor of Biomedical Informatics, Department of Health Sciences Research, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

abyzov commented 3 years ago

Hi, most likely there is some problem with the input file samp1.chr22.cnvnator.out. Is it in the correct format?

Alexej Abyzov, Ph.D. Senior Associate Consultant, Associate Professor of Biomedical Informatics, Department of Quantitative Health Sciences, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978 fax: +1-(507)-284-0745

Saeideh-Ashouri commented 3 years ago

Hi Alexej, I have this error instead, any help is appreciated.

cnvnator -pe ../KO.recal.bam -qual 20 -over 0.8 -root KO.root -f KO.CNVnator.txt

Break segmentation violation Generating stack trace... 0x000055a7c7839a00 in from ../../../software/CNVnator/cnvnator 0x000055a7c782d953 in from ../../../software/CNVnator/cnvnator 0x000055a7c77f9f5e in from ../../../software/CNVnator/cnvnator 0x000055a7c77fb06b in from ../../../software/CNVnator/cnvnator 0x000055a7c77c5412 in from ../../../software/CNVnator/cnvnator 0x00007fc7d859ab97 in __libc_start_main + 0xe7 from /lib/x86_64-linux-gnu/libc.so.6 0x000055a7c77c6eea in from ../../../software/CNVnator/cnvnator

Hello @Solyris83

I know almost one year has passed since you wrote this issue; however, as I have the same problem with -pe command (reported it as issue #245 ), would you please tell me how you solved this issue? And another question I have is that for -pe command, don't we need to redirect the results to an output file?

I hope you can reply my questions.

Thank you in advance, Best, Saeideh

abyzov commented 3 years ago

Hi, sorry for later response. Do you have the same issue if you run the command interactively?

Alexej Abyzov, Ph.D. Senior Associate Consultant, Associate Professor of Biomedical Informatics, Department of Quantitative Health Sciences, Center for Individualized Medicine, Mayo Clinic

Mayo Clinic, 200 1st street SW, Harwick 3-12 Rochester, MN 55905 www.abyzovlab.orghttp://www.abyzovlab.org tel: +1-(507)-538-0978