-
Notifications
You must be signed in to change notification settings - Fork 4
/
README
57 lines (33 loc) · 2.45 KB
/
README
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
################################################################################
betatree is a collection of python scripts to generate trees from the beta coalescent ensemble, calculate properties of those trees and gather statistics across many instances.
Authors: Taylor Kessinger and Richard Neher
Contact: [email protected]
Example output of the betatree generator are given as pdf files in example_trees. Example site frequency spectra for different parameters of the beta-coalescent ensemble are provided in example_SFS.
If you use betatree in a publication, please refer to
Neher, Kessinger, and Shraiman. "Coalescence and Genetic Diversity in Sexual Populations under Selection." PNAS 110: 15836-41. doi:10.1073/pnas.1309697110.
################################################################################
Tree generation
the script src/betatree.py provides a class that generates coalescent trees using pseudorandom numbers given an initial sample size n and a parameter alpha of the beta measure of the Lambda coalescent process. The following 3 lines will generate a tree of a sample size of 100 and draw it with the Biopython.Phylo package.
myT = betatree(100,alpha = 2)
myT.coalesce()
Phylo.draw(myT.BioTree)
the tree is internally stored as a biopython.phylo tree with all the associated functionality.
Sample code is appended to the definition of the class as will be exectuted if betatree.py is run as main.
################################################################################
Site frequency spectra
the script src/sfs.py uses the class betatree to generate many trees and calculate the SFS assuming that mutation are uniformly distributed on the tree. The following three lines will generate an SFS for a sample size 100 and alpha=1.5 by averaging 1000 trees.
mySFS = SFS(100,alpha=1.5)
mySFS.getSFS(ntrees=1000)
The sfs is accessible as mySFS.sfs and can be binned using different binning schemes or a user defined binning.
mySFS.binSFS(mode='logit', bins=20)
plt.plot(mySFS.bin_center, mySFS.binned_sfs)
The sfs can be saved and loaded by member functions.
Sample code is appended to the definition of the class as will be exectuted if sfs.py is run as main.
################################################################################
Dependencies:
python 2.7
numpy
scipy
biopython
matplotlib for plotting examples
################################################################################