search

ANGSD / angsd

Program for analysing NGS data.
227 stars 50 forks source link

Inbreeding, HWE and SFS for heterozygosity #580

Open irbinveliz opened 1 year ago

irbinveliz commented 1 year ago

Hello,

I want to obtain heterozygosity values for my individual samples. I have a doubt regarding on how to account for HWE when doing this.

From what I have seen, the standard procedure is to first calculate the genotype likelihoods and then the SFS. I have seen on the website that this is done as follows:

angsd -i my.bam -anc reference.fa -dosaf 1 -gl 1 realSFS angsdput.saf.idx > est.ml

Then, checking the section for SFS estimation, I saw that by using -doSaf 1, it assumes HWE. So what would happen if there is no HWE?

I have seen that there is the option to calculate per site posterior probabilites of the SFS (doSaf 2). However, if I understand it correctly, this wouldn't be the site allele frequency likelihood but the posterior probabilities of these values. So can we still calculate the realSFS with the produced output and the get the heterozygosity values as above?

Thank you in advance for any help!

irbinveliz commented 1 year ago

To continue the thread,

I have estimated the inbreeding coefficients using ngsF and now I want to incorporate this information to the estimation of the allele frequency spectrum. I am using therefore dosaf 2. However, I am encountering the following error:

[glfReader.cpp:fetch():34] Error reading full chunk: bytesRead:3840 expected:6400 will exit

What could be the reason behind this? I have searched on the issues section and I found one similar issue but there was no clear solution. This is what I am running in ANGSD:

angsd -glf nomaf.genolike.glf.gz -fai 90p-consensus.fasta.fai -anc 90p-consensus.fasta -nInd 80 -domajorminor 1 -domaf 1 -doPost 1 -doGeno 32 -indF 2_testF.indF -out gl.inbr

A question that I have here is also why do we need the fai file? When I run -dosaf 1 as a test I didn't include it and it run well. I also saw here that the argument `-issim 1´ can be added but I couldn't find its meaning.

Thank you for any advice!