brentp / slivar

genetic variant expressions, annotation, and filtering for great good.
MIT License
248 stars 23 forks source link

compound het variants on non trios #79

Closed atimms closed 3 years ago

atimms commented 3 years ago

Thanks for this great tool. I'm just transitioning over from gemini and trying to figure things out.

I was wondering if you had any suggestion on calling compound heterozygotes on singletons and duo.

Andrew

brentp commented 3 years ago

Hi, you can do this by calling variants that are het in the kid (and het or hom-ref in the parent) and then using --allow-non-trios to slivar compound-hets. You can see more about this here

atimms commented 3 years ago

thanks for getting back to me so quickly.

I'm having issues trying to get the correct variants from the VCF, I tried:

--family-expr 'comphet_side:fam.every(function(s) {return s.het == s.affected && (s.het == !s.affected || s.hom_ref == !s.affected) && s.GQ > 20 })'

and only got variants that were het in the child and hom_ref in the parent.

so I then tried just calling het in all samples i.e.

--family-expr 'comphet_one:fam.every(function(s) {return s.het == s.affected && s.het == !s.affected && s.GQ > 20})',

and again got no variants.

Any help would be great.

Andrew

brentp commented 3 years ago

I think you'd just want:

--family-expr 'comphet_side:fam.every(function(s) {return (s.het || s.hom_ref) && s.GQ >= 20 })'
atimms commented 3 years ago

That worked perfectly.

I don't really understand how that retrieved only heterozygous variants in the proband but it did.

Thanks again.

Andrew

brentp commented 3 years ago

that should get variants that are hom_ref in the proband and het in the parent. you can also require that it is het in the proband with:

--family-expr 'comphet_side:fam.every(function(s) {return (s.het || s.hom_ref) && s.GQ >= 20 }) && fam.some(function(s) { return s.het && s.affected })'
atimms commented 3 years ago

That makes sense, thanks for the explanation.

Andrew

atimms commented 3 years ago

Sorry for all the questions. The compound het function is working really well for singletons and duo (as well as trio's).

But I'm having an issue when there are multiple affected individuals but no parents. I get all het variants in the affected kids but slivar compound-hets doesn't find any compounds hets.

Any ideas?

brentp commented 3 years ago

can you find an example where there are 2 het variants in the same gene in an affected kid ? if you can share a sample VCF with that (you can change positions and id's to anonymize), then i can have a look. i just ran a test and it should find "compound hets" in a single affected sample.

atimms commented 3 years ago

I have some examples.

This VCF has 2 affected siblings, and the corresponding ped file.

I used the command --family-expr 'comphet_side:fam.every(function(s) {return (s.het || s.hom_ref) && s.GQ >= 20 })'

and then used slivar compound-hets on the resulting vcf.

If see multiple het variants in many genes including OR4C5, TAS2R19, TOE1, STIL and PAPPA2.

Thanks.

Andrew

LR09-416.txt

LR09-416.intersected.bcftools.GRCh37_87.vcf.gz

brentp commented 3 years ago

thanks for reporting and the test-case.

that file is not decomposed properly. There are lots of warnings from slivar you can look at any of those variants and see that the AD field is not correct. That said, this was a bug in compound hets where if parents were specified in pedigree file and not in vcf then it would error. (were you able to get it to run somehow?)

If you'd try the attached binary and let me know if it resolves your issue, I would appreciate it. slivar_dev.gz

atimms commented 3 years ago

Thanks for the response, and the updated binary file.

The binary resolved the issue, so now I have comp het calls for those cases where the parents are in the ped file.

As for the issue withe the AD field, I'm aware of the issue which was cause by the original VCF file being normalized without the header being adjusted. Unfortunately I don't have access to another version of the VCF, and I tried to correct using bcftools but couldn't, so it look like I'll have recall variants all my pedigrees.

Andrew

brentp commented 3 years ago

got it. thanks for testing it out. i'll get a new release out with this soon.