This is a test version of Biostars. For the public version, visit https://www.biostars.org.
calling variants between two strains

(newbie alert)

I'm trying to find the sites of variants between two strains. I've map them both to the same reference. However, I don't know what the next step should be. I've used freebayes to call the SNP between the strains and the reference, but I don't know how to find the variants between the two strains themselves. Should I manipulate the 2 vcf files produced by calling SNP between the strains and reference, or is there a tool that I can use to call variants between the 2 strains themselves?

Thanks a lot!

snp vcftools samtools freebayes
• 4,468 views • link
updated by apelin20 apelin20 • written by steven.lin.pur • steven.lin.pur

The final approach I used for this problem is simple and though the problem itself is likely elementary for many people here. I think I'd still post my solution:

  1. use samtools mpileup or GATK haplotypecaller (or other similar tool) to get the likelihood of every position of your strains.
  2. call variant between strain and the reference and merge your resulting vcf files. You should obtain a vcf file with the union of the all the variant sites across all strains. At this stage, you don't know if the ., ./. or .|. in the genotype field means (a) that there is a lack of variant between the reference and a particular strain or (b) there is a lack of data for a particular strain.
  3. call genotype again, but call genotype regardless of the presence of a variant or not for each position that exist in the merged file you obtained from step 2. Merge the resulting vcf files and you should get a vcf contains all the variating sites. And the lack of data for a particular strain should be presented as ., ./. or .|. in the genotype field while the lack of variant should be presented as 0, 0/0 or 0|0. Then, you'd know all the variants between the two strains.

This is the easiest solution I can think of. Hope it helps if you happen to be a newbie like me (:

• 0 views • link
updated by Ram Ram • written by steven.lin.pur • steven.lin.pur

2 answers

You need to merge (vcf-merge) the two VCF files you have. Take a look at vcf-compare also. It will tell you how many differences to reference your two strains have and how many differences are unique to each strain.

• 0 views • link
written by apelin20 apelin20

Thank you for the advice. I have some problems when practicing it however, can you please help me?

When I merge and compare the two VCF files, I can't tell if there is a difference between the two strains when the following case occur:

Reference has a genotype A at locus 100;

Strain i has a genotype C at locus 100;

Strain j does not have a variant found at locus 100;

in such a case I can't tell if strain j is different from strain i at locus 100 or not. Because the lack of variant in strain j might be due to either lack of data(in which case it might actually be the same as strain i) or lack of variation from the reference(in which case its genotype is A and is a SNP against strain i).

I'm planning on using the SNPs to perform QTL mapping. I thought that genome wide SNP calling for QTL mapping is a common method used and that there should be a tool that calls the SNP between two strains directly, but I've looked up Gatk, samtools, freebayes and I can't seem to find any of it doing this...

• 0 views • link
updated by Ram Ram • written by steven.lin.pur • steven.lin.pur

Why not use Mummer to align two strains directly?

• 0 views • link
updated by Ram Ram • written by poisonAlien poisonAlien

Thank you for the great advice!

it is my understanding that Mummer needs a reference too to perform alignment.

are you suggesting that I map one of the strain to the reference, creating a aligned genome of one strain, then used that aligned genome as a reference to map the other strain? Thanks again

• 0 views • link
updated by Ram Ram • written by steven.lin.pur • steven.lin.pur

Hi, Mummer is better option if you have reference genome for both strains. We have used it in past to align two starins of yeast (S288C and YJM789) to determine genome wide difference. In your case if you dont have reference genomes for both strains, then you can do as suggested in above answer (vcf compare).

Or if you dont mind trying out different options, assemble genomes for both strains using velvet and use them to compare using Mummer.

• 0 views • link
updated by Ram Ram • written by poisonAlien poisonAlien

Log in to answer this question.

More posts like this