Using NanoCaller with BED files and Indel Calling Gives Empty VCF files about nanocaller HOT 3 OPEN

Hi there,

I ran into the same issue with my ONT data. Does this issue resolve yet?

For my case I could see a few indel vcfs in the intermediate_indel_files but the all files are empty.

Tuan

from nanocaller.

I think I find the reason behind this, cross-ref with issue #1 (already closed)

The default NanoCaller workflow, which is the mode 'both', goes like this:

bam file > SNP calling by NanoCaller > SNP phasing by whatshap > Reads phasing by whatshap > Indel calling by NanoCaller on reads phased in previous step > Combine indel and SNP calls

If you choose 'indel' mode, then we do not have any SNPs that we can use to phase the reads in order to call indels. This option is allowed so that if a user prefers to use a different tool than whatshap to carry out SNPs and reads phasing, or if they already have phased reads that they believe are phased using some other reliable method, then they can bypass the whatshap subroutine by NanoCaller and supply their own phased bam file.

Another reason for this option is that if a user runs NanoCaller in 'both' mode, but the program crashes during indel calling stage for some likely scenarios such as memory exceeding system limits, then they can resume NanoCaller at the indel calling stage, instead of having to do the entire SNP calling and phasing all over again.

We have removed the option for bed files for now. We found for whole genome variant calling that using a bed file would get rid of potential variants near the boundaries of the interval in a bed file. So we found it better to leave it up to a user to filter variants with a bed file using a more sophisticated tool like rtg vcffilter. But I can see having an option like bed file can be useful if you you dont want to spend extra time in calling variants that are outside the intervals in a bed file, and submitting individual calls for each line in a bed file can be cumbersome. We will revisit this issue in our team. Do you have any preferences or thoughts about this?

We provide an example here: https://github.com/WGLab/NanoCaller/tree/master/sample You can see the command we used in the file: https://github.com/WGLab/NanoCaller/blob/master/sample/sample_call

from nanocaller.

Hi,

Regarding --bed, we changed how we used the bed file. Previously we were using it to filter out candidate sites after we generate a set of candidate sites from whole genome. But now it is meant to be used as a way to specify regions in which candidate sites are going to be generated and then we open the BAM file for each bed file entry. So if there are too many lines in the BED file then it will slow down NanoCaller significantly.

Whereas for indel calling failure, it might help to see if there are any phased SNP files created so we can narrow down where the failure happens. It seems like there is a problem with generating phased BAM file from whatshap.

Tuan, can you show me the error you received? Are you able to get phased SNP files? You can get more error details by replacing this line:
https://github.com/WGLab/NanoCaller/blob/master/nanocaller_src/indelCaller.py#L178

by

run_cmd("whatshap haplotag --ignore-read-groups --ignore-linked-read --reference %s %s %s --regions %s:%d-%d --tag-supplementary -o - | samtools view -b -1 --write-index -o %s" %(params['fasta_path'], output_contig_vcf, params['sam_path'], contig, start, end, phased_bam), verbose=True, output=True, error=True)

from nanocaller.