본문 바로가기
Tools

LDpop

by wycho 2020. 7. 12.

 

Coalescent theory에 관한 글에서 coalescent probability는 lineage, n에 비해 population, N이 충분히 크다면 2-lineage로 근사하여 계산한다. 즉, 한 세대에 1쌍의 coalescence가 일어난다고 생각할 수 있다. 따라서 2-locus linkage로 부터 mutation이나 recombination rate를 구할 수 있다.

 

LDpop은 recombination rate를 구하기 위해 two-locus의 likelihood를 계산하여 lookup table을 만들어주는 프로그램이다. 다른 프로그램과 다른 점은 population size가 변해도 정교하게 계산 가능하다는 것이다.

 

- Repository : https://github.com/popgenmethods/ldpop

- Poster : jackkamm.github.io/assets/ldpop_poster.pdf

- Paper : Two-Locus Likelihoods under Variable Population Size and Fine-Scale Recombination Rate Estimation, https://www.genetics.org/content/203/3/1381

 

INSTALLATION

$ git clone https://github.com/popgenmethods/ldpop

$ cd ldpop

$ pip install .

$ ./run/ldtable.py --help

더보기
usage: ldtable.py [-h] -n N [-th theta] -rh num_rh,max_rh [-s s0,s1,...,sD]
                  [-t t1,...,tD] [--approx] [--cores CORES] [--log logfile]
                  [-L total_length] [-S seg_sites]

Print a lookup table, as expected by ldhat or ldhelmet.

optional arguments:
  -h, --help         show this help message and exit
  -n N               sample size
  -th theta          twice the mutation rate
  -rh num_rh,max_rh  grid of rhos (twice the recomb rate). The grid has num_rh
                     uniformly spaced points from 0 to max_rh, inclusive.
                     (((Alternatively, to create a non-uniform grid, use "-rh
                     r0,step0,r1,step1,r2,...rK". This creates a grid
                     {r0,r0+step0,r0+2*step0,...,r1,r1+step1,...,rK} similar
                     to ldhelmet. Note that non-uniform grid is incompatible
                     with vanilla ldhat.)))
  -s s0,s1,...,sD    coalescent scaled population sizes (s0=present size,
                     sD=ancient size)
  -t t1,...,tD       times of size changes from present backwards. Must be
                     increasing positive reals.
  --approx           use finite moran model. A reasonable approximation that
                     is much faster than the exact formula. Accuracy of the
                     approximation can be improved by taking n larger than
                     needed, and using ldhat/lkgen.c to subsample. (Converges
                     to the "exact" diffusion as n->infty)
  --cores CORES      Number of parallel processes to use.
  --log logfile      Log extra info to logfile. If logfile=".", logs to
                     STDERR.
  -L total_length    Total length of DNA sequence
  -S seg_sites       Segregating sites

 

NOTATION

- θ/2 : the mutation rate per locus per unit time.

- P=(P_ij), i,j A : the transition probabilities between alleles given a mutation.

- A = {0,1} : the set of alleles.

- ρ/2 : the recombination rate per unit time.

- * : missing alleles.

- a haplotype : observed only at the left locus. ex) [(0, -1), (1, -1)]

- b haplotype : observed only at the right locus. ex) [(-1, 0), (-1, 1)]

- c haplotype : observed at both loci. ex) [(0, 0), (0, 1), (1, 0), (1, 1)]

- n = {n_i*,n*j,n_ij} : number of alleles.

- M_t : a collection of N two-locus haplotypes at time t (with no missing alleles). Then M_t is a Markov chain going forward in time that changes due to mutation, recombination, and copying events.

- Λ^d_(N) : the transition matrix of M_t within (t_d,t_d+1]. 

- η_d : the population size.

- 1/η_d : the coalescent rate within the time interval (t_d,t_d+1].

