Commit 05e317e2 authored by Ryan Gutenkunst's avatar Ryan Gutenkunst
Browse files

Add CUDA benchmarking code

parent 1419cf50
......@@ -253,8 +253,7 @@ def one_pop(phi, xx, T, nu=1, gamma=0, h=0.5, theta0=1.0, initial_t=0,
def two_pops(phi, xx, T, nu1=1, nu2=1, m12=0, m21=0, gamma1=0, gamma2=0,
h1=0.5, h2=0.5, theta0=1, initial_t=0, frozen1=False,
frozen2=False, nomut1=False, nomut2=False, enable_cuda_const=False,
transpose=False):
frozen2=False, nomut1=False, nomut2=False, enable_cuda_cached=False):
"""
Integrate a 2-dimensional phi foward.
......@@ -281,7 +280,7 @@ def two_pops(phi, xx, T, nu1=1, nu2=1, m12=0, m21=0, gamma1=0, gamma2=0,
nomut1,nomut2: If True, no new mutations will be introduced into the
given population.
enable_cuda_const: If True, enable CUDA integration with slower constant
enable_cuda_cached: If True, enable CUDA integration with slower constant
parameter method. Likely useful only for benchmarking.
Note: Generalizing to different grids in different phi directions is
......@@ -305,7 +304,7 @@ def two_pops(phi, xx, T, nu1=1, nu2=1, m12=0, m21=0, gamma1=0, gamma2=0,
if False and numpy.all([numpy.isscalar(var) for var in vars_to_check]):
# Constant integration with CUDA turns out to be slower,
# so we only use it in specific circumsances.
if not cuda_enabled or (cuda_enabled and enable_cuda_const):
if not cuda_enabled or (cuda_enabled and enable_cuda_cached):
return _two_pops_const_params(phi, xx, T, nu1, nu2, m12, m21,
gamma1, gamma2, h1, h2, theta0, initial_t,
frozen1, frozen2, nomut1, nomut2)
......@@ -356,22 +355,12 @@ def two_pops(phi, xx, T, nu1=1, nu2=1, m12=0, m21=0, gamma1=0, gamma2=0,
_inject_mutations_2D(phi, this_dt, xx, yy, theta0, frozen1, frozen2,
nomut1, nomut2)
if not transpose:
if not frozen1:
phi = int_c.implicit_2Dx(phi, xx, yy, nu1, m12, gamma1, h1,
this_dt, use_delj_trick)
if not frozen2:
phi = int_c.implicit_2Dy(phi, xx, yy, nu2, m21, gamma2, h2,
this_dt, use_delj_trick)
else:
phi = phi.transpose()
if not frozen1:
phi = int_c.implicit_2Dy(phi, xx, yy, nu1, m12, gamma1, h1,
this_dt, use_delj_trick)
phi = phi.transpose()
if not frozen2:
phi = int_c.implicit_2Dy(phi, xx, yy, nu2, m21, gamma2, h2,
this_dt, use_delj_trick)
if not frozen1:
phi = int_c.implicit_2Dx(phi, xx, yy, nu1, m12, gamma1, h1,
this_dt, use_delj_trick)
if not frozen2:
phi = int_c.implicit_2Dy(phi, xx, yy, nu2, m21, gamma2, h2,
this_dt, use_delj_trick)
current_t = next_t
return phi
......@@ -380,7 +369,7 @@ def three_pops(phi, xx, T, nu1=1, nu2=1, nu3=1,
m12=0, m13=0, m21=0, m23=0, m31=0, m32=0,
gamma1=0, gamma2=0, gamma3=0, h1=0.5, h2=0.5, h3=0.5,
theta0=1, initial_t=0, frozen1=False, frozen2=False,
frozen3=False, enable_cuda_const=False, transpose=False):
frozen3=False, enable_cuda_cached=False):
"""
Integrate a 3-dimensional phi foward.
......@@ -400,7 +389,7 @@ def three_pops(phi, xx, T, nu1=1, nu2=1, nu3=1,
initial_t: Time at which to start integration. (Note that this only matters
if one of the demographic parameters is a function of time.)
enable_cuda_const: If True, enable CUDA integration with slower constant
enable_cuda_cached: If True, enable CUDA integration with slower constant
parameter method. Likely useful only for benchmarking.
Note: Generalizing to different grids in different phi directions is
......@@ -426,7 +415,7 @@ def three_pops(phi, xx, T, nu1=1, nu2=1, nu3=1,
vars_to_check = [nu1,nu2,nu3,m12,m13,m21,m23,m31,m32,gamma1,gamma2,
gamma3,h1,h2,h3,theta0]
if False and numpy.all([numpy.isscalar(var) for var in vars_to_check]):
if not cuda_enabled or (cuda_enabled and enable_cuda_const):
if not cuda_enabled or (cuda_enabled and enable_cuda_cached):
return _three_pops_const_params(phi, xx, T, nu1, nu2, nu3,
m12, m13, m21, m23, m31, m32,
gamma1, gamma2, gamma3, h1, h2, h3,
......@@ -495,29 +484,15 @@ def three_pops(phi, xx, T, nu1=1, nu2=1, nu3=1,
_inject_mutations_3D(phi, this_dt, xx, yy, zz, theta0,
frozen1, frozen2, frozen3)
if not transpose:
if not frozen1:
phi = int_c.implicit_3Dx(phi, xx, yy, zz, nu1, m12, m13,
gamma1, h1, this_dt, use_delj_trick)
if not frozen2:
phi = int_c.implicit_3Dy(phi, xx, yy, zz, nu2, m21, m23,
gamma2, h2, this_dt, use_delj_trick)
if not frozen3:
phi = int_c.implicit_3Dz(phi, xx, yy, zz, nu3, m31, m32,
gamma3, h3, this_dt, use_delj_trick)
else:
phi = phi.transpose(1,2,0)
if not frozen1:
phi = int_c.implicit_3Dz(phi, xx, yy, zz, nu1, m12, m13,
gamma1, h1, this_dt, use_delj_trick)
phi = phi.transpose(2,1,0)
if not frozen2:
phi = int_c.implicit_3Dz(phi, xx, yy, zz, nu2, m21, m23,
gamma2, h2, this_dt, use_delj_trick)
phi = phi.transpose(0,2,1)
if not frozen3:
phi = int_c.implicit_3Dz(phi, xx, yy, zz, nu3, m31, m32,
gamma3, h3, this_dt, use_delj_trick)
if not frozen1:
phi = int_c.implicit_3Dx(phi, xx, yy, zz, nu1, m12, m13,
gamma1, h1, this_dt, use_delj_trick)
if not frozen2:
phi = int_c.implicit_3Dy(phi, xx, yy, zz, nu2, m21, m23,
gamma2, h2, this_dt, use_delj_trick)
if not frozen3:
phi = int_c.implicit_3Dz(phi, xx, yy, zz, nu3, m31, m32,
gamma3, h3, this_dt, use_delj_trick)
current_t = next_t
return phi
......@@ -528,7 +503,7 @@ def four_pops(phi, xx, T, nu1=1, nu2=1, nu3=1, nu4=1,
gamma1=0, gamma2=0, gamma3=0, gamma4=0,
h1=0.5, h2=0.5, h3=0.5, h4=0.5,
theta0=1, initial_t=0, frozen1=False, frozen2=False,
frozen3=False, frozen4=False, enable_cuda_const=False):
frozen3=False, frozen4=False):
"""
Integrate a 4-dimensional phi foward.
......@@ -653,8 +628,7 @@ def five_pops(phi, xx, T, nu1=1, nu2=1, nu3=1, nu4=1, nu5=1,
gamma1=0, gamma2=0, gamma3=0, gamma4=0, gamma5=0,
h1=0.5, h2=0.5, h3=0.5, h4=0.5, h5=0.5,
theta0=1, initial_t=0, frozen1=False, frozen2=False,
frozen3=False, frozen4=False, frozen5=False,
enable_cuda_const=False):
frozen3=False, frozen4=False, frozen5=False):
"""
Integrate a 5-dimensional phi foward.
......@@ -674,9 +648,6 @@ def five_pops(phi, xx, T, nu1=1, nu2=1, nu3=1, nu4=1, nu5=1,
initial_t: Time at which to start integration. (Note that this only matters
if one of the demographic parameters is a function of time.)
enable_cuda_const: If True, enable CUDA integration with slower constant
parameter method. Likely useful only for benchmarking.
Note: Generalizing to different grids in different phi directions is
straightforward. The tricky part will be later doing the extrapolation
correctly.
......
......@@ -45,4 +45,22 @@ def cuda_enabled(toggle=None):
Integration.cuda_enabled = False
return False
else:
raise ValueError("toggle must be True, False, or None")
\ No newline at end of file
raise ValueError("toggle must be True, False, or None")
def pts_to_RAM(pts, P):
"""
Approximate RAM usage for a given grid points and number of populations
pts: Grid points setting
P: Number of populations
"""
return 8*4*pts**P / 1024**3
def RAM_to_pts(RAM, P):
"""
Approximate maximum grid points given the number of populations and available RAM
pts: Grid points setting
P: Number of populations
"""
return int((RAM*1024**3/(8*4))**(1./P))
\ No newline at end of file
#!/bin/env python
import socket, sys
import numpy as np
import dadi
import models
import argparse
parser = argparse.ArgumentParser(description='Benchmark dadi integration')
parser.add_argument('--cuda', action='store_true', help='run using GPU')
parser.add_argument('--RAM', type=float, default=0.01, help='approximate amount of RAM to use (in GB)')
args = parser.parse_args()
if args.cuda:
if not dadi.cuda_enabled(True):
raise ValueError('Failed to initialize CUDA')
print('# Host: {0}'.format(socket.getfqdn()))
print('# CUDA enabled: {0}'.format(dadi.cuda_enabled()))
if dadi.cuda_enabled():
import pycuda.autoinit
device = pycuda.autoinit.device
print('# GPU: {0} ({1:.1f} GB RAM)'.format(device.name(),
device.total_memory()/1024**3))
print()
# Burn in to ensure GPU kernels are loaded
for ii in range(3):
models.OutOfAfrica_2L06(10, variant='original')
models.OutOfAfrica_2L06(10, variant='const_int')
maxpts = int(dadi.RAM_to_pts(args.RAM, 2))
pts_l_2D = 10 * 2**np.arange(0,np.log2(maxpts//10), dtype=int)
print('# Model: 2D Original')
for pts in pts_l_2D:
_,t = models.OutOfAfrica_2L06_timed(pts, variant='original')
print('pts: {0}, time: {1:.4g}'.format(pts, t), flush=True)
print()
print('# Model: 2D Cached integration')
for pts in pts_l_2D[:-1]:
_,t = models.OutOfAfrica_2L06_timed(pts, variant='cached_int')
print('pts: {0}, time: {1:.4g}'.format(pts, t), flush=True)
print()
maxpts = int(dadi.RAM_to_pts(args.RAM, 3))
pts_l_3D = 10 * 2**np.arange(0,np.log2(maxpts//10), dtype=int)
print('# Model: 3D Original')
for pts in pts_l_3D:
_, t = models.OutOfAfrica_3G09_timed(pts, variant='original')
print('pts: {0}, time: {1:.4g}'.format(pts, t), flush=True)
print()
print('# Model: 3D Constant')
for pts in pts_l_3D:
_, t = models.OutOfAfrica_3G09_timed(pts, variant='const')
print('pts: {0}, time: {1:.4g}'.format(pts, t), flush=True)
print()
print('# Model: 3D Cached integration')
for pts in pts_l_3D[:-1]:
_, t = models.OutOfAfrica_3G09_timed(pts, variant='cached_int')
print('pts: {0}, time: {1:.4g}'.format(pts, t), flush=True)
print()
maxpts = int(dadi.RAM_to_pts(args.RAM, 4))
pts_l_4D = 10 * 2**np.arange(0,np.log2(maxpts//10), dtype=int)
print('# Model: 4D Original')
for pts in pts_l_4D:
_,t = models.NewWorld_4G09_timed(pts, variant='original')
print('pts: {0}, time: {1:.4g}'.format(pts, t), flush=True)
print()
maxpts = int(dadi.RAM_to_pts(args.RAM, 5))
pts_l_5D = 10 * 2**np.arange(0,np.log2(maxpts//10), dtype=int)
print('# Model: 5D Original')
for pts in pts_l_5D:
_, t = models.OutOfAfricaArchaicAdmixture_5R19_timed(pts, variant='original')
print('pts: {0}, time: {1:.4g}'.format(pts, t), flush=True)
print()
\ No newline at end of file
import math, time
import dadi
def timeit(func):
def newfunc(*args, **kwargs):
start = time.time()
out = func(*args, **kwargs)
end = time.time()
return out, end-start
return newfunc
def OutOfAfricaArchaicAdmixture_5R19(pts, variant='original'):
"""
Integration of model from Ragsdale (2019) PLoS Genetics.
"""
# Output
# Pop 1: African
# Pop 2: CEU
# Pop 3: CHB
# Parameter values from stdpopsim GitHub
# From here ***
# First we set out the maximum likelihood values of the various parameters
# given in Table 1 (under archaic admixture).
N_0 = 3600
N_YRI = 13900
N_B = 880
N_CEU0 = 2300
N_CHB0 = 650
# Times are provided in years, so we convert into generations.
# In the published model, the authors used a generation time of 29 years to
# convert from genetic to physical units
generation_time = 29
T_AF = 300e3 / generation_time
T_B = 60.7e3 / generation_time
T_EU_AS = 36.0e3 / generation_time
T_arch_afr_split = 499e3 / generation_time
T_arch_afr_mig = 125e3 / generation_time
T_nean_split = 559e3 / generation_time
T_arch_adm_end = 18.7e3 / generation_time
# We need to work out the starting (diploid) population sizes based on
# the growth rates provided for these two populations
r_CEU = 0.00125
r_CHB = 0.00372
N_CEU = N_CEU0 / math.exp(-r_CEU * T_EU_AS)
N_CHB = N_CHB0 / math.exp(-r_CHB * T_EU_AS)
# Migration rates during the various epochs.
m_AF_B = 52.2e-5
m_YRI_CEU = 2.48e-5
m_YRI_CHB = 0e-5
m_CEU_CHB = 11.3e-5
m_AF_arch_af = 1.98e-5
m_OOA_nean = 0.825e-5
# To here ***
xx = dadi.Numerics.default_grid(pts)
phi = dadi.PhiManip.phi_1D(xx)
# Split off Neanderthal pop
phi = dadi.PhiManip.phi_1D_to_2D(xx, phi)
phi = dadi.Integration.two_pops(phi, xx, (T_nean_split-T_arch_afr_split)/(2*N_0))
# Split off archaic African pop
phi = dadi.PhiManip.phi_2D_to_3D(phi, 1, xx,xx,xx)
phi = dadi.Integration.three_pops(phi, xx, (T_arch_afr_split-T_AF)/(2*N_0))
# African population growth
phi = dadi.Integration.three_pops(phi, xx, (T_AF-T_arch_afr_mig)/(2*N_0), nu1=N_YRI/N_0)
# Archaic African migration begins
phi = dadi.Integration.three_pops(phi, xx, (T_arch_afr_mig-T_B)/(2*N_0), nu1=N_YRI/N_0,
m13=m_AF_arch_af*2*N_0, m31=m_AF_arch_af*2*N_0)
# Split of Eurasian ancestral pop
phi = dadi.PhiManip.phi_3D_to_4D(phi, 1,0, xx,xx,xx,xx)
phi = dadi.Integration.four_pops(phi, xx, (T_B-T_EU_AS)/(2*N_0),
nu1=N_YRI/N_0, nu4=N_B/N_0,
m13=m_AF_arch_af*2*N_0, m31=m_AF_arch_af*2*N_0,
m14=m_AF_B*2*N_0, m41=m_AF_B*2*N_0,
m24=m_OOA_nean*2*N_0, m42=m_OOA_nean*2*N_0)
# Split of European and Asian ancestral pops
phi = dadi.PhiManip.phi_4D_to_5D(phi, 0,0,0, xx,xx,xx,xx,xx)
nuEu_func = lambda t: N_CEU0/N_0*(N_CEU/N_CEU0)**(t * 2*N_0/(T_EU_AS))
nuAs_func = lambda t: N_CHB0/N_0*(N_CHB/N_CHB0)**(t * 2*N_0/(T_EU_AS))
phi = dadi.Integration.five_pops(phi, xx, (T_EU_AS-T_arch_adm_end)/(2*N_0),
nu1=N_YRI/N_0, nu4=nuEu_func, nu5=nuAs_func,
m13=m_AF_arch_af*2*N_0, m31=m_AF_arch_af*2*N_0,
m14=m_YRI_CEU*2*N_0, m41=m_YRI_CEU*2*N_0,
m15=m_YRI_CHB*2*N_0, m51=m_YRI_CHB*2*N_0,
m24=m_OOA_nean*2*N_0, m42=m_OOA_nean*2*N_0,
m25=m_OOA_nean*2*N_0, m52=m_OOA_nean*2*N_0,
m45=m_CEU_CHB*2*N_0, m54=m_CEU_CHB*2*N_0,
)
# End of archaic migration. Remove archaic pops for efficiency.
phi = dadi.PhiManip.filter_pops(phi, xx, [1,4,5])
# We use initial_t argument so we can reuse nuEu_func and nuAs_func
phi = dadi.Integration.three_pops(phi, xx, T_arch_adm_end/(2*N_0), initial_t=(T_EU_AS-T_arch_adm_end)/(2*N_0),
nu1=N_YRI/N_0, nu2=nuEu_func, nu3=nuAs_func,
m12=m_YRI_CEU*2*N_0, m21=m_YRI_CEU*2*N_0,
m13=m_YRI_CHB*2*N_0, m31=m_YRI_CHB*2*N_0,
m23=m_CEU_CHB*2*N_0, m32=m_CEU_CHB*2*N_0,
)
return phi
OutOfAfricaArchaicAdmixture_5R19_timed = timeit(OutOfAfricaArchaicAdmixture_5R19)
def NewWorld_4G09(pts, variant='original'):
"""
Integration of New World model from Gutenkunst (2009) PLoS Genetics.
Note that in the original publication, the African population was removed late
in the simulation, to keep the total number of populations to three. Now
we can run with that fourth population.
"""
# From Table S3
nuAf, nuB = 1.68, 0.287
mAfB, mAfEu, mAfAs = 3.65, 0.44, 0.28
TAf, TB = 0.607, 0.396
# From Table S6
nuEu0, nuEu = 0.208, 2.42
nuAs0, nuAs = 0.081, 4.186
mEuAs = 1.98
TEuAs = 0.072
nuMx0, nuMx = 0.103, 7.94
TMx = 0.059
fEuMx = 0.48
xx = dadi.Numerics.default_grid(pts)
phi = dadi.PhiManip.phi_1D(xx)
phi = dadi.Integration.one_pop(phi, xx, TAf-TB-TEuAs-TMx, nu=nuAf)
phi = dadi.PhiManip.phi_1D_to_2D(xx, phi)
phi = dadi.Integration.two_pops(phi, xx, TB-TEuAs-TMx, nu1=nuAf, nu2=nuB,
m12=mAfB, m21=mAfB)
phi = dadi.PhiManip.phi_2D_to_3D_split_2(xx, phi)
nuEu_func = lambda t: nuEu0*(nuEu/nuEu0)**(t/(TEuAs+TMx))
nuAs_func = lambda t: nuAs0*(nuAs/nuAs0)**(t/(TEuAs+TMx))
phi = dadi.Integration.three_pops(phi, xx, TEuAs,
nu1=nuAf,
nu2=nuEu_func, nu3=nuAs_func,
m12=mAfEu, m21=mAfEu, m13=mAfAs, m31=mAfAs,
m23=mEuAs, m32=mEuAs)
# Split off New World population, with contribution entirely from
# population 3 (East Asia)
phi = dadi.PhiManip.phi_3D_to_4D(phi, 0,0, xx,xx,xx,xx)
nuEu0 = nuEu_func(TEuAs)
nuAs0 = nuAs_func(TEuAs)
nuEu_func = lambda t: nuEu0*(nuEu/nuEu0)**(t/TMx)
nuAs_func = lambda t: nuAs0*(nuAs/nuAs0)**(t/TMx)
nuMx_func = lambda t: nuMx0*(nuMx/nuMx0)**(t/TMx)
phi = dadi.Integration.four_pops(phi, xx, TMx,
nu1=nuAf, nu2=nuEu_func,
nu3=nuAs_func, nu4=nuMx_func,
m12=mAfEu, m21=mAfEu, m13=mAfAs, m31=mAfAs,
m23=mEuAs, m32=mEuAs)
phi = dadi.PhiManip.phi_4D_admix_into_4(phi, 0,fEuMx,0, xx,xx,xx,xx)
return phi
NewWorld_4G09_timed = timeit(NewWorld_4G09)
def OutOfAfrica_3G09(pts, variant='original'):
"""
Integration of Out-of-Africa model from Gutenkunst (2009) PLoS Genetics.
Potential variants:
'original': Model from the original publication
'const': Model in which European and Asian population have constant size, for testing
'cached_int': As in 'const', but forcing CUDA code to use cached matrices
"""
# Parameter values from https://github.com/popsim-consortium/stdpopsim/blob/master/stdpopsim/catalog/homo_sapiens.py
generation_time = 25
# First we set out the maximum likelihood values of the various parameters
# given in Table 1.
N_A = 7300
N_B = 2100
N_AF = 12300
N_EU0 = 1000
N_AS0 = 510
# Times are provided in years, so we convert into generations.
T_AF = 220e3 / generation_time
T_B = 140e3 / generation_time
T_EU_AS = 21.2e3 / generation_time
# We need to work out the starting (diploid) population sizes based on
# the growth rates provided for these two populations
r_EU = 0.004
r_AS = 0.0055
# Migration rates during the various epochs.
m_AF_B = 25e-5
m_AF_EU = 3e-5
m_AF_AS = 1.9e-5
m_EU_AS = 9.6e-5
xx = dadi.Numerics.default_grid(pts)
phi = dadi.PhiManip.phi_1D(xx)
phi = dadi.Integration.one_pop(phi, xx, (T_AF-T_B)/(2*N_A), nu=N_AF/N_A)
phi = dadi.PhiManip.phi_1D_to_2D(xx, phi)
phi = dadi.Integration.two_pops(phi, xx, (T_B-T_EU_AS)/(2*N_A), nu1=N_AF/N_A, nu2=N_B/N_A,
m12=m_AF_B * 2*N_A, m21=m_AF_B * 2*N_A)
phi = dadi.PhiManip.phi_2D_to_3D_split_2(xx, phi)
enable_cuda_cached = False
if variant=='original':
def nuEu_func(t): return N_EU0/N_A * math.exp(r_EU * t*2*N_A)
def nuAs_func(t): return N_AS0/N_A * math.exp(r_AS * t*2*N_A)
elif variant=='const':
nuEu_func = N_EU0/N_A
nuAs_func = N_AS0/N_A
elif variant=='cached_int':
nuEu_func = N_EU0/N_A
nuAs_func = N_AS0/N_A
enable_cuda_cached = True
phi = dadi.Integration.three_pops(phi, xx, T_EU_AS/(2*N_A), nu1=N_AF/N_A, nu2=nuEu_func, nu3=nuAs_func,
m12=m_AF_EU *2*N_A, m13=m_AF_AS *2*N_A, m21=m_AF_EU *2*N_A,
m23=m_EU_AS *2*N_A, m31=m_AF_AS *2*N_A, m32=m_EU_AS *2*N_A,
enable_cuda_cached=enable_cuda_cached)
return phi
OutOfAfrica_3G09_timed = timeit(OutOfAfrica_3G09)
def OutOfAfrica_2L06(pts, variant='original'):
"""
Integration of Out-of-Africa model from Li (2006) PLoS Genetics.
Potential variants:
'original': Model from the original publication
'cached_int': Forcing CUDA code to use cached matrices
"""
# Parameter values from https://github.com/popsim-consortium/stdpopsim/blob/master/stdpopsim/catalog/drosophila_melanogaster.py
# African Parameter values from "Demographic History of the African
# Population" section
N_A0 = 8.603e06
t_A0 = 600000 # assuming 10 generations / year
N_A1 = N_A0 / 5.0 # This is the ancestral size
# European Parameter values from "Demo History of Euro Population"
N_E0 = 1.075e06
N_E1 = 2200
t_AE = 158000 # generations
t_E1 = t_AE - 3400
# Ancestral size is N_A1
xx = dadi.Numerics.default_grid(pts)
phi = dadi.PhiManip.phi_1D(xx)
phi = dadi.Integration.one_pop(phi, xx, (t_A0 - t_AE)/(2*N_A1), nu=N_A0/N_A1)
phi = dadi.PhiManip.phi_1D_to_2D(xx, phi)
enable_cuda_cached = False
if variant == 'cached_int':
enable_cuda_cached = True
phi = dadi.Integration.two_pops(phi, xx, (t_AE - t_E1)/(2*N_A1), nu1=N_A0/N_A1, nu2=N_E1/N_A1,
enable_cuda_cached=enable_cuda_cached)
phi = dadi.Integration.two_pops(phi, xx, t_E1/(2*N_A1), nu1=N_A0/N_A1, nu2=N_E0/N_A1,
enable_cuda_cached=enable_cuda_cached)
return phi
OutOfAfrica_2L06_timed = timeit(OutOfAfrica_2L06)
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment