Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

shorah sampling status #87

Open
snaketron opened this issue Mar 15, 2021 · 1 comment
Open

shorah sampling status #87

snaketron opened this issue Mar 15, 2021 · 1 comment

Comments

@snaketron
Copy link

I am performing local haplotype reconstruction on more than 100 samples with shorah.

Is there a way to check whether the MCMC sampling has converged for a particular region or overall?

@DrYak
Copy link
Member

DrYak commented May 17, 2021

You can inspect the shorah.log file.

  • The MCMC sampling is done by the sub-program diri_sampler, and ShoRAH starts it from run_dmp()
  • For each window, once the sampling is finished, you should see a line like:
DEBUG 2020/09/21 16:49:35 shotgun.py: run_dpm() 223:    run  -i w-NC_045512.2-7840-8040.reads.fas -j 2535 -t 507 -a 0.100000 -K 20 -R 42 finished

(In the above example: NC_045512.2 is the reference sequence/chromosome, 7840-8040 is the region covered by this window, 2535 is the total number of iteration to sample over, 507 is the number of iterations after which that window will start recording history to look for convergence)

Overall, when all sampling is finished, the code will move onto the next part and will switch to output log lines such as:

INFO 2020/09/21 17:02:56 shotgun.py: main() 502:        reading windows for start position 1

To find how many windows there are in total, you can look at this line right before the block of run_dmp():

INFO 2020/09/21 16:49:35 shotgun.py: main() 468:        will run on 375 windows

To get details about the windows, you can look inside file coverage.txt.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants