Can't find a variant which suppose must have in a vcf file

Question

[Chris]

Hi all,

I get a vcf file after running nf-core pipeline then I annotated ID column using bcftools. But when I search for a variant which know must have because the sample from a person has this mutation, I didn't find it. Would you please suggest what wrong in this case? Thank you so much!

[Ram]

There are a whole lot of reasons this could be happening, including the fact that you could be misdiagnosing it. What are you comparing to decide that your expected variant locus is indeed being called a hom-ref locus?

[Chris]

hom-ref locus means homozygous reference locus? Sorry for my understanding. I am not sure your question.

[Ram]

Yes, it means the locus is REF-REF (if you're dealing with a biallelic organism).

[GenoMax]

Just because a pipeline ran does not mean it produced correct results. If you are sure about the variant you are looking for then you should back track and check intermediate result files. If you have not done this already start by looking at the alignment to make sure the expected base changes are there and they made it to the VCF file.

[Chris]

Thanks GenoMax! I know the chromosome and the gene has this mutation so how can I check the alignment?

[GenoMax]

You must have started with fastq data? There should be BAM alignment files that you can start checking in a genome browser like IGV.

[Chris]

Seem I don't have bam file with nf-core sarek. Thanks GenoMax for your instruction.

[GenoMax]

You should have "cram" files which are the equivalent alignment files. I see them in example for sarek.

https://nf-co.re/sarek/results#sarek/results-ed1cc8499366dcefea216fe37e36c6189537d57b/germline_test/preprocessing/recalibrated/NA12878/

[Chris]

Yes, I have 2 cram files:

./preprocessing/recalibrated/sample_1/sample_1.recal.cram
./preprocessing/markduplicates/sample_1/sample_1.md.cram

Not sure which file I should view on IGV.

[Ram]

View the one with the more recent time stamp - it's from later in the pipeline.

[Chris]

Thanks @RAM! Which a gene length about 50k bases, I am not sure how to identify the variant. Is that a variant A -> G? Would you explain why other reads are grey but this is red?

[GenoMax]

Yes that seems to be a variant.

Color of reads explained in: Meaning of read color on IGV

[Chris]

Thanks Max! So red mean large insertion but not a point mutation in this case?

[Chris]

I used other tools such as samtools and some variants show in both vcf file but some only in one file. How can I verify which tools results the correct vcf?

[Pierre Lindenbaum]

look at the bam with e.g. IGV