Again on vcf format

[linsson]

Dear @singing-scientist,

Thank you for developing and sharing SNPGenie.

I would like to ask for some clarification regarding the appropriate input format, and I apologize if this has already been discussed in the past.

We are studying the domestication of a diploid plant species. Specifically, we are interested in the genetic basis of certain fruit traits that vary across groups of plants cultivated in different geographic regions. After performing association studies, we identified a set of candidate genes, and we are now aiming to calculate Ka/Ks ratios for these genes within each group.

For genotyping, we used GATK on sequencing data from N separate individuals, so we do not have pooled sequencing data. Instead, we generated a multi-sample VCF with individual-level genotypes.

Our plan was to split the main VCF by group using bcftools, recalculating AF and AN using fill-tags, and then run SNPGenie with --vcfformat=1 on each group-specific VCF. Each input vcf has N_i individuals (SUM(N_i)=N).

We are using the product_results.txt file and calculating Ka/Ks as mean_dN_vs_ref/mean_dS_vs_ref

Could you please confirm whether this is the correct approach for computing Ka/Ks with SNPGenie in our case?

Thank you very much

[singing-scientist]

Thanks very much @linsson!

Which is the 'correct' approach will depend on the hypothesis, e.g., do you hypothesize that whole-gene dN/dS will differ between the two groups?

Note that the statistic you mention, "mean_dN_vs_ref/mean_dS_vs_ref", is specifically computing the mean dN and dS divergences from the reference sequence (i.e., "vs_ref"), thus it's quantifying divergence of each group from the reference rather than diversity within each group. Is that your hypothesis, i.e., does the reference sequence have a special evolutionary significance for you? If not, and if you're concerned with the amount of diversity/divergence within each group (which seems to be the case), I suggest just using piN and piS instead — this is mathematically the same as mean dN and dS among all sequences in each group; it's "pi" just because it's not corrected for multiple hits (e.g., not JC-corrected), which is fine for close sequences. (Very often, reference sequences are arbitrary and have no special evolutionary significance; one notable exception is when the reference is the known ancestral genome and you want to look for patterns of divergence from that ancestor.)

Assuming you generate a VCF file that somehow provides the AF of every SNV in a group, you should be able to straightforwardly run SNPGenie once for each group (2x). vcfformat 4 seems to be the commonest format people use, but any format will suffice as long as it matches the input specifications.

Let me know if that helps, and feel free to ask any follow-ups!
Chase