- Type : {n00, n01, n10, n11}의 가능한 조합의 수. (let h=n//2, type= 1 + h + h(h-1)(h+4)//6 + (h-1)(h+2)//2

- Tridiagonal_matrix : a band matrix that has nonzero elements on the main diagonal.

- The importance sampler yields a posterioi sample of two-locus ARGs.

 

The state of M_t changes due to

- Copying event with rate 1/2η_d.

- Mutation with rate θ/2.

- Recombination (c → a,b) with rate ρ/2.

- Recoalescence event (a,b  c) with rate 1/η_d.

 

n_t is a backward-in-time Markov chain :

P(n_s1|n_s2,n_s3)=P(n_s1|n_s2) for -∞<s1<s2<s3≤0

 

EXECUTION

- Lookup table을 구하는 방법으로는 exact하게 계산 하여 구하는 방법과 N이 충분히 크다고 가정하고 근사하는 방법이 있다. 근사하여 구하는 방법은 multiprocessing이 되어 계산이 빠르다.

 

$ ldtable.py -n 10 -th 0.001 -rh 11,10

더보기

$ cat lk_n10_th0.001_rh11_10.lk

10 50
1 0.001
11 10.0


1 # 9 0 0 1 : -17.674165163339033 -18.11665728607374 -18.37499013353314 -18.546993444456426 -18.671091609254574 -18.76566262526687 -18.840635516359647 -18.90186011565102 -18.953017628764186 -18.99654829342412 -19.03413932069867
2 # 8 1 1 0 : -19.708043132709598 -19.690740080458937 -19.677024376197814 -19.66692041474687 -19.65928974893083 -19.65333706221456 -19.64855847312643 -19.644630404253192 -19.641338192838596 -19.638534477576158 -19.636114852440762
3 # 8 1 0 1 : -21.78156318282099 -21.525183186401947 -21.496238928748433 -21.504655775141813 -21.52056772598629 -21.536824526264276 -21.551584304549273 -21.564506450873942 -21.57569569556118 -21.585372187589286 -21.593764310783172
4 # 7 2 1 0 : -21.786715118115428 -21.762917445613088 -21.745722297729024 -21.7334478645285 -21.72442452377522 -21.717587623953666 -21.712270145611882 -21.70804256637116 -21.704618931897013 -21.701802777034217 -21.699455267653715
5 # 8 0 0 2 : -19.383305057491032 -20.007386174024997 -20.413368533885045 -20.706691315539164 -20.93232346256279 -21.113489506316963 -21.26358443541618 -21.390939299227004 -21.501036718319437 -21.597649303650893 -21.683469142588653
6 # 7 1 1 1 : -30.383045732624304 -25.84496913932525 -25.23208888895376 -24.93563324627234 -24.755486385405597 -24.632135839331752 -24.541116347573684 -24.470436056987992 -24.413495016575432 -24.36634311040092 -24.326459699039653
7 # 6 2 2 0 : -23.731899423231283 -23.70217675067733 -23.686452967884904 -23.67769235902812 -23.672846534733765 -23.670358452500544 -23.66935310039161 -23.66930617521522 -23.669888659745574 -23.670886811983628 -23.672158052163297
8 # 7 2 0 1 : -23.037191108022473 -22.908068426041535 -22.88090061380508 -22.876653024981184 -22.87911707506993 -22.883533506870226 -22.888281106166875 -22.89278645158421 -22.896865118496603 -22.900484005831256 -22.903668091797666
9 # 6 3 1 0 : -23.039105990884607 -23.008134939162183 -22.988794369741925 -22.97581340238414 -22.966658086656523 -22.95996118710668 -22.954922798117025 -22.951046027013184 -22.948007569415694 -22.94558912913157 -22.943638765477893
10 # 7 1 0 2 : -23.72464259022971 -23.350280259658486 -23.37216786495492 -23.446612324478618 -23.527078621022557 -23.603085265451345 -23.672341291002244 -23.73483221222514 -23.791181261450777 -23.842136115982935 -23.888404590378283
11 # 6 2 1 1 : -31.921267240102242 -27.373870350204324 -26.72583155898064 -26.393828689057674 -26.184135475667073 -26.03693350917254 -25.926634546162678 -25.840223121478044 -25.770298419562902 -25.712306556609906 -25.663274271718944
12 # 5 3 2 0 : -24.830100957930473 -24.79099810330199 -24.779060544722608 -24.777004359699074 -24.779419678793143 -24.784033916550445 -24.789738672026946 -24.795944602194812 -24.80232186097665 -24.808680724691865 -24.81491132934575
13 # 7 0 0 3 : -20.456387021492713 -21.196636984661204 -21.706587664804573 -22.091676336674652 -22.3981117146893 -22.65064148131129 -22.864098637700295 -23.048055279473722 -23.209026881290264 -23.35163378977949 -23.4792617237907
14 # 6 1 1 2 : -32.216319828564046 -27.260925895965332 -26.690780683340282 -26.447689975508165 -26.31960461854008 -26.2445253649386 -26.19757914229449 -26.166968037941157 -26.146461010001506 -26.13250639481429 -26.122959748242724
15 # 5 2 2 1 : -33.43316489743534 -28.762646568390615 -28.088995101729466 -27.73446551863796 -27.50666506305476 -27.34509398217093 -27.22334908016135 -27.127744984206238 -27.05036931034485 -26.986283898006896 -26.932226526244946
16 # 4 3 3 0 : -25.745810340205193 -25.68857248615876 -25.68714816550794 -25.701352842565033 -25.720830499698867 -25.741876800321247 -25.76295126471435 -25.783372629775304 -25.802844023408657 -25.821253526598166 -25.83858180456946
17 # 6 3 0 1 : -23.73122199751227 -23.650859305190217 -23.62620895488314 -23.61698798478473 -23.613497878103644 -23.61244498796976 -23.612494569898146 -23.613046091675024 -23.613811544655384 -23.614647850779328 -23.615483262838445
18 # 5 4 1 0 : -23.731968754400864 -23.6912862080864 -23.670241512550508 -23.65746738218061 -23.64906094067236 -23.643243500395865 -23.63907452072441 -23.636008299818723 -23.633707260475905 -23.631952451242906 -23.630596562767696
19 # 6 2 0 2 : -24.827690316553024 -24.657275253430356 -24.661510006461135 -24.699625494686615 -24.744990929514366 -24.790132474179615 -24.83272236550716 -24.87213112550095 -24.908342947833795 -24.941563375413658 -24.972065848587015
20 # 5 3 1 1 : -32.708378777751506 -28.162257075212175 -27.50449742426639 -27.159891227907416 -26.939005781557796 -26.782512521477074 -26.664616220314173 -26.57199232743767 -26.49696029426157 -26.43473903010456 -26.382176905675855
21 # 4 4 2 0 : -25.340505539010376 -25.283362398431994 -25.276078307261653 -25.283234267099747 -25.295570276100936 -25.30973127729995 -25.32429686884466 -25.338619665220644 -25.35240066929789 -25.36551027328056 -25.377905049192886
22 # 6 1 0 3 : -24.82156434154419 -24.37830358819717 -24.44257482355402 -24.567942174397068 -24.69836033753638 -24.82121281320576 -24.933733248922103 -25.035953415353163 -25.12872392968611 -25.21307576993064 -25.29000520823643
23 # 5 2 1 2 : -33.580428310728934 -28.706470371125654 -28.073650438942934 -27.76807133940261 -27.58610033219086 -27.46546858453705 -27.37991310163718 -27.316304024127604 -27.267322355561255 -27.228562189694188 -27.197213536130903
24 # 4 3 2 1 : -34.09871992932527 -29.354258085982615 -28.683126024103174 -28.331623084953108 -28.10699906392459 -27.94861693181719 -27.830001984244433 -27.7374292555701 -27.662964952506705 -27.60166054750313 -27.550250328050197
25 # 3 4 3 0 : -26.03267518106679 -25.934139654375986 -25.944921839763357 -25.980835240033034 -26.02324324663063 -26.066225631961917 -26.107630119253937 -26.146678996938387 -26.183156192416295 -26.21708755769175 -26.24860261448747
26 # 6 0 0 4 : -21.06079936218715 -21.86859776952598 -22.442758871812135 -22.887564954666242 -23.248856754534 -23.551437237514325 -23.810402839253147 -24.035677405230945 -24.23415246015679 -24.410816068742644 -24.569397935005494
27 # 5 1 1 3 : -33.15020768335405 -27.992299485908838 -27.45066816520753 -27.241120310395473 -27.146085793473368 -27.102236678479358 -27.084342453711685 -27.080574377264426 -27.084798095835758 -27.093608445655825 -27.105022825641043
28 # 4 2 2 2 : -34.88424444721014 -29.99647788077746 -29.306944302551216 -28.950386025980457 -28.725128970563645 -28.56789282895805 -28.45118505782399 -28.360830011899914 -28.28867802024518 -28.229673046560173 -28.180495638239186
29 # 3 3 3 1 : -34.56275080869954 -29.677775447376135 -29.022573209731426 -28.691996118631977 -28.487767618544694 -28.348133061235032 -28.246475542665955 -28.169187023262246 -28.108512972631573 -28.05968579164577 -28.01960392923982
30 # 2 4 4 0 : -26.03078781415632 -25.830253722974724 -25.864915662878502 -25.942358351262968 -26.027841454301992 -26.11159049072308 -26.190519872164792 -26.263790924841473 -26.331403746991157 -26.393676961381907 -26.4510372136087
31 # 5 4 0 1 : -23.95480805368079 -23.899367339539882 -23.876657327438206 -23.86500534945362 -23.858327468517626 -23.854244804136716 -23.851645315706865 -23.84994721265235 -23.8488215448051 -23.848071455973756 -23.84757393261481
32 # 4 5 1 0 : -23.95480805368079 -23.899367339539882 -23.876657327438206 -23.86500534945362 -23.858327468517626 -23.854244804136716 -23.851645315706865 -23.84994721265235 -23.8488215448051 -23.848071455973756 -23.84757393261481
33 # 5 3 0 2 : -25.33986289532523 -25.246748056005924 -25.244560368671227 -25.26415696136341 -25.28991703571838 -25.31683880794493 -25.343081798398995 -25.367957569193084 -25.391247175605987 -25.412933751419256 -25.433087590003833
34 # 4 4 1 1 : -32.95612432116707 -28.41063935126291 -27.750572801843703 -27.402580794057236 -27.17860454203106 -27.019516308021046 -26.89948845753294 -26.805120443070003 -26.728657485523833 -26.66525607788239 -26.611714561043165
35 # 3 5 2 0 : -25.33986289532523 -25.246748056005924 -25.244560368671227 -25.26415696136341 -25.28991703571838 -25.31683880794493 -25.343081798398995 -25.367957569193084 -25.391247175605987 -25.412933751419256 -25.433087590003833
36 # 5 2 0 3 : -25.743317170677468 -25.550113673608497 -25.57679185403195 -25.643989677607127 -25.718972843623103 -25.792556478156182 -25.861866154715205 -25.926135316897298 -25.985375706690018 -26.039889005190478 -26.090071187028418
37 # 4 3 1 2 : -34.143906032176645 -29.337736459070754 -28.67958946272954 -28.343827808870998 -28.13402275928525 -27.98891905480989 -27.88206358437708 -27.799898941753717 -27.73467697445762 -27.681618078722252 -27.637601255358895
38 # 3 4 2 1 : -34.143906032176645 -29.337736459070754 -28.67958946272954 -28.343827808870998 -28.13402275928525 -27.98891905480989 -27.88206358437708 -27.799898941753717 -27.73467697445762 -27.681618078722252 -27.637601255358895
39 # 2 5 3 0 : -25.743317170677464 -25.550113673608493 -25.576791854031946 -25.643989677607124 -25.7189728436231 -25.79255647815618 -25.8618661547152 -25.926135316897295 -25.985375706690014 -26.039889005190474 -26.090071187028414
40 # 5 1 0 4 : -25.331586642152192 -24.85659600999635 -24.94416712048194 -25.097010732039227 -25.254702434571083 -25.40351211733604 -25.540317728168855 -25.665062742075794 -25.778624491230023 -25.88211403365648 -25.97663187449567
41 # 4 2 1 3 : -34.29254048734296 -29.27653289516612 -28.655366936245002 -28.365747151399482 -28.200132440503708 -28.095229256436852 -28.024455359525533 -27.97460141507215 -27.93836358420743 -27.911390062339926 -27.89093974305883
42 # 3 3 2 2 : -35.2608233924436 -30.373172224209007 -29.671787568038475 -29.303082529890098 -29.066956733826814 -28.90025551802121 -28.77534071291354 -28.677850751566872 -28.599464542921908 -28.534980482755365 -28.48095886148686
43 # 2 4 3 1 : -34.29254048734296 -29.276532895166117 -28.655366936245 -28.36574715139948 -28.200132440503705 -28.09522925643685 -28.02445535952553 -27.974601415072147 -27.938363584207426 -27.911390062339922 -27.890939743058826
44 # 1 5 4 0 : -25.331586642152192 -24.85659600999635 -24.94416712048194 -25.097010732039227 -25.254702434571083 -25.40351211733604 -25.540317728168855 -25.665062742075794 -25.778624491230023 -25.88211403365648 -25.97663187449567
45 # 5 0 0 5 : -21.256777724044703 -22.086893829843525 -22.682932562383343 -23.148691245386473 -23.52973686312638 -23.85073788925365 -24.126745720422768 -24.36769547087625 -24.580524465450864 -24.77028683635222 -24.94078928336218
46 # 4 1 1 4 : -33.44292014830852 -28.223599300826123 -27.691751234785507 -27.493585957757475 -27.40985563094509 -27.37679601481365 -27.36903328073057 -27.374700492567356 -27.387671224850962 -27.404570158522084 -27.42345477794246
47 # 3 2 2 3 : -35.30774048158377 -30.35028517070117 -29.659506727152365 -29.305237897875408 -29.08328787412062 -28.92964154619174 -28.816527218640587 -28.72965107885456 -28.66081308522928 -28.60493847049341 -28.558704676367654
48 # 2 3 3 2 : -35.30774048158377 -30.350285170701174 -29.659506727152365 -29.305237897875404 -29.083287874120618 -28.929641546191736 -28.816527218640584 -28.729651078854555 -28.660813085229275 -28.604938470493405 -28.55870467636765
49 # 1 4 4 1 : -33.44292014830852 -28.223599300826123 -27.691751234785507 -27.493585957757475 -27.40985563094509 -27.37679601481365 -27.36903328073057 -27.374700492567356 -27.387671224850962 -27.404570158522084 -27.42345477794246
50 # 0 5 5 0 : -21.256777724044703 -22.086893829843525 -22.682932562383343 -23.148691245386473 -23.52973686312638 -23.85073788925365 -24.126745720422768 -24.36769547087625 -24.580524465450864 -24.77028683635222 -24.94078928336218

$ ldtable.py -n 10 -th 0.001 -rh 11,10 --approx

더보기

$ cat a_lk_n10_th0.001_rh11_10.lk

10 50
1 0.001
11 10.0


1 # 9 0 0 1 : -17.674165164700412 -18.040527234510023 -18.263108571790525 -18.415962993629243 -18.528984326382677 -18.616860680804898 -18.687724106008382 -18.74646635957096 -18.796218781690712 -18.839086098972437 -18.87654006031636
2 # 8 1 1 0 : -19.70804313270812 -19.672339921132753 -19.65442569976029 -19.64314443842394 -19.635156109585786 -19.629085690665093 -19.624252456729913 -19.62027572336303 -19.616923290436134 -19.614044021073724 -19.61153446348321
3 # 8 1 0 1 : -21.781563182826886 -21.656544332049336 -21.627760714184937 -21.621989837422184 -21.62329209462191 -21.626819554982617 -21.630865994679056 -21.63480074767461 -21.63840232524231 -21.641613279391358 -21.64444214452221
4 # 7 2 1 0 : -21.78671511811353 -21.750409911003576 -21.731068907484524 -21.71875066881128 -21.710149346100398 -21.703788003872514 -21.6988910329144 -21.69500709231527 -21.691854051084906 -21.689245911494254 -21.687054867239137
5 # 8 0 0 2 : -19.383305059162655 -19.903422262020197 -20.254757554642904 -20.516877733913457 -20.724045693669677 -20.894221728127174 -21.037960299558005 -21.161956271008837 -21.270699764415628 -21.36733484232771 -21.45414067798116
6 # 7 1 1 1 : -30.38304573262679 -25.14487345257676 -24.750205179150008 -24.56559370353076 -24.45272788121696 -24.374100566631995 -24.314951850257252 -24.268181579035286 -24.229897793106417 -24.197759854832825 -24.170259206180646
7 # 6 2 2 0 : -23.731899423229322 -23.69224333167182 -23.674572442582487 -23.665477063705218 -23.6606257278305 -23.658143173101855 -23.65707917248016 -23.656902795299885 -23.65729705204788 -23.6580634506834 -23.659073208224264
8 # 7 2 0 1 : -23.037191108023055 -22.950925170824636 -22.925569606700268 -22.91662345959463 -22.91350607504013 -22.912737501796695 -22.91296326163143 -22.913601447417843 -22.914387424310444 -22.91519794462012 -22.915976337923265
9 # 6 3 1 0 : -23.03910599088246 -22.99999173561128 -22.97939180301931 -22.966504819006232 -22.95769959701795 -22.95134242204208 -22.946572533896813 -22.942888949193012 -22.939979234125524 -22.937638313939637 -22.935726127738945
10 # 7 1 0 2 : -23.724642590241746 -23.57712328983502 -23.591926215526325 -23.63929366037667 -23.69340435303887 -23.746872594901408 -23.797386245316527 -23.844330928299176 -23.88771098302436 -23.927754823591417 -23.964759301377693
11 # 6 2 1 1 : -31.921267240102818 -26.98669377581034 -26.468409522369694 -26.201496709081564 -26.03026689580918 -25.908253613605925 -25.81566368281806 -25.742366991101274 -25.682546126873415 -25.632582398915044 -25.5900878628614
12 # 5 3 2 0 : -24.830100957928675 -24.78928840914448 -24.777824916303747 -24.77617258531158 -24.778674373183907 -24.783108060158423 -24.78844716227697 -24.794169974781955 -24.79999580608828 -24.805769101622257 -24.811403227754024
13 # 7 0 0 3 : -20.456387023326087 -21.07421386104646 -21.514701973779196 -21.857610638991545 -22.137890719770283 -22.374333149738966 -22.57830970257153 -22.757258061753937 -22.91631925511036 -23.059197934650285 -23.188652510901626
14 # 6 1 1 2 : -32.216319828570256 -26.472239657198987 -26.169381326054957 -26.0628134228558 -26.016801501697113 -25.996026257635474 -25.987430721321335 -25.985232017676044 -25.986563749441817 -25.989891687598696 -25.994348269720266
15 # 5 2 2 1 : -33.43316489743556 -28.403592396341104 -27.82567614889313 -27.52412024136246 -27.330531486281956 -27.19310957971883 -27.089418122508782 -27.00786401635547 -26.941753650070666 -26.8869116564476 -26.84057985148175
16 # 4 3 3 0 : -25.74581034020419 -25.702704309955994 -25.705308324087756 -25.720229137873428 -25.73895251842308 -25.758542998313974 -25.777839190385738 -25.796361754407602 -25.813925198172388 -25.830479372212512 -25.846038307403173
17 # 6 3 0 1 : -23.731221997511113 -23.666689652659926 -23.64311040808358 -23.63208616648535 -23.62621140279422 -23.62281485592456 -23.62074024869625 -23.619423665711693 -23.618565702695264 -23.61799690675037 -23.61761646508662
18 # 5 4 1 0 : -23.731968754398682 -23.687980387474905 -23.66639442099277 -23.653528702442994 -23.64506193729105 -23.639140398356044 -23.63482261608593 -23.631576304490654 -23.629077286416166 -23.627117004445136 -23.625555494645415
19 # 6 2 0 2 : -24.827690316555973 -24.73972469378161 -24.750102605971207 -24.782091472880055 -24.8189565075388 -24.85567874989878 -24.890601784497186 -24.9232310374026 -24.953511374216397 -24.98155653515331 -25.007539409113967
20 # 5 3 1 1 : -32.708378777751236 -27.865163125619496 -27.29610607168241 -26.998516799794402 -26.806820402586695 -26.670249292203962 -26.566833297579564 -26.48521728518749 -26.4188397047581 -26.363603078165653 -26.316797317156464
21 # 4 4 2 0 : -25.340505539009147 -25.294824448892275 -25.29070828320594 -25.298386513579697 -25.31010524265587 -25.323124188611562 -25.33631873823553 -25.349194869450454 -25.361535586160112 -25.373255016190715 -25.38433207439564
22 # 6 1 0 3 : -24.821564341560403 -24.66130668380474 -24.708471487652346 -24.794786004683722 -24.888115981860878 -24.979296058710375 -25.065444311287678 -25.14580366811547 -25.220408711701978 -25.289591284126814 -25.353780808121627
23 # 5 2 1 2 : -33.58042831073169 -28.269189945383168 -27.803756218984084 -27.579723745965204 -27.44450224280935 -27.35344511810972 -27.287907689864554 -27.2385427727424 -27.200097452544135 -27.169376824172097 -27.14432076662942
24 # 4 3 2 1 : -34.09871992932563 -29.03198417242239 -28.452152038699722 -28.150541498635143 -27.957773927727935 -27.821628035955037 -27.71945172461635 -27.63953044021249 -27.575098498464428 -27.521936575510395 -27.477259441481223
25 # 3 4 3 0 : -26.032675181067457 -25.97791528750709 -25.996393524089736 -26.031647549364823 -26.070625168348396 -26.10923773062446 -26.146069736911915 -26.180665821979897 -26.212951059113447 -26.24300644945603 -26.27097412965468
26 # 6 0 0 4 : -21.06079936410463 -21.735149138548348 -22.229957451863317 -22.62439282364456 -22.953145796259296 -23.234945896715292 -23.481220093881976 -23.69952685508991 -23.895162203298387 -24.072000866378684 -24.232974905473174
27 # 5 1 1 3 : -33.15020768336262 -27.179859935309757 -26.926422830019696 -26.864195580080356 -26.857857773569293 -26.87260417027896 -26.89587355755647 -26.92233298497961 -26.949513678737937 -26.97622811597509 -27.00190897901861
28 # 4 2 2 2 : -34.884244447211955 -29.55632333371561 -28.991641059920497 -28.70321322323509 -28.521736012178692 -28.395301010732794 -28.301539197587907 -28.228965832671108 -28.170998028564092 -28.12356183437715 -28.083988332304152
29 # 3 3 3 1 : -34.56275080870113 -29.331749959746716 -28.79309012997242 -28.52195599643826 -28.3534875951665 -28.23754455200556 -28.152614833208013 -28.087687486074103 -28.036467906104455 -27.995070910718034 -27.960958552181708
30 # 2 4 4 0 : -26.030787814161243 -25.943830087335904 -25.98625419141336 -26.054752504792212 -26.12767560878453 -26.198656717501088 -26.2656670616326 -26.328174052588977 -26.38621070723088 -26.44002576855719 -26.489940729704095
31 # 5 4 0 1 : -23.954808053678885 -23.9029930960025 -23.880500636333462 -23.868228522242017 -23.860683322158568 -23.855689251309126 -23.8522136659688 -23.849705355369213 -23.847844705635634 -23.846434627580383 -23.845347658960957
32 # 4 5 1 0 : -23.954808053678885 -23.9029930960025 -23.880500636333462 -23.868228522242017 -23.860683322158568 -23.855689251309126 -23.8522136659688 -23.849705355369213 -23.847844705635634 -23.846434627580383 -23.845347658960957
33 # 5 3 0 2 : -25.339862895325272 -25.281082730920712 -25.284569960339038 -25.303377024625465 -25.326352451699744 -25.349904848412887 -25.372727480033628 -25.394345612779727 -25.41462097043053 -25.43355825485313 -25.451222589118164
34 # 4 4 1 1 : -32.95612432116654 -28.134686620462915 -27.552421362938873 -27.24696347900646 -27.050041668792993 -26.909762771454147 -26.803596701683503 -26.719875047315625 -26.65184594783448 -26.59528914185183 -26.54741161665568
35 # 3 5 2 0 : -25.339862895325272 -25.281082730920712 -25.284569960339038 -25.303377024625465 -25.326352451699744 -25.349904848412887 -25.372727480033628 -25.394345612779727 -25.41462097043053 -25.43355825485313 -25.451222589118164
36 # 5 2 0 3 : -25.743317170681895 -25.655983482063338 -25.690185783941676 -25.74931141191045 -25.81293649593069 -25.875041325243654 -25.933702745644272 -25.98840880262721 -26.039179312400538 -26.086235863816132 -26.129867482582984
37 # 4 3 1 2 : -34.143906032177746 -29.0009467432225 -28.45389466375061 -28.17524220250504 -28.000029618690434 -27.878001205386873 -27.787551327437612 -27.717593431739388 -27.66177093242225 -27.61614710461067 -27.578140633741697
38 # 3 4 2 1 : -34.143906032177746 -29.0009467432225 -28.45389466375061 -28.17524220250504 -28.000029618690434 -27.878001205386873 -27.787551327437612 -27.717593431739388 -27.66177093242225 -27.61614710461067 -27.578140633741697
39 # 2 5 3 0 : -25.74331717068189 -25.655983482063334 -25.690185783941672 -25.749311411910448 -25.812936495930686 -25.87504132524365 -25.93370274564427 -25.988408802627205 -26.039179312400535 -26.08623586381613 -26.12986748258298
40 # 5 1 0 4 : -25.331586642170517 -25.16538635672111 -25.229461199657212 -25.33613464675064 -25.450207743069637 -25.56159631185714 -25.667072417708287 -25.765744694591543 -25.857607769922925 -25.943000429737392 -26.022383671239734
41 # 4 2 1 3 : -34.29254048734688 -28.826178203593486 -28.38615489890133 -28.183973326695714 -28.06791905775748 -27.993897105608234 -27.943635957137637 -27.90806140358079 -27.882131197870663 -27.862818912399256 -27.848202474299853
42 # 3 3 2 2 : -35.26082339244509 -29.964249471459127 -29.375859324337107 -29.06941794001051 -28.87361423759956 -28.735463132561527 -28.631929309973913 -28.55108027308851 -28.48601627316087 -28.432430986411063 -28.387480499220786
43 # 2 4 3 1 : -34.29254048734688 -28.826178203593482 -28.386154898901328 -28.18397332669571 -28.067919057757475 -27.99389710560823 -27.943635957137634 -27.908061403580785 -27.88213119787066 -27.862818912399252 -27.84820247429985
44 # 1 5 4 0 : -25.331586642170517 -25.16538635672111 -25.229461199657212 -25.33613464675064 -25.450207743069637 -25.56159631185714 -25.667072417708287 -25.765744694591543 -25.857607769922925 -25.943000429737392 -26.022383671239734
45 # 5 0 0 5 : -21.256777725988584 -21.949750915191473 -22.46292645786226 -22.875223321508116 -23.22116802668468 -23.519385852272208 -23.781241827985617 -24.014262731861535 -24.223734669681036 -24.413538159823574 -24.586622352960056
46 # 4 1 1 4 : -33.442920148317896 -27.40606719657937 -27.16862621300563 -27.1210441724704 -27.12803318360568 -27.15487251967614 -27.18910171809251 -27.22549203964848 -27.261675964811978 -27.296560982574043 -27.329667035500457
47 # 3 2 2 3 : -35.307740481586144 -29.88861088784387 -29.33001383525542 -29.047765786211976 -28.872121535823013 -28.75107744229613 -28.66226162211739 -28.594216005131628 -28.540393197232525 -28.496755834109187 -28.460669691229164
48 # 2 3 3 2 : -35.307740481586144 -29.88861088784387 -29.330013835255418 -29.047765786211972 -28.87212153582301 -28.751077442296125 -28.662261622117388 -28.594216005131624 -28.54039319723252 -28.496755834109184 -28.46066969122916
49 # 1 4 4 1 : -33.442920148317896 -27.40606719657937 -27.16862621300563 -27.1210441724704 -27.12803318360568 -27.15487251967614 -27.18910171809251 -27.22549203964848 -27.261675964811978 -27.296560982574043 -27.329667035500457
50 # 0 5 5 0 : -21.256777725988584 -21.949750915191473 -22.46292645786226 -22.875223321508116 -23.22116802668468 -23.519385852272208 -23.781241827985617 -24.014262731861535 -24.223734669681036 -24.413538159823574 -24.586622352960056

[ s: population size, t: times of size changes ]

$ ldtable.py -n 10 -th 0.001 -rh 11,10 -s 10,30,20,50 -t 1,3,5

더보기
10 50
1 0.001
11 10.0


1 # 9 0 0 1 : -10.370719622677827 -12.2900734519261 -12.503506561105914 -12.608623987921407 -12.672706643385647 -12.715517292911942 -12.745840936999132 -12.768283113849352 -12.78548123018527 -12.799039859321207 -12.809982611014322
2 # 8 1 1 0 : -13.179439158809531 -12.968578763686173 -12.947310551801811 -12.93726405301962 -12.931285008091399 -12.927356224248674 -12.924606933782128 -12.922590857165202 -12.921056864405857 -12.91985426622154 -12.918888001960438
3 # 8 1 0 1 : -14.57463107442854 -14.556574467239182 -14.606678907782145 -14.632275795730477 -14.648402357862524 -14.659335906629616 -14.667099544911812 -14.67282135593099 -14.67717393249262 -14.680575890549855 -14.683297228166706
4 # 7 2 1 0 : -14.963235240431468 -14.747043227928751 -14.7310801427521 -14.72393649015475 -14.719760474646746 -14.717062283118723 -14.715209944939291 -14.713878545181595 -14.712885094089964 -14.7121203242037 -14.711515954138287
5 # 8 0 0 2 : -11.977612003809817 -14.732915927764871 -15.188483954798421 -15.449012850110863 -15.623827864459615 -15.749146729572077 -15.843056618819753 -15.915913838909626 -15.974052991959747 -16.021532312789294 -16.061054363621178
6 # 7 1 1 1 : -20.028998586521087 -16.887277863688375 -16.740660023173227 -16.675136346599757 -16.63781881676098 -16.613841057268658 -16.597168335428375 -16.584903629851304 -16.575496115598586 -16.56804544714508 -16.561994415341395
7 # 6 2 2 0 : -16.61159694113873 -16.452788458890197 -16.46295253777152 -16.470352584816695 -16.47547180143134 -16.47919567998401 -16.48203029297254 -16.484263936430125 -16.48607102265639 -16.487563508299964 -16.48881688903326
8 # 7 2 0 1 : -15.887943076733494 -15.790564574903685 -15.810158143318768 -15.819644969186045 -15.825524608533778 -15.829491022082886 -15.832300786502069 -15.83436780035172 -15.835937561570242 -15.837162636186859 -15.838141332667961
9 # 6 3 1 0 : -16.081754803088685 -15.872187968111923 -15.861472691945579 -15.856983590805868 -15.854407868732828 -15.852760010242491 -15.851637836727662 -15.850836933005027 -15.850242873291222 -15.849787746373472 -15.849429401006539
10 # 7 1 0 2 : -16.195666087307426 -16.677867227255252 -16.921805125300562 -17.066573313972576 -17.164760174880623 -17.235200004698964 -17.2878198973409 -17.328459790892357 -17.360734028163566 -17.38696733049872 -17.40870794047196
11 # 6 2 1 1 : -21.456172204482055 -18.125128438935814 -17.945136829104282 -17.862583221944416 -17.81466930664732 -17.78361367179647 -17.76197624497887 -17.746084468733738 -17.73393366006567 -17.724346275296295 -17.716589235230014
12 # 5 3 2 0 : -17.572666204531814 -17.49331284851714 -17.532768846469768 -17.55647753192311 -17.572084768274923 -17.583034206002633 -17.59109750264769 -17.597266236290537 -17.602131533179627 -17.60606450517241 -17.609308614474667
13 # 7 0 0 3 : -13.018210911162125 -16.439560184226814 -17.07014568495005 -17.434065736430135 -17.677742469758805 -17.850801111864595 -17.978836031679858 -18.076796271319253 -18.153897800377976 -18.216044081487574 -18.267147796753903
14 # 6 1 1 2 : -21.540771185708525 -18.63035111750182 -18.62526277273989 -18.642369704263626 -18.66000389814369 -18.674935318048522 -18.687085705223726 -18.69695846852554 -18.70506279747467 -18.711804412715235 -18.71748755638855
15 # 5 2 2 1 : -22.76107629399832 -19.30160904663659 -19.108791463522692 -19.019383456903228 -18.967077286803956 -18.933116205738422 -18.909499767975987 -18.892215205948887 -18.87905160939128 -18.868705077738902 -18.860363442647326
16 # 4 3 3 0 : -18.346267482843494 -18.432485480278213 -18.53008208478288 -18.586785116012766 -18.623864036215974 -18.649672756001678 -18.668494099700652 -18.68274883436975 -18.693886140137774 -18.70281399453071 -18.71012450577138
17 # 6 3 0 1 : -16.61750913212031 -16.472290282073644 -16.478928225846403 -16.482178984820933 -16.484183247286737 -16.485534338113 -16.486492763458507 -16.487199168803627 -16.48773660007381 -16.48815665950307 -16.488492650631834
18 # 5 4 1 0 : -16.700771408084638 -16.50452986698806 -16.49865528797123 -16.49639408913242 -16.49512968889967 -16.49433170128562 -16.4937936082825 -16.493412461432815 -16.493131355758457 -16.49291687814348 -16.492748489177256
19 # 6 2 0 2 : -17.37972307225597 -17.66417440842789 -17.828452479147863 -17.9244973085897 -17.98870803047759 -18.03427643223965 -18.06803374162208 -18.093931601975953 -18.1143848896453 -18.130932054733343 -18.144590025502723
20 # 5 3 1 1 : -22.21900378499455 -18.802166584266587 -18.61051963631457 -18.522437703633795 -18.471122886051795 -18.43780390211083 -18.41458592315717 -18.3975456481245 -18.38453095552364 -18.374274096630202 -18.36598485892285
21 # 4 4 2 0 : -18.003623775886712 -18.021804951556042 -18.094388766045686 -18.136684487449745 -18.164430316253597 -18.18384123799425 -18.19808045665593 -18.208928384563524 -18.217450487372762 -18.224315817270867 -18.229962266708554
22 # 6 1 0 3 : -17.16049200049263 -18.02206601357919 -18.39871681227843 -18.620152523674896 -18.768486714158428 -18.873367166120964 -18.95055758299072 -19.009338871778475 -19.055424145755403 -19.092455400834993 -19.122832769370497
23 # 5 2 1 2 : -22.80757560619405 -19.59722396980184 -19.50403999089238 -19.4707112798977 -19.454901033288287 -19.446307321999885 -19.441177695148593 -19.437890853968625 -19.435664943064204 -19.43408984395247 -19.43293517156931
24 # 4 3 2 1 : -23.35368025584886 -19.89176645110762 -19.714521800563404 -19.63393660823006 -19.587256127129596 -19.55717871072963 -19.536398940375854 -19.52127497123549 -19.509810914464268 -19.50083575341216 -19.493623738030568
25 # 3 4 3 0 : -18.544691873311553 -18.83245044498279 -18.998871070014136 -19.09507772969535 -19.158231115081733 -19.202254607537025 -19.234344175243976 -19.258615121061002 -19.277548431107952 -19.29270314690987 -19.30509636892189
26 # 6 0 0 4 : -13.610421560087484 -17.490059925747776 -18.22768663909486 -18.644415733891467 -18.91821408599588 -19.109206286417045 -19.248252465751495 -19.353179309127565 -19.434817794058194 -19.49999819286121 -19.553178850829017
27 # 5 1 1 3 : -22.35157162896833 -19.60629408763316 -19.692556368127498 -19.76240095172593 -19.814399535187484 -19.852991254358802 -19.882149773759583 -19.904713911547237 -19.922597148124712 -19.937080636534112 -19.94903391958495
28 # 4 2 2 2 : -23.913864999276612 -20.482872693368662 -20.323468864208976 -20.25140580154189 -20.209913755729303 -20.183385472420483 -20.16519556468682 -20.152043440422506 -20.14212853546848 -20.134400937505827 -20.128214131775422
29 # 3 3 3 1 : -23.717841725787117 -20.35808469447123 -20.232426120349174 -20.180958815223494 -20.153122172887105 -20.136070784414716 -20.12473920574103 -20.116740274295957 -20.11082413595477 -20.106284411288847 -20.102696717501395
30 # 2 4 4 0 : -18.436760987072333 -19.06207365640144 -19.345388470240497 -19.50902968846037 -19.616786134718662 -19.691924251453365 -19.746613372838777 -19.787890776566773 -19.820022201339786 -19.845692068840023 -19.86665027380036
31 # 5 4 0 1 : -16.87662364067997 -16.700767586591173 -16.700135547384853 -16.70015559520898 -16.70022543880007 -16.700295428426795 -16.700356978057226 -16.70040898229698 -16.700452366184262 -16.70048851953238 -16.700518780001172
32 # 4 5 1 0 : -16.87662364067997 -16.700767586591173 -16.700135547384853 -16.70015559520898 -16.70022543880007 -16.700295428426795 -16.700356978057226 -16.70040898229698 -16.700452366184262 -16.70048851953238 -16.700518780001172
33 # 5 3 0 2 : -17.946156477418963 -18.082803982954317 -18.19513545658477 -18.26027440064766 -18.303285247539204 -18.3335262605702 -18.355779084625294 -18.37276525944885 -18.386127493952813 -18.39690303196028 -18.405773037124824
34 # 4 4 1 1 : -22.462235222222485 -19.019700908140507 -18.824969672931708 -18.735504206302053 -18.68335287927473 -18.649481116803926 -18.625878801251403 -18.608560178217417 -18.59533667343495 -18.584918335647643 -18.576500941883044
35 # 3 5 2 0 : -17.946156477418963 -18.082803982954317 -18.19513545658477 -18.26027440064766 -18.303285247539204 -18.3335262605702 -18.355779084625294 -18.37276525944885 -18.386127493952813 -18.39690303196028 -18.405773037124824
36 # 5 2 0 3 : -18.18086338234885 -18.716914723518098 -18.970164950749574 -19.117380437892248 -19.21486949776466 -19.28322441117808 -19.333229487686356 -19.371138875114912 -19.400759068469508 -19.424496511135388 -19.44392716040802
37 # 4 3 1 2 : -23.368721992092517 -19.98950362908256 -19.84466642770449 -19.782188953236783 -19.747228996194146 -19.725261622517174 -19.710369446959433 -19.699688734968717 -19.691686840601058 -19.68548149177768 -19.680534298362765
38 # 3 4 2 1 : -23.368721992092517 -19.98950362908256 -19.84466642770449 -19.782188953236783 -19.747228996194146 -19.725261622517174 -19.710369446959433 -19.699688734968717 -19.691686840601058 -19.68548149177768 -19.680534298362765
39 # 2 5 3 0 : -18.18086338234885 -18.716914723518098 -18.970164950749574 -19.117380437892248 -19.21486949776466 -19.28322441117808 -19.333229487686356 -19.371138875114912 -19.400759068469508 -19.424496511135388 -19.44392716040802
40 # 5 1 0 4 : -17.618420698138113 -18.691190071424533 -19.135827040706022 -19.392743055658602 -19.562482636640006 -19.68104809089733 -19.767400797496773 -19.83258433019829 -19.883320925015436 -19.923849395508064 -19.956935056767442
41 # 4 2 1 3 : -23.4155144538949 -20.279234306001747 -20.23541574229725 -20.230730479251182 -20.233589349846117 -20.23787831210098 -20.24204396026545 -20.24572327500804 -20.248888387292197 -20.251599041249023 -20.25392872512635
42 # 3 3 2 2 : -24.249010805492144 -20.756456022857574 -20.57959386626736 -20.497226341156114 -20.44883808647228 -20.417474750640682 -20.395761102497183 -20.37994909597759 -20.367963679569783 -20.35858163616944 -20.351043503167315
43 # 2 4 3 1 : -23.4155144538949 -20.279234306001747 -20.23541574229725 -20.230730479251182 -20.233589349846117 -20.23787831210098 -20.24204396026545 -20.24572327500804 -20.248888387292197 -20.251599041249023 -20.25392872512635
44 # 1 5 4 0 : -17.618420698138113 -18.691190071424533 -19.135827040706022 -19.392743055658602 -19.562482636640006 -19.68104809089733 -19.767400797496773 -19.83258433019829 -19.883320925015436 -19.923849395508064 -19.956935056767442
45 # 5 0 0 5 : -13.803163956583148 -17.850055903907112 -18.62326077559763 -19.054322536769696 -19.3347503057033 -19.528849755215887 -19.669315119079386 -19.774838941551472 -19.856671097581263 -19.921848301037517 -19.97493235995905
46 # 4 1 1 4 : -22.610470101634526 -19.92546711497649 -20.043107708127796 -20.130479746912428 -20.19358547337479 -20.239630315527524 -20.2740251800213 -20.300420883428227 -20.321208642844752 -20.337960998908127 -20.35173178783654
47 # 3 2 2 3 : -24.263860750674873 -20.848534823820895 -20.704554218464114 -20.6412666276238 -20.605427613616776 -20.5827933177372 -20.567420828834905 -20.55638911664577 -20.548122326891686 -20.541710170870203 -20.536596573812098
48 # 2 3 3 2 : -24.263860750674873 -20.848534823820895 -20.704554218464114 -20.6412666276238 -20.605427613616776 -20.5827933177372 -20.567420828834905 -20.55638911664577 -20.548122326891686 -20.541710170870203 -20.536596573812098
49 # 1 4 4 1 : -22.610470101634526 -19.92546711497649 -20.043107708127796 -20.130479746912428 -20.19358547337479 -20.239630315527524 -20.2740251800213 -20.300420883428227 -20.321208642844752 -20.337960998908127 -20.35173178783654
50 # 0 5 5 0 : -13.803163956583148 -17.850055903907112 -18.62326077559763 -19.054322536769696 -19.3347503057033 -19.528849755215887 -19.669315119079386 -19.774838941551472 -19.856671097581263 -19.921848301037517 -19.97493235995905

 

$ lkgen -lk lk_n10_th0.001_rh11_10.lk -nseq 10

더보기

$ cat new_lk_n10_th0.001_rh11_10.lk

10 50
1  0.001000
11 10.000000


Type    00  01  10  11 Rho    0.0     1.0     2.0     3.0     4.0     5.0     6.0     7.0     8.0     9.0    10.0

   1 #   9   0   0   1  :   -17.67  -18.12  -18.37  -18.55  -18.67  -18.77  -18.84  -18.90  -18.95  -19.00  -19.03
   2 #   8   1   1   0  :   -19.71  -19.69  -19.68  -19.67  -19.66  -19.65  -19.65  -19.64  -19.64  -19.64  -19.64
   3 #   8   1   0   1  :   -21.78  -21.53  -21.50  -21.50  -21.52  -21.54  -21.55  -21.56  -21.58  -21.59  -21.59
   4 #   7   2   1   0  :   -21.79  -21.76  -21.75  -21.73  -21.72  -21.72  -21.71  -21.71  -21.70  -21.70  -21.70
   5 #   8   0   0   2  :   -19.38  -20.01  -20.41  -20.71  -20.93  -21.11  -21.26  -21.39  -21.50  -21.60  -21.68
   6 #   7   1   1   1  :   -30.38  -25.84  -25.23  -24.94  -24.76  -24.63  -24.54  -24.47  -24.41  -24.37  -24.33
   7 #   6   2   2   0  :   -23.73  -23.70  -23.69  -23.68  -23.67  -23.67  -23.67  -23.67  -23.67  -23.67  -23.67
   8 #   7   2   0   1  :   -23.04  -22.91  -22.88  -22.88  -22.88  -22.88  -22.89  -22.89  -22.90  -22.90  -22.90
   9 #   6   3   1   0  :   -23.04  -23.01  -22.99  -22.98  -22.97  -22.96  -22.95  -22.95  -22.95  -22.95  -22.94
  10 #   7   1   0   2  :   -23.72  -23.35  -23.37  -23.45  -23.53  -23.60  -23.67  -23.73  -23.79  -23.84  -23.89
  11 #   6   2   1   1  :   -31.92  -27.37  -26.73  -26.39  -26.18  -26.04  -25.93  -25.84  -25.77  -25.71  -25.66
  12 #   5   3   2   0  :   -24.83  -24.79  -24.78  -24.78  -24.78  -24.78  -24.79  -24.80  -24.80  -24.81  -24.81
  13 #   7   0   0   3  :   -20.46  -21.20  -21.71  -22.09  -22.40  -22.65  -22.86  -23.05  -23.21  -23.35  -23.48
  14 #   6   1   1   2  :   -32.22  -27.26  -26.69  -26.45  -26.32  -26.24  -26.20  -26.17  -26.15  -26.13  -26.12
  15 #   5   2   2   1  :   -33.43  -28.76  -28.09  -27.73  -27.51  -27.35  -27.22  -27.13  -27.05  -26.99  -26.93
  16 #   4   3   3   0  :   -25.75  -25.69  -25.69  -25.70  -25.72  -25.74  -25.76  -25.78  -25.80  -25.82  -25.84
  17 #   6   3   0   1  :   -23.73  -23.65  -23.63  -23.62  -23.61  -23.61  -23.61  -23.61  -23.61  -23.61  -23.62
  18 #   5   4   1   0  :   -23.73  -23.69  -23.67  -23.66  -23.65  -23.64  -23.64  -23.64  -23.63  -23.63  -23.63
  19 #   6   2   0   2  :   -24.83  -24.66  -24.66  -24.70  -24.74  -24.79  -24.83  -24.87  -24.91  -24.94  -24.97
  20 #   5   3   1   1  :   -32.71  -28.16  -27.50  -27.16  -26.94  -26.78  -26.66  -26.57  -26.50  -26.43  -26.38
  21 #   4   4   2   0  :   -25.34  -25.28  -25.28  -25.28  -25.30  -25.31  -25.32  -25.34  -25.35  -25.37  -25.38
  22 #   6   1   0   3  :   -24.82  -24.38  -24.44  -24.57  -24.70  -24.82  -24.93  -25.04  -25.13  -25.21  -25.29
  23 #   5   2   1   2  :   -33.58  -28.71  -28.07  -27.77  -27.59  -27.47  -27.38  -27.32  -27.27  -27.23  -27.20
  24 #   4   3   2   1  :   -34.10  -29.35  -28.68  -28.33  -28.11  -27.95  -27.83  -27.74  -27.66  -27.60  -27.55
  25 #   3   4   3   0  :   -26.03  -25.93  -25.94  -25.98  -26.02  -26.07  -26.11  -26.15  -26.18  -26.22  -26.25
  26 #   6   0   0   4  :   -21.06  -21.87  -22.44  -22.89  -23.25  -23.55  -23.81  -24.04  -24.23  -24.41  -24.57
  27 #   5   1   1   3  :   -33.15  -27.99  -27.45  -27.24  -27.15  -27.10  -27.08  -27.08  -27.08  -27.09  -27.11
  28 #   4   2   2   2  :   -34.88  -30.00  -29.31  -28.95  -28.73  -28.57  -28.45  -28.36  -28.29  -28.23  -28.18
  29 #   3   3   3   1  :   -34.56  -29.68  -29.02  -28.69  -28.49  -28.35  -28.25  -28.17  -28.11  -28.06  -28.02
  30 #   2   4   4   0  :   -26.03  -25.83  -25.86  -25.94  -26.03  -26.11  -26.19  -26.26  -26.33  -26.39  -26.45
  31 #   5   4   0   1  :   -23.95  -23.90  -23.88  -23.87  -23.86  -23.85  -23.85  -23.85  -23.85  -23.85  -23.85
  32 #   4   5   1   0  :   -23.95  -23.90  -23.88  -23.87  -23.86  -23.85  -23.85  -23.85  -23.85  -23.85  -23.85
  33 #   5   3   0   2  :   -25.34  -25.25  -25.24  -25.26  -25.29  -25.32  -25.34  -25.37  -25.39  -25.41  -25.43
  34 #   4   4   1   1  :   -32.96  -28.41  -27.75  -27.40  -27.18  -27.02  -26.90  -26.81  -26.73  -26.67  -26.61
  35 #   3   5   2   0  :   -25.34  -25.25  -25.24  -25.26  -25.29  -25.32  -25.34  -25.37  -25.39  -25.41  -25.43
  36 #   5   2   0   3  :   -25.74  -25.55  -25.58  -25.64  -25.72  -25.79  -25.86  -25.93  -25.99  -26.04  -26.09
  37 #   4   3   1   2  :   -34.14  -29.34  -28.68  -28.34  -28.13  -27.99  -27.88  -27.80  -27.73  -27.68  -27.64
  38 #   3   4   2   1  :   -34.14  -29.34  -28.68  -28.34  -28.13  -27.99  -27.88  -27.80  -27.73  -27.68  -27.64
  39 #   2   5   3   0  :   -25.74  -25.55  -25.58  -25.64  -25.72  -25.79  -25.86  -25.93  -25.99  -26.04  -26.09
  40 #   5   1   0   4  :   -25.33  -24.86  -24.94  -25.10  -25.25  -25.40  -25.54  -25.67  -25.78  -25.88  -25.98
  41 #   4   2   1   3  :   -34.29  -29.28  -28.66  -28.37  -28.20  -28.10  -28.02  -27.97  -27.94  -27.91  -27.89
  42 #   3   3   2   2  :   -35.26  -30.37  -29.67  -29.30  -29.07  -28.90  -28.78  -28.68  -28.60  -28.53  -28.48
  43 #   2   4   3   1  :   -34.29  -29.28  -28.66  -28.37  -28.20  -28.10  -28.02  -27.97  -27.94  -27.91  -27.89
  44 #   1   5   4   0  :   -25.33  -24.86  -24.94  -25.10  -25.25  -25.40  -25.54  -25.67  -25.78  -25.88  -25.98
  45 #   5   0   0   5  :   -21.26  -22.09  -22.68  -23.15  -23.53  -23.85  -24.13  -24.37  -24.58  -24.77  -24.94
  46 #   4   1   1   4  :   -33.44  -28.22  -27.69  -27.49  -27.41  -27.38  -27.37  -27.37  -27.39  -27.40  -27.42
  47 #   3   2   2   3  :   -35.31  -30.35  -29.66  -29.31  -29.08  -28.93  -28.82  -28.73  -28.66  -28.60  -28.56
  48 #   2   3   3   2  :   -35.31  -30.35  -29.66  -29.31  -29.08  -28.93  -28.82  -28.73  -28.66  -28.60  -28.56
  49 #   1   4   4   1  :   -33.44  -28.22  -27.69  -27.49  -27.41  -27.38  -27.37  -27.37  -27.39  -27.40  -27.42
  50 #   0   5   5   0  :   -21.26  -22.09  -22.68  -23.15  -23.53  -23.85  -24.13  -24.37  -24.58  -24.77  -24.94

 

CODE : ldtable.py

더보기
#! /usr/bin/env python
from __future__ import print_function
from ldpop import LookupTable, rhos_from_string
import argparse
import logging
import multiprocessing as mp
import math
import subprocess


if __name__ == '__main__':
    parser = argparse.ArgumentParser(
        description='Print a lookup table, as expected by ldhat or ldhelmet.')
    parser.add_argument('-n', type=int, metavar='N', required=True,
                        help='sample size')
    parser.add_argument('-th', type=float, metavar='theta',
                        help='twice the mutation rate')
    parser.add_argument('-rh', type=str,
                        required=False, metavar='num_rh,max_rh',
                        help='grid of rhos (twice the recomb rate). The grid '
                             'has num_rh uniformly spaced points from 0 to '
                             'max_rh, inclusive. (((Alternatively, to create '
                             'a non-uniform grid, use '
                             '"-rh r0,step0,r1,step1,r2,...rK". '
                             'This creates a grid '
                             '{r0,r0+step0,r0+2*step0,...,r1,r1+step1,...,rK} '
                             'similar to ldhelmet. Note that non-uniform grid '
                             'is incompatible with vanilla ldhat.)))',
                        default='11,10')
    parser.add_argument('-s', type=str, metavar='s0,s1,...,sD',
                        help='coalescent scaled population sizes '
                             '(s0=present size, sD=ancient size)')
    parser.add_argument('-t', type=str, metavar='t1,...,tD',
                        help='times of size changes from present backwards. '
                             'Must be increasing positive reals.')
    parser.add_argument('--approx', action='store_true',
                        help='use finite moran model. A reasonable '
                             'approximation that is much faster than '
                             'the exact formula. Accuracy of the  '
                             'approximation can be improved by taking n '
                             'larger than needed, and using ldhat/lkgen.c '
                             'to subsample. (Converges to the "exact" '
                             'diffusion as n->infty)')
    parser.add_argument('-c','--cores', type=int, default=mp.cpu_count(),
                        help='Number of parallel processes to use.')
    parser.add_argument('--log', type=str, metavar='logfile',
                        help='Log extra info to logfile. '
                             'If logfile=".", logs to STDERR.')

    parser.add_argument('-L', type=int, metavar='total_length',
                        help='Total length of DNA sequence')

    parser.add_argument('-S', type=int, metavar='seg_sites', 
                        help='Segregating sites')

    parser.add_argument('-add',type=str,metavar='suffix')
    
    args = parser.parse_args()

    if args.log == '.':
        logging.basicConfig(level=logging.INFO)
    elif args.log is not None:
        logging.basicConfig(filename=args.log, level=logging.INFO)

    rho='_'.join((args.rh).split(','))
    
    #theta=0
    if args.L and args.S:
        for i in range(args.n-1):
            theta = theta + 1/(i+1)
        theta = theta*math.log(args.L/(args.L-args.S))
    else:
        theta=args.th
        
    exact = not args.approx

    apprx=''
    if exact == False:
        apprx='a_'

    numCores = args.cores
    
    rhos = rhos_from_string(args.rh)
    rho='_'.join((args.rh).split(','))
    
    popSizes, times = args.s, args.t
    assert (popSizes is None) == (times is None)
    if popSizes is None:
        popSizes = [1.]
        times = []

        pt=''
    else:
        popSizes = [float(_) for _ in popSizes.split(',')]
        times = [float(_) for _ in times.split(',')]

        pt='_s'+str('_'.join(map(str,popSizes)))+'_t'+str('_'.join(map(str,times)))
        
    assert len(popSizes) == len(times)+1
    
    lk=LookupTable(args.n, theta, rhos, popSizes, times, exact, numCores)

    #print(args.n, theta, rhos, popSizes, times, exact, numCores)

    filename=apprx+'n'+str(args.n)+'_th'+str(round(theta,8))+'_rh'+str(rho)+pt+'.lk'
    
    with open(filename,'w') as f:
        f.write(str(lk)+'\n')

    cmd='/ssd/program/LDhat/lkgen -lk '+filename+' -nseq '+str(args.n)
    subprocess.run(cmd ,shell=True)
ldtable.py
0.00MB

 

'Tools' 카테고리의 다른 글

checkVCF  (0) 2020.09.15
Windows Terminal  (0) 2020.08.17
Slurm - Workload manager  (0) 2020.08.15
GATK  (0) 2020.08.07
Shapeit4  (0) 2020.07.09

댓글