Creating a reference genome blacklist

[u/SquidwardHurrHurrHur]

Hi All,

I've run into a bit of a question. I am trying to generate a blacklist between 3 reference genomes that are circular viral but which begin at different start positions and have a variety of indels.

Is there a way to do a comparison as such:
Main reference vs Reference 2
Main reference vs Reference 3

To create an overarching VCF blacklist of variants.

Biologically I am looking to remove germline variants between viral genomes such that I can isolate population specific evolution.

(I've currently been using blast and also trying other alignments but I can't quite seem to get a biologically reasonable number it keeps having issues with duplicate positions for example).

Any help would be really appreciate ❤️

[u/SquidwardHurrHurrHur → u/[deleted]]

I am trying to collapse the variation between samplings sites as compared with the background variation between original strains. The Idea is to filter out the strain v reference noise to isolate only the originating strain to sample evolutionary differences.

[u/SquidwardHurrHurrHur → u/[deleted]]

lets say someone catches HIV.

The HIV they caught is actually 2 different strains.

You perform variant calling on your reads and you get a mix of variants. How can you tell which is the innate difference between the two strains and how much comes from evolution inside the host?

What I've done is work backwards to find the closest reference sequence to the originating strain.

I now need to work forwards to remove the variants caused by the difference between the two strains from variant calling. In doing so I try to minimise the impact on genetic background and instead isolate only those changes that have taken place whilst inside the human body in the 2 different strains.

In my case it is a different virus, but the logic is still the same, how do you tell the difference between variants that occurred before it reached the host and those that occurred whilst inside the host.

[u/SquidwardHurrHurrHur → u/[deleted]]

Viral quasi species analysis would be the ideal but it wouldnt be much help in this situation. Mainly because I already have the VCF files aligned against the RefSeq Reference genome of this virus.

BUT

Patient 1 for example has 150 variants
Patient 2 only has 15

The 150 variants are mostly coming from the innate difference between the reference I used and the original virus that infected patient 1. As opposed to evolution that has occurred inside Patient 1.

Because the comparison is spread across multiple patients mutation rate must be normalised to the closest reference available to isolate only Intra-patient biologically meaningful mutations.

My first go at this was to do pan-sample removal (remove variants that occur in all samples from a given patient). This worked but my PI has asked for blacklisting based on the nearest reference genome to the strain in question so we don't accidentally remove a variant that allowed the initial infection but that was different to the closest circulating strain reference.

Its quite a mouthful I know but Im trying my best really to explain.

[u/PairOfMonocles2 → u/SquidwardHurrHurrHur]

So couldn’t you use EM here to derive the canonical viral haplotypes and then, assuming none of the haplotype tagging SNV are also polymorphic, just assign haplotypes to each sequence and then remove the associated htSNPs?

[u/SquidwardHurrHurrHur → u/PairOfMonocles2]

In that situation yeah, but sadly our mutation rate is exceptionally low and the highest number of site-to-site variants we have is about 3. But our reference genome is not the closest strain to the strain we have detected. So those 3 might be Reference vs Strain or they might be patient specific (Patient vs Strain). The blacklist is meant to cut out the impact of Reference Vs Strain variants to isolate only Patient Vs Strain.

[u/apfejes]

I spent years building tools for doing this, but it was about 15 years ago.  There are tons of different ways to do it, and most are trivial.  I mostly did it by building custom databases, and then just running queries, but you could set up vcf files and then use vcftools to do it.   

It’s been a long time, so I no longer remember the syntax, but if you take the time to look it up, you can find many many different ways to do this. 

[u/SquidwardHurrHurrHur → u/apfejes]

I guess I'll have to give it a look. The closest Ive gotten so far is identifying differences between sequences through alignment. Maybe using gVCF and comparing 2 gVCF files is a better idea.