kaks.py script produces inflated SNP counts

[jfaberha]

Thank you for the useful script, but I noticed an error in my output for kaks.py. Because the script totals the numbers of missense and sense SNPs per gene, I found the counts to be redundant for genes that have alternatively spliced transcripts annotated in the gff/gtf file used for SnpEff. For example, two missense SNPs in one gene were each counted 12 times (once per isoform that included the exon in question) giving me a kA of 24 and a kS of 0. I was able to modify the script on my own to give me kA/kS for each feature (mRNA accession) as opposed to each gene, but it might also be useful to output kA/kS per gene while accounting for SNP redundancy.

[naborlozada]

Hi @jfaberha,

I noticed his issue too and I was trying to modify the script (with no success so far). I was wondering if you could share your solution?
Thanks.

[jfaberha]

My solution was a quick and dirty one where I just swapped out "geneid_index" for "feature_id_index", which gave me ka/ks values for each unique RNA isoform instead of each gene. I would imagine one possible solution for users to figure out ka/ks per gene would be something like a bedtools merge type of step to modify the gff prior to running SNPeff, but to check for redundancy the output of your script would be more complicated. I'm not sure the best way to do that within your script though. -Josh

…

On Wed, Jul 7, 2021 at 12:13 AM naborlozada ***@***.***> wrote: Hi @jfaberha <https://github.com/jfaberha>, I noticed his issue too and I was trying to modify the script (with no success so far). I was wondering if you could share your solution? Thanks. — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#31 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AEHUDFNMFS4AQJXCNJUGYFDTWP5CHANCNFSM4ZUMSM7Q> .

[naborlozada]

I see. Thanks for your reply, Josh, worth it to know how you solved it, even if was "quick and dirty". If I come with a solution, I'll share it here.
Best,
Nabor