Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

octopus appears to hang #4

Closed
jaredo opened this issue Oct 14, 2016 · 5 comments
Closed

octopus appears to hang #4

jaredo opened this issue Oct 14, 2016 · 5 comments
Assignees

Comments

@jaredo
Copy link
Contributor

jaredo commented Oct 14, 2016

Hello! I have been trying to run octopus (Ubuntu 16.04) on a small test example. I sampled ~100kb either side of an interesting gene on each autosome from NA12878 to construct a small test bam. The tool appears to start calling variants but then stop, CPU usage remains at 100% (and there is plenty of memory available). Also ran into the same issue with 30X WGS bams.

The mini-bam (hg19) is available here: https://dl.dropboxusercontent.com/u/15887058/NA12878.bam

$ software/octopus/bin/octopus  -I NA12878.bam -R genome.fa -o test.vcf.gz
[2016-10-14 09:38:37] <INFO> ------------------------------------------------------------------------
[2016-10-14 09:38:37] <INFO> octopus v0.1 alpha
[2016-10-14 09:38:37] <INFO> Copyright (c) 2016 University of Oxford
[2016-10-14 09:38:37] <INFO> ------------------------------------------------------------------------
Not a BGZF file: test.vcf.gz
[2016-10-14 09:38:37] <INFO> Done initialising calling components in 184ms
[2016-10-14 09:38:37] <INFO> Invoked calling model: individual
[2016-10-14 09:38:37] <INFO> Detected 1 sample: "NA12878"
[2016-10-14 09:38:37] <INFO> Writing calls to "test.vcf.gz"
[2016-10-14 09:38:37] <INFO> ------------------------------------------------------------------------
[2016-10-14 09:38:37] <INFO>      current      |                   |     time      |     estimated
[2016-10-14 09:38:37] <INFO>      position     |     completed     |     taken     |     ttc
[2016-10-14 09:38:37] <INFO> ------------------------------------------------------------------------
[2016-10-14 09:38:42] <INFO>  chr1:197001111             6.4%              5s           1.33m

it has been stuck in this state for ~6 hours.

@dancooke dancooke self-assigned this Oct 14, 2016
@jaredo
Copy link
Contributor Author

jaredo commented Oct 14, 2016

PS. I just tried on octopus-0.1.1-alpha and it appears to still be there.

@dancooke dancooke added the bug label Oct 14, 2016
@dancooke
Copy link
Member

dancooke commented Oct 14, 2016

Hi Jared!

Many thanks for flagging this. I've found the problem! It's because the reference genome hg19 contains a mixture of lower & upper case bases while the reads are upper case only. Octopus was using the raw input, and considers mixed case versions of the same base as being distinct, causing high quality candidate variants at almost every position!

I've added sequence capitalisation transforms to both reference and read sequences. They can be turned off with the command line options --disable-reference-base-capitalisation and --disable-read-base-capitalisation respectively. The patch is available in release 0.1.2-alpha. The sample you posted now runs in about 5 minutes on my machine and in less than 2 minutes with --fast enabled. I'd also recommend taking a look at the command --threads if you have multiple cores available.

I'm also a bit concerned by the spurious message Not a BGZF file: test.vcf.gz you're getting. Does this appear every time you try to run octopus? I'd be grateful if you'd consider making a new issue detailing this.

Best wishes,

Dan

@gerton
Copy link
Collaborator

gerton commented Oct 14, 2016

Hi Dan,

Good find! However, I would consider removing the new options —disable-reference-base-capitalisation and likewise for reads. I see no situation where you want to consider upper vs. lower case to be non-matching.

How do you treat Ns in the reads? In Stampy I set the quality to 0 to effectively make N match anything. (And for the reference I set it to G. One has to do something…)

Stuck in Copenhagen — flight delayed. Having fun implementing a bit more of the BWT algorithm (in the company of a Carlsberg..!)

BW
Gerton

On 14 Oct 2016, at 20:03, Daniel Cooke <notifications@github.com mailto:notifications@github.com> wrote:

Hi!

Many thanks for flagging this. I've found the problem! It's because the reference genome hg19 contains a mixture of lower & upper case bases while the reads are upper case only. Octopus was using the raw input, and considers mixed case versions of the same base as being distinct, causing high quality candidate variants at almost every position!

I've added sequence capitalisation transforms to both reference and read sequences. They can be turned off with the command line options --disable-reference-base-capitalisation and --disable-read-base-capitalisation respectively. The patch is available in release 0.1.2-alpha. The sample you posted now runs in about 5 minutes on my machine and in less than 2 minutes with --fast enabled. I'd also recommend taking a look at the command --threads if you have multiple cores available.

I'm also a bit concerned by the spurious message Not a BGZF file: test.vcf.gz you're getting. Does this appear every time you try to run octopus? I'd be grateful if you'd consider making a new issue detailing this.

Best wishes,

Dan


You are receiving this because you are subscribed to this thread.
Reply to this email directly, view it on GitHub #4 (comment), or mute the thread https://github.com/notifications/unsubscribe-auth/AD3zPjHX4RqDs5MVnzncO0whzpBeEtbHks5qz8QNgaJpZM4KXB3U.

@dancooke
Copy link
Member

dancooke commented Oct 14, 2016

Hi Gerton,

Yes I think you're right, I actually just read the mixed-case is only intended to mark low complexity regions. I'll remove them from the next release.

The raw cigar generator won't propose candidate base substitutions where either the reference or read base is N. Indels are allowed to contain Ns however. The local reassembler ignores any read k-mer with an N, creating a break in the k-mer path, but doesn't throw away the entire read. Reference sequence containing Ns is not allowed by the assembler as it requires a single contiguous path. Unfortunately this means the entire active region (~1000bp) won't be assembled.

I'll expect a flurry of commits!

Best,
Dan

@jaredo
Copy link
Contributor Author

jaredo commented Oct 15, 2016

Hi guys,

Thanks for the fast response. I can confirm this is resolved (takes 3 minutes on my AMI).

I am still getting this Not a BGZF file and will raise another issue once/if I can confirm it is not something weird to do with my environment.

Running on some full WGS data now - very excited to see the results!

@jaredo jaredo closed this as completed Oct 15, 2016
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

3 participants