luntergroup / octopus

Bayesian haplotype-based mutation calling
MIT License
302 stars 38 forks source link

`std::out_of_range` during call set refinement (after 3 days of runtime 😖) #114

Closed DaGaMs closed 4 years ago

DaGaMs commented 4 years ago

Describe the bug I ran Octopus (0.6.3-beta) on some whole-genome sequencing cancer data from TCGA as a test. I re-aligned the data to the GATK-provided hg38 full reference, including all alt alleles and unassembled contigs. The first step finished ok:

[2020-03-18 05:26:44] <INFO> chr16_KI270728v1_random:1090052         96.0%           2d 11h            2h 30m
[2020-03-18 05:32:08] <INFO>               chr22:15369970            96.1%           2d 11h            2h 26m
[2020-03-18 07:00:31] <WARN> Skipping region chr22_KI270733v1_random:138388-138438 as there are too many haplotypes
[2020-03-18 07:30:33] <INFO> chr22_KI270736v1_random:129515          96.2%           2d 13h            2h 27m
[2020-03-18 10:32:38] <INFO>       chrUn_KI270749v1:71896            96.3%           2d 16h            2h 30m
[2020-03-18 11:41:03] <INFO>   chr1_KI270763v1_alt:839577            96.4%           2d 17h            2h 28m
[2020-03-18 11:43:09] <INFO>   chr3_KI270777v1_alt:173649            96.5%           2d 17h            2h 24m
[2020-03-18 11:43:50] <INFO>   chr5_KI270793v1_alt:126136            96.6%           2d 17h            2h 20m
[2020-03-18 11:44:23] <INFO>   chr6_GL000250v2_alt:904697            96.7%           2d 17h            2h 11m
[2020-03-18 11:44:23] <INFO>  chr6_GL000250v2_alt:4326264            96.8%           2d 17h             2h 7m
[2020-03-18 11:44:38] <INFO>   chr7_KI270803v1_alt:196626            96.9%           2d 17h             2h 3m
[2020-03-18 11:46:19] <INFO>   chr8_KI270821v1_alt:985506            97.0%           2d 17h            1h 59m
[2020-03-18 11:47:02] <INFO>  chr10_GL383545v1_alt:179254            97.1%           2d 17h            1h 54m
[2020-03-18 11:47:49] <INFO>  chr12_KI270835v1_alt:238139            97.2%           2d 17h            1h 50m
[2020-03-18 11:48:30] <INFO> chr14_KI270847v1_alt:1416150            97.3%           2d 17h            1h 46m
[2020-03-18 11:49:38] <INFO>  chr15_KI270849v1_alt:244917            97.4%           2d 17h            1h 42m
[2020-03-18 11:50:03] <INFO> chr16_KI270853v1_alt:1967662            97.5%           2d 17h            1h 38m
[2020-03-18 11:50:30] <INFO> chr17_KI270857v1_alt:1428024            97.6%           2d 17h            1h 33m
[2020-03-18 11:50:42] <INFO> chr17_GL000258v2_alt:1457260            97.8%           2d 17h            1h 29m
[2020-03-18 11:51:33] <INFO>  chr19_GL383573v1_alt:385657            97.9%           2d 17h            1h 25m
[2020-03-18 11:52:45] <INFO>  chr22_KI270875v1_alt:101103            98.0%           2d 17h            1h 21m
[2020-03-18 11:53:18] <INFO>   chr4_KI270896v1_alt:147480            98.1%           2d 17h            1h 17m
[2020-03-18 11:53:47] <INFO>  chr6_GL000251v2_alt:1835983            98.2%           2d 17h            1h 13m
[2020-03-18 11:54:07] <INFO>    chr8_KI270901v1_alt:75515            98.3%           2d 17h             1h 9m
[2020-03-18 11:55:16] <INFO> chr15_KI270905v1_alt:1557596            98.4%           2d 17h             1h 4m
[2020-03-18 11:56:18] <INFO> chr15_KI270905v1_alt:4771153            98.5%           2d 17h                1h
[2020-03-18 11:56:29] <INFO>  chr19_GL949747v2_alt:524942            98.6%           2d 17h           56m 46s
[2020-03-18 11:56:47] <INFO>  chr6_GL000252v2_alt:1530335            98.7%           2d 17h           52m 40s
[2020-03-18 11:56:50] <INFO>   chr8_KI270926v1_alt:143032            98.8%           2d 17h           48m 34s
[2020-03-18 11:57:11] <INFO>  chr6_GL000253v2_alt:2369531            98.9%           2d 17h           44m 28s
[2020-03-18 11:57:19] <INFO>  chr6_GL000254v2_alt:2686249            99.0%           2d 17h           40m 23s
[2020-03-18 11:57:32] <INFO>  chr6_GL000253v2_alt:3727571            99.1%           2d 17h           36m 19s
[2020-03-18 11:57:32] <INFO>  chr6_GL000253v2_alt:4677643            99.2%           2d 17h           32m 15s
[2020-03-18 11:57:35] <INFO>  chr19_GL949752v1_alt:987100            99.3%           2d 17h           28m 11s
[2020-03-18 11:57:36] <INFO>  chr6_GL000256v2_alt:2356368            99.4%           2d 17h            24m 8s
[2020-03-18 11:57:40] <INFO>  chr6_GL000256v2_alt:4443203            99.5%           2d 17h            16m 4s
[2020-03-18 14:58:15] <INFO> chrUn_JTFH01000785v1_decoy:306          99.6%           2d 20h           12m 35s
[2020-03-18 17:27:27] <INFO>                            -             100%           2d 23h                 -
[2020-03-18 17:30:10] <INFO> Starting Call Set Refinement (CSR) filtering
[2020-03-18 17:30:27] <INFO> CSR: Starting registration pass
[2020-03-18 17:30:27] <INFO> -------------------------------------------------------------------------------------
[2020-03-18 17:30:27] <INFO>            current             |                   |     time      |     estimated   
[2020-03-18 17:30:27] <INFO>            position            |     completed     |     taken     |     ttc         
[2020-03-18 17:30:27] <INFO> -------------------------------------------------------------------------------------

but then the "Call Set Refinement" dies:

[2020-03-18 21:54:11] <INFO>   chr8_KI270813v1_alt:300230            97.0%           4h 23m            7m 55s
[2020-03-18 21:54:44] <INFO>  chr5_KI270897v1_alt:1144418            98.0%           4h 24m             5m 8s
[2020-03-18 21:54:45] <INFO>  chr6_GL000251v2_alt:3248255            98.1%           4h 24m            4m 52s
[2020-03-18 21:54:52] <INFO>  chr15_KI270906v1_alt:196384            98.2%           4h 24m            4m 35s
[2020-03-18 21:54:52] <INFO> chr15_KI270905v1_alt:3224785            98.3%           4h 24m            4m 19s
[2020-03-18 21:54:55] <INFO>  chr17_KI270908v1_alt:635380            98.4%           4h 24m             4m 2s
[2020-03-18 21:54:58] <INFO>   chr6_GL000252v2_alt:462437            98.5%           4h 24m            3m 46s
[2020-03-18 21:54:58] <INFO>  chr6_GL000252v2_alt:3690848            98.6%           4h 24m            3m 30s
terminate called after throwing an instance of 'std::out_of_range'
  what():  ContigRegion: expanded past contig start

Version

octopus version 0.6.3-beta (HEAD fd981298)
Target: x86_64 Linux 3.10.0-1062.1.2.el7.x86_64
Compiler: GNU 8.3.0
Boost: 1_72

Command Command line to run octopus:

octopus -X 2G -R Data/genomes/hg38_full_gatk_HPV_HBV_HCV_EBV.fa -I tumour.bam control.bam -N control --threads 4 --somatics-only --somatic-forest-file octopus/resources/forests/somatic.v0.6.3-beta.forest | bcftools view -O z -o calls.vcf.gz
dancooke commented 4 years ago

Any way I can access this data to reproduce?

Also, do you really need to be calling all contigs in hg38? Probably a large amount of time is spent processing alt, unplaced, and the decoy contig that most people don't use. It's also where the error occurs. I'd highly recommend restricting calling to the standard chromosomes unless you have plans for the other contigs.

dancooke commented 4 years ago

Also, you don't need to use bcftools to get bgzipped output; Octopus can do this automatically - it recognises the provided file extension to --output. So you can just use:

octopus -X 2G -R Data/genomes/hg38_full_gatk_HPV_HBV_HCV_EBV.fa -I tumour.bam control.bam -N control --threads 4 --somatics-only --somatic-forest-file octopus/resources/forests/somatic.v0.6.3-beta.forest -o calls.vcf.gz
DaGaMs commented 4 years ago

Well, I just ran it on everything because... Why not? Who knows, it might be interesting at some point 😉

Anyway, the data are on our private cluster. What file(s) would you want to have access to? I can move some of it to the BDI cluster.

DaGaMs commented 4 years ago

Also, is there a way to rescue the process, or will I have to start from 0? I have some 364MB of temp.bcf files in octopus-temp...

dancooke commented 4 years ago

Well, I just ran it on everything because... Why not?

Well it could be a third of your runtime is spent in these contigs so that's one reason... I'd at-least skip chrUn_JTFH01000785v1_decoy as all you'll see nothing useful there.

What file(s) would you want to have access to?

I'll just need the BAM files.

Also, is there a way to rescue the process, or will I have to start from 0?

You can merge all these bcfs to get an unfiltered vcf, then you just need to filter this; you can use the --filter-vcf option.

dancooke commented 4 years ago

Closing due to inactivity. Please re-open if needed.