Menu
▾
▴
simupop-devel — Development dicussions
You can subscribe to this list here.
| 2011 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
(111) |
Jul
(77) |
Aug
(25) |
Sep
(12) |
Oct
(12) |
Nov
(24) |
Dec
(12) |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 2012 |
Jan
(7) |
Feb
(38) |
Mar
(152) |
Apr
(53) |
May
(4) |
Jun
(7) |
Jul
(13) |
Aug
(13) |
Sep
(51) |
Oct
(20) |
Nov
(4) |
Dec
(16) |
| 2013 |
Jan
(2) |
Feb
(19) |
Mar
(2) |
Apr
(12) |
May
(11) |
Jun
|
Jul
(1) |
Aug
|
Sep
|
Oct
|
Nov
|
Dec
(3) |
| 2014 |
Jan
(2) |
Feb
(46) |
Mar
(7) |
Apr
(8) |
May
(16) |
Jun
(25) |
Jul
(33) |
Aug
(4) |
Sep
(7) |
Oct
(23) |
Nov
(2) |
Dec
(1) |
| 2015 |
Jan
(2) |
Feb
(5) |
Mar
|
Apr
(17) |
May
(2) |
Jun
|
Jul
(1) |
Aug
(11) |
Sep
(1) |
Oct
(3) |
Nov
(1) |
Dec
|
| 2016 |
Jan
(11) |
Feb
(8) |
Mar
|
Apr
|
May
|
Jun
|
Jul
|
Aug
|
Sep
|
Oct
|
Nov
|
Dec
|
|
From: <bpe...@us...> - 2016-02-15 00:45:42
|
Revision: 5015
http://sourceforge.net/p/simupop/code/5015
Author: bpeng2000
Date: 2016-02-15 00:45:39 +0000 (Mon, 15 Feb 2016)
Log Message:
-----------
Add option -fPIC for linux system
Modified Paths:
--------------
trunk/setup.py
Modified: trunk/setup.py
===================================================================
--- trunk/setup.py 2016-02-10 22:30:25 UTC (rev 5014)
+++ trunk/setup.py 2016-02-15 00:45:39 UTC (rev 5015)
@@ -749,8 +749,10 @@
#distutils.log.set_verbosity(distutils.log.FATAL)
if os.name == 'nt':
NO_WARNING_ARG = [] # ['/nowarn'] is for vc 2015 only
+ SHLIB_ARG = []
else:
NO_WARNING_ARG = ['-w']
+ SHLIB_ARG = ['-fPIC']
try:
if not os.path.isfile('src/boost_pch.h.pch') or \
os.path.getmtime('src/boost_pch.h.pch') > os.path.getmtime('src/boost_pch.h'):
@@ -779,7 +781,7 @@
objects = c.compile(LIB_FILES,
include_dirs=['gsl', 'gsl/specfunc', 'build', '.', boost_include_dir] + common_extra_include_dirs,
output_dir='build',
- extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
+ extra_preargs = common_extra_compile_args + NO_WARNING_ARG + SHLIB_ARG,
macros = COMMON_MACROS + STD_MACROS
)
c.create_static_lib(objects, "simuPOP_shstd", output_dir='build')
@@ -787,7 +789,7 @@
objects = c.compile(LIB_FILES,
include_dirs=['gsl', 'gsl/specfunc', 'build', '.', boost_include_dir] + common_extra_include_dirs,
output_dir='build',
- extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
+ extra_preargs = common_extra_compile_args + NO_WARNING_ARG + SHLIB_ARG,
macros = COMMON_MACROS + OPT_MACROS
)
c.create_static_lib(objects, "simuPOP_shop", output_dir='build')
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-02-10 22:30:27
|
Revision: 5014
http://sourceforge.net/p/simupop/code/5014
Author: bpeng2000
Date: 2016-02-10 22:30:25 +0000 (Wed, 10 Feb 2016)
Log Message:
-----------
Fix setup.py for windows system
Modified Paths:
--------------
trunk/setup.py
Modified: trunk/setup.py
===================================================================
--- trunk/setup.py 2016-02-10 16:26:29 UTC (rev 5013)
+++ trunk/setup.py 2016-02-10 22:30:25 UTC (rev 5014)
@@ -187,9 +187,6 @@
if distutils.sysconfig.get_config_var('CC') is not None:
USE_ICC = 'icc' in distutils.sysconfig.get_config_var('CC')
-if not USE_ICC and USE_OPENMP:
- COMMON_MACROS.append(('_GLIBCXX_PARALLEL', None))
-
# simuPOP works with these boost versions. Newer versions will be used if these
# versions are not available, and will most likely work just fine.
boost_versions = ['1_35_0', '1_36_0', '1_37_0', '1_38_0', '1_39_0', '1_40_0',
@@ -557,9 +554,12 @@
COMMON_MACROS = [
('BOOST_UBLAS_NDEBUG', None),
- ('_HAS_ITERATOR_DEBUGGING', 0)
+ ('_HAS_ITERATOR_DEBUGGING', 0),
]
-
+
+if not USE_ICC and USE_OPENMP:
+ COMMON_MACROS.append(('_GLIBCXX_PARALLEL', None))
+
if os.name == 'nt':
COMMON_MACROS.extend([('BOOST_ALL_NO_LIB', None),
('NO_ZLIB', 0), ('NO_BZIP' , 1),
@@ -567,28 +567,20 @@
#('_SCL_SECURE_NO_WARNINGS', None)
])
+STD_MACROS = [('_SECURE_SCL', 1)]
+OPT_MACROS = [('NDEBUG', None), ('_SECURE_SCL', 0), ('OPTIMIZED', None)]
+
MACROS = {
- 'std': [('SIMUPOP_MODULE', 'simuPOP_std'),
- ('_SECURE_SCL', 1)],
- 'op': [('SIMUPOP_MODULE', 'simuPOP_op'), ('OPTIMIZED', None),
- ('NDEBUG', None), ('_SECURE_SCL', 0)],
- 'la': [('SIMUPOP_MODULE', 'simuPOP_la'), ('LONGALLELE', None),
- ('_SECURE_SCL', 1)],
- 'laop': [('SIMUPOP_MODULE', 'simuPOP_laop'), ('LONGALLELE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
- 'ba': [('SIMUPOP_MODULE', 'simuPOP_ba'), ('BINARYALLELE', None),
- ('_SECURE_SCL', 1)],
- 'baop': [('SIMUPOP_MODULE', 'simuPOP_baop'), ('BINARYALLELE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
- 'mu': [('SIMUPOP_MODULE', 'simuPOP_mu'), ('MUTANTALLELE', None),
- ('_SECURE_SCL', 1)],
- 'muop': [('SIMUPOP_MODULE', 'simuPOP_muop'), ('MUTANTALLELE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
- 'lin': [('SIMUPOP_MODULE', 'simuPOP_lin'), ('LINEAGE', None),
- ('_SECURE_SCL', 1)],
- 'linop': [('SIMUPOP_MODULE', 'simuPOP_linop'), ('LINEAGE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
-
+ 'std': [('SIMUPOP_MODULE', 'simuPOP_std')] + STD_MACROS,
+ 'op': [('SIMUPOP_MODULE', 'simuPOP_op') ] + OPT_MACROS,
+ 'la': [('SIMUPOP_MODULE', 'simuPOP_la'), ('LONGALLELE', None)] + STD_MACROS,
+ 'laop': [('SIMUPOP_MODULE', 'simuPOP_laop'), ('LONGALLELE', None)] + OPT_MACROS,
+ 'ba': [('SIMUPOP_MODULE', 'simuPOP_ba'), ('BINARYALLELE', None)] + STD_MACROS,
+ 'baop': [('SIMUPOP_MODULE', 'simuPOP_baop'), ('BINARYALLELE', None)] + OPT_MACROS,
+ 'mu': [('SIMUPOP_MODULE', 'simuPOP_mu'), ('MUTANTALLELE', None)] + STD_MACROS,
+ 'muop': [('SIMUPOP_MODULE', 'simuPOP_muop'), ('MUTANTALLELE', None)] + OPT_MACROS,
+ 'lin': [('SIMUPOP_MODULE', 'simuPOP_lin'), ('LINEAGE', None)] + STD_MACROS,
+ 'linop': [('SIMUPOP_MODULE', 'simuPOP_linop'),('LINEAGE', None)] + OPT_MACROS
}
WRAP_INFO = {
@@ -648,19 +640,52 @@
# copyGenotype for binary modules. Fortunately, the compiler provides options to
# continue to use libstdc++. Hopefully this quick fix can work a little bit longer.
#
-if is_maverick():
- common_extra_link_args = ['-stdlib=libstdc++']
- common_extra_include_dirs = ['/usr/include/c++/4.2.1']
- #common_extra_compile_args = ['-Wno-error=unused-command-line-argument-hard-error-in-future']
- common_extra_compile_args = ['-Wno-error=unused-command-line-argument']
-elif os.name == 'nt':
+
+common_library_dirs = ['build']
+
+if os.name == 'nt':
+ # I have a portable stdint.h for msvc, to avoid distributing
+ # zdll1.dll, I also build zlib from source
+ # zdll.lib is under win32
+ common_library_dirs.append('win32')
common_extra_link_args = []
- common_extra_include_dirs = []
- common_extra_compile_args = []
+ if PY3:
+ # Python3 uses VC 2010, which has stdint
+ common_extra_include_dirs = ['win32/zlib-1.2.3']
+ else:
+ common_extra_include_dirs = ['win32', 'win32/zlib-1.2.3']
+ # msvc does not have O3 option, /GR is to fix a C4541 warning
+ # /EHsc is for VC exception handling,
+ # /wd4819 disables warning messages for non-unicode character in boost/uitlity/enable_if.hpp
+ # /wd4996 disables warning messages for unsafe function call in boost/serialization
+ # /wd4068 disables warning messages for unknown pragma set by gcc
+ common_extra_compile_args = ['/O2', '/GR', '/EHsc', '/wd4819', '/wd4996', '/wd4068']
+ # Enable openMP if USE_OPENMP = True
+ if USE_OPENMP:
+ if USE_ICC:
+ common_extra_compile_args.append('/Qopenmp')
+ else:
+ common_extra_compile_args.append('/openmp')
else:
- common_extra_link_args = ['-static-libgcc', '-static-libstdc++']
+ common_extra_link_args = []
common_extra_include_dirs = []
- common_extra_compile_args = []
+ common_extra_compile_args = ['-O3', '-Wall', '-Wno-unknown-pragmas', '-Wno-unused-parameter']
+ if is_maverick():
+ common_extra_link_args.append('-stdlib=libstdc++')
+ common_extra_include_dirs.append('/usr/include/c++/4.2.1')
+ common_extra_compile_args.append('-Wno-error=unused-command-line-argument')
+ else:
+ common_extra_link_args.extend(['-static-libgcc', '-static-libstdc++'])
+ if not USE_ICC: # for gcc, turn on extra warning message
+ common_extra_compile_args.append('-Wextra')
+ if USE_OPENMP:
+ if USE_ICC:
+ common_extra_compile_args.append('-openmp')
+ else:
+ common_extra_compile_args.append('-fopenmp')
+ # if Intel ICC is used, turn off remark 981
+ if USE_ICC:
+ common_extra_compile_args.extend(['-wd981', '-wd191'])
#
def ModuInfo(modu, SIMUPOP_VER, SIMUPOP_REV):
@@ -675,9 +700,9 @@
#
# lib
if os.name == 'nt': # Windows, build zlib from source
- res['libraries'] = ['simuPOP_shared']
+ res['libraries'] = []
else:
- res['libraries'] = ['z', 'simuPOP_shared']
+ res['libraries'] = ['z']
if USE_OPENMP:
if USE_ICC:
res['libraries'].append('iomp5')
@@ -686,44 +711,12 @@
res['libraries'].append('gomp')
res['libraries'].extend(boost_lib_names)
res['include_dirs'] = ['.', 'gsl', boost_inc_path] + common_extra_include_dirs
- res['library_dirs'] = ['build']
- if os.name == 'nt':
- # I have a portable stdint.h for msvc, to avoid distributing
- # zdll1.dll, I also build zlib from source
- if PY3:
- # Python3 uses VC 2010, which has stdint
- res['include_dirs'].extend(['win32/zlib-1.2.3'])
- else:
- res['include_dirs'].extend(['win32', 'win32/zlib-1.2.3'])
- # zdll.lib is under win32
- res['library_dirs'].append('win32')
- if os.name == 'nt':
- # msvc does not have O3 option, /GR is to fix a C4541 warning
- # /EHsc is for VC exception handling,
- # /wd4819 disables warning messages for non-unicode character in boost/uitlity/enable_if.hpp
- # /wd4996 disables warning messages for unsafe function call in boost/serialization
- # /wd4068 disables warning messages for unknown pragma set by gcc
- res['extra_compile_args'] = ['/O2', '/GR', '/EHsc', '/wd4819', '/wd4996', '/wd4068'] + common_extra_compile_args
- # Enable openMP if USE_OPENMP = True
- if USE_OPENMP:
- if USE_ICC:
- res['extra_compile_args'].append('/Qopenmp')
- else:
- res['extra_compile_args'].append('/openmp')
- else:
- res['extra_compile_args'] = ['-O3', '-Wall', '-Wno-unknown-pragmas', '-Wno-unused-parameter'] + common_extra_compile_args
- if not USE_ICC: # for gcc, turn on extra warning message
- res['extra_compile_args'].append('-Wextra')
- if USE_OPENMP:
- if USE_ICC:
- res['extra_compile_args'].append('-openmp')
- else:
- res['extra_compile_args'].append('-fopenmp')
- # if Intel ICC is used, turn off remark 981
- if USE_ICC:
- res['extra_compile_args'].extend(['-wd981', '-wd191'])
# define_macros (deep copy)
res['define_macros'] = COMMON_MACROS + MACROS[modu]
+ if 'op' in modu:
+ res['libraries'].append('simuPOP_shop')
+ else:
+ res['libraries'].append('simuPOP_shstd')
res['undef_macros'] = []
return res
@@ -753,40 +746,51 @@
('SIMUPOP_REV', SIMUPOP_REV)
])
- if os.name != 'nt':
+ #distutils.log.set_verbosity(distutils.log.FATAL)
+ if os.name == 'nt':
+ NO_WARNING_ARG = [] # ['/nowarn'] is for vc 2015 only
+ else:
NO_WARNING_ARG = ['-w']
- else:
- NO_WARNING_ARG = []
try:
- distutils.log.set_verbosity(distutils.log.INFO)
- if os.name != 'nt' and (not os.path.isfile('src/boost_pch.h.pch') or \
- os.path.getmtime('src/boost_pch.h.pch') > os.path.getmtime('src/boost_pch.h')):
- c = new_compiler(verbose=1)
- # allow compiling .h file
- c.src_extensions.append('.hpp')
- c.obj_extension = '.hpp.gch'
- c.compile(['src/boost_pch.hpp'], output_dir='.',
- include_dirs=['.', 'gsr', boost_include_dir] + common_extra_include_dirs,
- extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
- macros=COMMON_MACROS)
- if not os.path.isfile('src/boost_pch.hpp.pch'):
- raise RuntimeError('Failed to pre-compile boost headers')
+ if not os.path.isfile('src/boost_pch.h.pch') or \
+ os.path.getmtime('src/boost_pch.h.pch') > os.path.getmtime('src/boost_pch.h'):
+ if os.name == 'nt':
+ pass
+ else:
+ c = new_compiler(verbose=1)
+ # allow compiling .h file
+ c.src_extensions.append('.hpp')
+ c.obj_extension = '.hpp.gch'
+ c.compile(['src/boost_pch.hpp'], output_dir='.',
+ include_dirs=['.', 'gsr', boost_include_dir] + common_extra_include_dirs,
+ extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
+ macros=COMMON_MACROS)
+ if not os.path.isfile('src/boost_pch.hpp.pch'):
+ raise RuntimeError('Failed to pre-compile boost headers')
except Exception as e:
# it is ok if boost_pch cannot be precompiled.
print('Failed to pre-compile boost headers: {}'.format(e))
try:
# try to get
- print('Building a static library')
+ print('Building static libraries')
c = new_compiler(verbose=1)
# -w suppress all warnings caused by the use of boost libraries
objects = c.compile(LIB_FILES,
include_dirs=['gsl', 'gsl/specfunc', 'build', '.', boost_include_dir] + common_extra_include_dirs,
output_dir='build',
extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
- macros = COMMON_MACROS
+ macros = COMMON_MACROS + STD_MACROS
)
- c.create_static_lib(objects, "simuPOP_shared", output_dir='build')
+ c.create_static_lib(objects, "simuPOP_shstd", output_dir='build')
+ #
+ objects = c.compile(LIB_FILES,
+ include_dirs=['gsl', 'gsl/specfunc', 'build', '.', boost_include_dir] + common_extra_include_dirs,
+ output_dir='build',
+ extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
+ macros = COMMON_MACROS + OPT_MACROS
+ )
+ c.create_static_lib(objects, "simuPOP_shop", output_dir='build')
except Exception as e:
sys.exit("Failed to build a shared supporting library: {}".format(e))
@@ -846,7 +850,7 @@
Extension('simuPOP._gsl',
sources = GSL_FILES + ['src/gsl_wrap.c'],
include_dirs = ['gsl', 'gsl/specfunc', 'build', '.'] + common_extra_include_dirs,
- extra_compile_args = common_extra_compile_args + ['-w'],
+ extra_compile_args = common_extra_compile_args + NO_WARNING_ARG,
extra_link_args = common_extra_link_args,
)
]
@@ -855,10 +859,10 @@
EXT_MODULES.append(
Extension('simuPOP._simuPOP_%s' % modu,
sources = info['src'],
- extra_compile_args = info['extra_compile_args'],
+ extra_compile_args = common_extra_compile_args,
# src for config.h etc, build for swigpyrun.h
include_dirs = info['include_dirs'] + ['src', 'build'],
- library_dirs = info['library_dirs'],
+ library_dirs = common_library_dirs,
libraries = info['libraries'],
define_macros = info['define_macros'],
undef_macros = info['undef_macros'],
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-02-10 16:26:31
|
Revision: 5013
http://sourceforge.net/p/simupop/code/5013
Author: bpeng2000
Date: 2016-02-10 16:26:29 +0000 (Wed, 10 Feb 2016)
Log Message:
-----------
Fix pre-compile header code
Modified Paths:
--------------
trunk/setup.py
Modified: trunk/setup.py
===================================================================
--- trunk/setup.py 2016-02-10 15:43:45 UTC (rev 5012)
+++ trunk/setup.py 2016-02-10 16:26:29 UTC (rev 5013)
@@ -44,8 +44,8 @@
from distutils.ccompiler import new_compiler
from distutils.errors import CompileError
from distutils.errors import DistutilsExecError
+import distutils.log
-
USE_SETUPTOOLS = False
try:
from setuptools import setup, find_packages, Extension
@@ -758,15 +758,19 @@
else:
NO_WARNING_ARG = []
try:
+ distutils.log.set_verbosity(distutils.log.INFO)
if os.name != 'nt' and (not os.path.isfile('src/boost_pch.h.pch') or \
os.path.getmtime('src/boost_pch.h.pch') > os.path.getmtime('src/boost_pch.h')):
c = new_compiler(verbose=1)
# allow compiling .h file
c.src_extensions.append('.hpp')
- c.compile(['src/boost_pch.hpp'], output_dir='src',
+ c.obj_extension = '.hpp.gch'
+ c.compile(['src/boost_pch.hpp'], output_dir='.',
include_dirs=['.', 'gsr', boost_include_dir] + common_extra_include_dirs,
extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
macros=COMMON_MACROS)
+ if not os.path.isfile('src/boost_pch.hpp.pch'):
+ raise RuntimeError('Failed to pre-compile boost headers')
except Exception as e:
# it is ok if boost_pch cannot be precompiled.
print('Failed to pre-compile boost headers: {}'.format(e))
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-02-10 15:43:48
|
Revision: 5012
http://sourceforge.net/p/simupop/code/5012
Author: bpeng2000
Date: 2016-02-10 15:43:45 +0000 (Wed, 10 Feb 2016)
Log Message:
-----------
Minor fixes of manifest file and suppress the last warning message
Modified Paths:
--------------
trunk/MANIFEST.in
trunk/simuPOP_version.py
trunk/src/stator.cpp
Modified: trunk/MANIFEST.in
===================================================================
--- trunk/MANIFEST.in 2016-02-10 04:18:10 UTC (rev 5011)
+++ trunk/MANIFEST.in 2016-02-10 15:43:45 UTC (rev 5012)
@@ -18,7 +18,7 @@
include config_linux.h
include config_mac.h
include src/simuPOP_cfg.h
-include src/boost_pch.h
+include src/boost_pch.hpp
include src/mutant_vector.h
include src/utility.h
include src/genoStru.h
@@ -258,23 +258,33 @@
include test/test_16_sampling.py
include test/test_17_utils.py
include test/test_18_plotter.py
-include test/test_18_demography.py
+include test/test_19_demography.py
include test/sample_32_ba.pop
include test/sample_32_la.pop
include test/sample_32_std.pop
+include test/sample_32_ba_v1.pop
+include test/sample_32_la_v1.pop
+include test/sample_32_lin_v1.pop
+include test/sample_32_mu_v1.pop
+include test/sample_32_std_v1.pop
+include test/sample_32_ba_v2.pop
+include test/sample_32_la_v2.pop
+include test/sample_32_lin_v2.pop
+include test/sample_32_mu_v2.pop
+include test/sample_32_std_v2.pop
include test/sample_64_ba.pop
include test/sample_64_la.pop
include test/sample_64_std.pop
-include test/sample_32_ba_v1.pop
-include test/sample_32_la_v1.pop
-include test/sample_32_std_v1.pop
include test/sample_64_ba_v1.pop
include test/sample_64_la_v1.pop
+include test/sample_64_lin_v1.pop
+include test/sample_64_mu_v1.pop
include test/sample_64_std_v1.pop
-include test/sample_64_mu_v1.pop
-include test/sample_64_muop_v1.pop
-include test/sample_64_lin_v1.pop
-include test/sample_64_linop_v1.pop
+include test/sample_64_ba_v2.pop
+include test/sample_64_la_v2.pop
+include test/sample_64_lin_v2.pop
+include test/sample_64_mu_v2.pop
+include test/sample_64_std_v2.pop
# documentations
include doc/userGuide.pdf
Modified: trunk/simuPOP_version.py
===================================================================
--- trunk/simuPOP_version.py 2016-02-10 04:18:10 UTC (rev 5011)
+++ trunk/simuPOP_version.py 2016-02-10 15:43:45 UTC (rev 5012)
@@ -1,2 +1,2 @@
-SIMUPOP_VER="1.1.8"
+SIMUPOP_VER="1.1.8svn"
SIMUPOP_REV="5003"
Modified: trunk/src/stator.cpp
===================================================================
--- trunk/src/stator.cpp 2016-02-10 04:18:10 UTC (rev 5011)
+++ trunk/src/stator.cpp 2016-02-10 15:43:45 UTC (rev 5012)
@@ -724,9 +724,11 @@
if (m_loci.empty())
return true;
- // get actual list of loci
+#ifndef MUTANTALLELE
+ // get actual list of loci, only used for non-mutant modules
const vectoru & loci = m_loci.elems(&pop);
DBG_DO(DBG_STATOR, cerr << "Count number of segregating sites for " << loci.size() << " loci " << endl);
+#endif
std::set<size_t> allFixedSites;
std::set<size_t> allSegSites;
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-02-10 04:18:13
|
Revision: 5011
http://sourceforge.net/p/simupop/code/5011
Author: bpeng2000
Date: 2016-02-10 04:18:10 +0000 (Wed, 10 Feb 2016)
Log Message:
-----------
Completely remove all warning messages for gcc
Modified Paths:
--------------
trunk/gsl/gsl_inline.h
trunk/setup.py
trunk/src/customizedTemplates.cpp
trunk/src/genoStru.cpp
trunk/src/individual.cpp
trunk/src/mating.h
trunk/src/migrator.cpp
trunk/src/operator.h
trunk/src/population.cpp
trunk/src/transmitter.cpp
trunk/src/utility.cpp
Modified: trunk/gsl/gsl_inline.h
===================================================================
--- trunk/gsl/gsl_inline.h 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/gsl/gsl_inline.h 2016-02-10 04:18:10 UTC (rev 5011)
@@ -48,8 +48,8 @@
#ifdef HAVE_INLINE
# if defined(__GNUC_STDC_INLINE__) || defined(GSL_C99_INLINE) || defined(HAVE_C99_INLINE)
-# define INLINE_DECL inline /* use C99 inline */
-# define INLINE_FUN inline
+# define INLINE_DECL static inline /* use C99 inline */
+# define INLINE_FUN static inline
# else
# define INLINE_DECL /* use GNU extern inline */
# define INLINE_FUN extern inline
Modified: trunk/setup.py
===================================================================
--- trunk/setup.py 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/setup.py 2016-02-10 04:18:10 UTC (rev 5011)
@@ -389,7 +389,7 @@
# since it is troublesome to link to external gsl library,
# I embed some GSL files with simuPOP.
-GSL_LIB_FILES = [
+LIB_FILES = [
'gsl/sys/infnan.c',
'gsl/sys/coerce.c',
'gsl/sys/fdiv.c',
@@ -474,9 +474,7 @@
'gsl/cdf/gamma.c',
'gsl/cdf/poisson.c',
'gsl/error.c'
-]
-
-LIB_FILES = [x for x in glob.glob(os.path.join(boost_serialization_dir, '*.cpp')) if 'xml' not in x and 'binary' not in x]\
+] + [x for x in glob.glob(os.path.join(boost_serialization_dir, '*.cpp')) if 'xml' not in x and 'binary' not in x]\
+ [x for x in glob.glob(os.path.join(boost_iostreams_dir, '*.cpp')) if 'bzip' not in x]\
+ glob.glob(os.path.join(boost_regex_dir, '*.cpp'))
@@ -675,12 +673,6 @@
for src in SOURCE_FILES:
res['src'].append('build/%s/%s' % (modu, src))
#
- # For some reason that I have not got a change to investigate, including GSL_LIB_FILES
- # to LIB_FILES would lead to symbol not found and some error like that. I suspect
- # that this is related to the order the library is linked, or some MACRO that is only
- # used when the GSL files are compiled with simuPOP modules.
- #
- res['src'].extend(GSL_LIB_FILES)
# lib
if os.name == 'nt': # Windows, build zlib from source
res['libraries'] = ['simuPOP_shared']
@@ -761,6 +753,10 @@
('SIMUPOP_REV', SIMUPOP_REV)
])
+ if os.name != 'nt':
+ NO_WARNING_ARG = ['-w']
+ else:
+ NO_WARNING_ARG = []
try:
if os.name != 'nt' and (not os.path.isfile('src/boost_pch.h.pch') or \
os.path.getmtime('src/boost_pch.h.pch') > os.path.getmtime('src/boost_pch.h')):
@@ -769,7 +765,7 @@
c.src_extensions.append('.hpp')
c.compile(['src/boost_pch.hpp'], output_dir='src',
include_dirs=['.', 'gsr', boost_include_dir] + common_extra_include_dirs,
- extra_preargs = common_extra_compile_args,
+ extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
macros=COMMON_MACROS)
except Exception as e:
# it is ok if boost_pch cannot be precompiled.
@@ -780,10 +776,10 @@
print('Building a static library')
c = new_compiler(verbose=1)
# -w suppress all warnings caused by the use of boost libraries
- objects = c.compile(LIB_FILES, extra_postargs=['-w', '-unknownnn'],
+ objects = c.compile(LIB_FILES,
include_dirs=['gsl', 'gsl/specfunc', 'build', '.', boost_include_dir] + common_extra_include_dirs,
output_dir='build',
- extra_preargs = common_extra_compile_args,
+ extra_preargs = common_extra_compile_args + NO_WARNING_ARG,
macros = COMMON_MACROS
)
c.create_static_lib(objects, "simuPOP_shared", output_dir='build')
Modified: trunk/src/customizedTemplates.cpp
===================================================================
--- trunk/src/customizedTemplates.cpp 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/customizedTemplates.cpp 2016-02-10 04:18:10 UTC (rev 5011)
@@ -91,7 +91,7 @@
/// CPPONLY
template <typename T>
int
-setarrayitem_template(struct arrayobject_template<T> * ap, Py_ssize_t i, PyObject * v)
+setarrayitem_template(struct arrayobject_template<T> * /* ap */, Py_ssize_t /* i */, PyObject * /* v */)
{
return -1;
}
@@ -641,7 +641,7 @@
/// CPPONLY
template<typename T>
-bool is_carrayobject_template(PyObject * op)
+bool is_carrayobject_template(PyObject * /* op */)
{
return false;
}
@@ -665,7 +665,7 @@
/// CPPONLY
template<typename T>
-PyObject * newcarrayobject_template(T begin, T end)
+PyObject * newcarrayobject_template(T /* begin */, T /* end */)
{
return(NULL);
}
Modified: trunk/src/genoStru.cpp
===================================================================
--- trunk/src/genoStru.cpp 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/genoStru.cpp 2016-02-10 04:18:10 UTC (rev 5011)
@@ -714,7 +714,6 @@
vectorstr::iterator pch1 = find(gs1.m_chromNames.begin(), gs1.m_chromNames.end(), gs2.m_chromNames[ch2]);
// if there is no matching chromosome, add the chromosome from gs2
if (pch1 == gs1.m_chromNames.end()) {
- size_t ch1 = pch1 - gs1.m_chromNames.begin();
loci.push_back(gs2.m_numLoci[ch2]);
chromNames.push_back(gs2.m_chromNames[ch2]);
chromTypes.push_back(gs2.m_chromTypes[ch2]);
Modified: trunk/src/individual.cpp
===================================================================
--- trunk/src/individual.cpp 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/individual.cpp 2016-02-10 04:18:10 UTC (rev 5011)
@@ -423,7 +423,6 @@
PyObject * Individual::mutAtLoci(const lociList & lociList)
{
ssize_t ply = ploidy();
- bool autosome_only = chromX() == -1 && chromY() == -1;
PyObject * mutDict = PyDefDict_New();
PyObject * dkey = NULL;
@@ -434,6 +433,7 @@
if (lociList.allAvail()) {
#ifdef MUTANTALLELE
+ bool autosome_only = chromX() == -1 && chromY() == -1;
if (autosome_only) {
vectorm::const_val_iterator m_ptr = m_genoPtr.get_val_iterator();
@@ -501,6 +501,7 @@
#endif
} else {
#ifdef MUTANTALLELE
+ bool autosome_only = chromX() == -1 && chromY() == -1;
if (autosome_only) {
const vectoru & loci = lociList.elems(this);
Modified: trunk/src/mating.h
===================================================================
--- trunk/src/mating.h 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/mating.h 2016-02-10 04:18:10 UTC (rev 5011)
@@ -791,7 +791,7 @@
/// Initialize a parent chooser for subpopulation \e subPop of \e population pop
- virtual void initialize(Population & pop, size_t subPop)
+ virtual void initialize(Population & /* pop */, size_t /* subPop */)
{
}
Modified: trunk/src/migrator.cpp
===================================================================
--- trunk/src/migrator.cpp 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/migrator.cpp 2016-02-10 04:18:10 UTC (rev 5011)
@@ -269,7 +269,7 @@
int begin, int end, int step, const intList & at,
const intList & reps, const subPopList & subPops, const stringList & infoFields)
: BaseOperator("", begin, end, step, at, reps, subPops, infoFields),
- m_rate(rate.elems()), m_inverse_rate(rate.elems().size(), rate.elems().size()), m_mode(mode), m_symmetric_matrix(true)
+ m_rate(rate.elems()), m_inverse_rate(rate.elems().size(), rate.elems().size()), m_symmetric_matrix(true), m_mode(mode)
{
DBG_FAILIF(!subPops.empty() && subPops.size() != m_rate.size(),
ValueError, "Length of param subPop must match rows of rate matrix.");
Modified: trunk/src/operator.h
===================================================================
--- trunk/src/operator.h 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/operator.h 2016-02-10 04:18:10 UTC (rev 5011)
@@ -209,6 +209,8 @@
/// CPPONLY
string infoField(size_t idx, const GenoStruTrait * trait = NULL) const
{
+ // suppress a warning about unused variable
+ (void)trait;
DBG_ASSERT(idx < m_infoFields.elems(trait).size(), IndexError,
(boost::format("Given info index %1% is out of range of 0 ~ %2%") % idx % m_infoFields.elems().size()).str());
return m_infoFields.elems()[idx];
Modified: trunk/src/population.cpp
===================================================================
--- trunk/src/population.cpp 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/population.cpp 2016-02-10 04:18:10 UTC (rev 5011)
@@ -1045,7 +1045,7 @@
for (size_t i = 0; i < infoFields.size(); ++i)
fields[i] = infoIdx(infoFields[i]);
for (size_t sp = 0; sp < numSubPop(); ++sp)
- parallelSort(rawIndBegin(sp), rawIndEnd(sp), indCompare(fields, reverse=reverse));
+ parallelSort(rawIndBegin(sp), rawIndEnd(sp), indCompare(fields, reverse));
setIndOrdered(false);
}
@@ -1552,7 +1552,7 @@
// if sp is being merged
if (find(sps.begin(), sps.end(), sp) != sps.end()) {
merged_size += subPopSize(sp);
- if (sp == toSubPop) {
+ if (sp == static_cast<size_t>(toSubPop)) {
merged_idx = new_size.size();
new_size.push_back(0);
if (!m_subPopNames.empty())
@@ -1578,7 +1578,7 @@
// find the new subpop order
vectoru sp_order;
// subpopulations before toSubPop
- for (size_t sp = 0; sp < toSubPop; ++sp)
+ for (size_t sp = 0; sp < static_cast<size_t>(toSubPop); ++sp)
if (find(sps.begin(), sps.end(), sp) == sps.end())
sp_order.push_back(sp);
// all merged subpopulations
Modified: trunk/src/transmitter.cpp
===================================================================
--- trunk/src/transmitter.cpp 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/transmitter.cpp 2016-02-10 04:18:10 UTC (rev 5011)
@@ -620,7 +620,10 @@
// prepare m_bt
vectorf vecP;
//
+#ifndef OPTIMIZED
+ // this is only used in non-optimized mode for checking
const vectoru & loci = m_loci.elems(&ind);
+#endif
DBG_FAILIF(m_intensity < 0 && (m_rates.empty() && !loci.empty()), ValueError,
"You should specify m_intensity, or m_rates "
Modified: trunk/src/utility.cpp
===================================================================
--- trunk/src/utility.cpp 2016-02-10 04:04:51 UTC (rev 5010)
+++ trunk/src/utility.cpp 2016-02-10 04:18:10 UTC (rev 5011)
@@ -2409,7 +2409,9 @@
// here we use version 2 because this is the latest version that supported by
// both python 2 and python 3, also because it is the one that handles simuPOP's
// defdict type using its __reduce__ interface.
- PyObject * pres = PyObject_CallMethod(pickle, "dumps", "(Oi)", m_dict, 2);
+ char dumps[] = "dumps";
+ char Oi[] = "(Oi)";
+ PyObject * pres = PyObject_CallMethod(pickle, dumps, Oi, m_dict, 2);
if (pres == NULL) {
PyErr_Clear();
/* If the dictionary is not pickleable, we have to go a longer way and test if
@@ -2422,18 +2424,19 @@
while (PyDict_Next(m_dict, &pos, &key, &value)) {
// try to pickle each item
// if the key is not pickleable, remove the key from the copied dictionary
- if (PyObject_CallMethod(pickle, "dumps", "(Oi)", key, 2) == NULL) {
+ // to avoid a warning
+ if (PyObject_CallMethod(pickle, dumps, Oi, key, 2) == NULL) {
PyErr_Clear();
PyDict_DelItem(new_dict, key);
}
- if (PyObject_CallMethod(pickle, "dumps", "(Oi)", value, 2) == NULL) {
+ if (PyObject_CallMethod(pickle, dumps, Oi, value, 2) == NULL) {
PyErr_Clear();
PyDict_DelItem(new_dict, key);
}
}
// now try to pickale the whole new object, in any case, the original
// dictionary is not touched.
- pres = PyObject_CallMethod(pickle, "dumps", "(Oi)", new_dict, 2);
+ pres = PyObject_CallMethod(pickle, dumps, Oi, new_dict, 2);
Py_DECREF(new_dict);
if (pres == NULL) {
PyErr_Print();
@@ -2508,7 +2511,9 @@
Py_XDECREF(m_dict);
}
m_ownVars = true;
- m_dict = PyObject_CallMethod(pickle, "loads", "(O)", args);
+ char loads[] = "loads";
+ char O[] = "(O)";
+ m_dict = PyObject_CallMethod(pickle, loads, O, args);
if (m_dict == NULL) {
PyErr_Print();
PyErr_Clear();
@@ -3059,7 +3064,7 @@
// all flags will be cleared to 0
StreamProvider::StreamProvider(const string & output, const pyFunc & func, const string & mode)
- : m_filename(output), m_filenameExpr(), m_func(func), m_flags(0), m_filePtr(NULL), m_mode(mode)
+ : m_filename(output), m_filenameExpr(), m_func(func), m_flags(0), m_mode(mode), m_filePtr(NULL)
{
if (!m_filename.empty() && m_filename[0] == '!') {
m_filenameExpr.setExpr(m_filename.substr(1));
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-02-10 04:04:53
|
Revision: 5010
http://sourceforge.net/p/simupop/code/5010
Author: bpeng2000
Date: 2016-02-10 04:04:51 +0000 (Wed, 10 Feb 2016)
Log Message:
-----------
Remove tests for rpy based plotter
Modified Paths:
--------------
trunk/test/test_18_plotter.py
Modified: trunk/test/test_18_plotter.py
===================================================================
--- trunk/test/test_18_plotter.py 2016-02-09 23:36:44 UTC (rev 5009)
+++ trunk/test/test_18_plotter.py 2016-02-10 04:04:51 UTC (rev 5010)
@@ -24,18 +24,9 @@
from simuPOP import *
from time import sleep
-hasRPy = True
-try:
- from simuPOP.plotter import *
-except (ImportError, RuntimeError, AttributeError):
- print("simuRPy can not be imported. Either rpy or r is not installed properly.")
- hasRPy = False
-
class TestPlotter(unittest.TestCase):
def testDerivedArgs(self):
'Testing class derivedARgs'
- if not hasRPy:
- return True
pop = Population(0)
pop.dvars().gen = 100
args = DerivedArgs(
@@ -93,8 +84,6 @@
def testVarPlotterBase(self):
'Testing byRep parameter of VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[50,50,100], ploidy=2, loci=[3,4], infoFields = ['migrate_to']),
rep=3)
@@ -113,8 +102,6 @@
def testVarPlotterByRep(self):
'Testing byRep parameter of VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[50,50,100], ploidy=2, loci=[3,4], infoFields = ['migrate_to']),
rep=3)
@@ -138,8 +125,6 @@
def testVarPlotterByDim(self):
'Testing byDim paramter of VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[200, 100], ploidy=2, loci=[3,4], infoFields = ['migrate_to']),
rep=3)
@@ -162,8 +147,6 @@
def testVarPlotterSelectedRep(self):
'Testing byDim paramter of VarPlotter using selected replicates'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[50, 50, 100], ploidy=2, loci=[3,4], infoFields = ['migrate_to']),
rep=5)
@@ -184,8 +167,6 @@
def testVarPlotterTogether(self):
'Testing plotting all lines together using VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[50, 50, 100], ploidy=2, loci=[3,4], infoFields = ['migrate_to']),
rep=3)
@@ -209,8 +190,6 @@
def testVarPlotterByRepDim(self):
'Testing byDim paramter of VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[200, 100], ploidy=2, loci=[3, 4], infoFields = ['migrate_to']),
rep=2)
@@ -234,8 +213,6 @@
def testVarPlotterSaveAs(self):
'Testing saveAs parameter of VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[50,50,100], ploidy=2, loci=[3,4], infoFields = ['migrate_to']),
rep=3)
@@ -293,8 +270,6 @@
def testVarPlotterPar(self):
'Testing parameter passing of VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[50,50,100], ploidy=2, loci=[3,4], infoFields=['migrate_to']),
rep=2)
@@ -317,8 +292,6 @@
def testVarPlotterHook(self):
'Testing ylim parameter of VarPlotter'
- if not hasRPy:
- return True
simu = Simulator(
Population(size=[50,50,100], ploidy=2, loci=[3,4], infoFields=['migrate_to']),
rep=3)
@@ -351,8 +324,6 @@
def testScatterPlotter(self):
'Testing ScatterPlotter'
- if not hasRPy:
- return True
import random
pop = Population([100, 200], infoFields=['x', 'y'])
initSex(pop)
@@ -373,8 +344,6 @@
def testScatterPlotterSP(self):
'Testing ScatterPlotter with multiple virtual subpopulations'
- if not hasRPy:
- return True
import random
pop = Population([100, 200], infoFields=['x', 'y'])
initSex(pop)
@@ -399,8 +368,6 @@
def testScatterPlotterSubSet(self):
'Testing ScatterPlotter with partial individuals'
- if not hasRPy:
- return True
import random
pop = Population([100, 200], infoFields=['x', 'y'])
initSex(pop)
@@ -426,274 +393,5 @@
sleep(1)
r.dev_off()
- def testInfoPlotterBase(self):
- 'Testing basic histogram using InfoPlotter'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 1000], infoFields=['x'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x')
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x')]),
- postOps = HistPlotter(infoFields='x', main='histogram of x in green',
- angle=60, col='green', step=2),
- #Pause(),
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testInfoPlotterFields(self):
- 'Testing Stat plotter with multiple fields and subpopulations'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 1000], infoFields=['x', 'y'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = [
- HistPlotter(infoFields=['x', 'y'],
- subPops = [(0, 0), (0, 1)],
- col_fld=['red', 'green'],
- density_sp=[5, 20], step=2,
- main_spfld=['Field x, Male', 'Field y, Male',
- 'Field x, Female', 'Field y, Female'],
- ),
- #Pause(),
- ],
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testInfoPlotterQQplot(self):
- 'Testing barplotter with multiple fields and subpopulations'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x', 'y'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- def qqline(r, data=None, **kwargs):
- r.qqline(data)
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = [
- QQPlotter(infoFields=['x', 'y'],
- subPops = [(0, 0), (0, 1)],
- col_fld=['red', 'green'],
- pch_sp=[1, 2],
- main_spfld=['Field x, Male', 'Field y, Male',
- 'Field x, Female', 'Field y, Female'],
- plotHook=qqline, step=2
- ),
- #Pause(),
- ],
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testInfoPlotterNoFunc(self):
- 'Testing the Stat plotter when no function is specified'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x', 'y'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- def qqplot(r, data=None, field=None, subPop=None, **kwargs):
- r.qqnorm(data, main='QQ drawn in plotHook fld: %s sp: %s' % (field, subPop))
- r.qqline(data)
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = [
- InfoPlotter(infoFields=['x', 'y'],
- subPops = [(0, 0), (0, 1)],
- plotHook=qqplot, step=2
- ),
- #Pause(),
- ],
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
-
- def testBoxPlotterBase(self):
- 'Testing the base boxplotter'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x')
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x')]),
- postOps = [
- BoxPlotter(infoFields='x', xlab='x', main='boxplot for field x',
- step=2),
- #Pause(),
- ],
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testBoxPlotterFields(self):
- 'Testing barplotter with multiple fields and subpopulations'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x', 'y'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = BoxPlotter(infoFields=['x', 'y'], step=2),
- #Pause(),
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testBoxPlotterFieldsAndSubPop(self):
- 'Testing BoxPlotter with both fields and subpopulation'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x', 'y'],
- subPopNames=['sp1', 'sp2'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = [
- BoxPlotter(infoFields=['x', 'y'],
- subPops=[(0, 0), (0, 1)], style='quantile',
- step=2, main='Boxplot with 2 fields x 2 subpops',
- col='green', horizontal=True,
- ),
- #Pause(),
- ],
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testBoxPlotterByField(self):
- 'Testing BoxPlotter separated by information field'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x', 'y'],
- subPopNames=['sp1', 'sp2'])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = BoxPlotter(infoFields=['x', 'y'], byField=True,
- subPops=[(0,0), (0,1)],
- step=2, col_fld=['yellow', 'green'],
- ),
- #Pause(),
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testBoxPlotterBySubPop(self):
- 'Testing BoxPlotter separated by subpopulation'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x', 'y'], subPopNames=['', ''])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = BoxPlotter(infoFields=['x', 'y'], bySubPop=True,
- subPops=[(0,0), (0,1)],
- step=2, horizontal_sp=[True, False],
- ),
- #Pause(),
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
- def testBoxPlotterByFieldSubPop(self):
- 'Testing BoxPlotter separated by both field and subpopulation'
- if not hasRPy:
- return True
- import random
- pop = Population([500, 100], infoFields=['x', 'y'], subPopNames=['', ''])
- initSex(pop)
- pop.setIndInfo([random.random() for i in range(100)], 'x', 0)
- pop.setIndInfo([1 + random.random() for i in range(100)], 'x', 1)
- pop.setIndInfo([random.random() for i in range(300)], 'y')
- pop.setVirtualSplitter(SexSplitter())
- pop.evolve(
- matingScheme=RandomMating(ops=[MendelianGenoTransmitter(),
- InheritTagger(PATERNAL, infoFields='x'),
- InheritTagger(MATERNAL, infoFields='y'),
- ]),
- postOps = BoxPlotter(infoFields=['x', 'y'], bySubPop=True,
- subPops=[(0,0), (0,1)], byField=True,
- step=2, col_spfld = ['red', 'green', 'yellow', 'blue'],
- ),
- #Pause(),
- gen = 5,
- )
- sleep(1)
- r.dev_off()
-
if __name__ == '__main__':
unittest.main()
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-02-09 23:36:46
|
Revision: 5009
http://sourceforge.net/p/simupop/code/5009
Author: bpeng2000
Date: 2016-02-09 23:36:44 +0000 (Tue, 09 Feb 2016)
Log Message:
-----------
Add pre-compile boost header
Added Paths:
-----------
trunk/src/boost_pch.hpp
Added: trunk/src/boost_pch.hpp
===================================================================
--- trunk/src/boost_pch.hpp (rev 0)
+++ trunk/src/boost_pch.hpp 2016-02-09 23:36:44 UTC (rev 5009)
@@ -0,0 +1,71 @@
+/**
+ * $File: boost_pch.h $
+ * $LastChangedDate: 2015-08-29 11:10:32 -0500 (Sat, 29 Aug 2015) $
+ * $Rev: 4983 $
+ *
+ * This file is part of simuPOP, a forward-time population genetics
+ * simulation environment. Please visit http://simupop.sourceforge.net
+ * for details.
+ *
+ * Copyright (C) 2004 - 2010 Bo Peng (bp...@md...)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#ifndef _BOOST_PCH_H
+#define _BOOST_PCH_H
+
+#if TR1_SUPPORT == 0
+# include <map>
+#elif TR1_SUPPORT == 1
+# include <unordered_map>
+#else
+# include <tr1/unordered_map>
+#endif
+
+
+#include <boost/format.hpp>
+
+#include <boost/serialization/vector.hpp>
+#include <boost/serialization/version.hpp>
+#include <boost/serialization/split_member.hpp>
+#include <boost/serialization/split_free.hpp>
+
+#include "boost/tuple/tuple.hpp"
+#include <boost/numeric/ublas/matrix.hpp>
+
+#include <boost/archive/text_iarchive.hpp>
+#include <boost/archive/text_oarchive.hpp>
+
+#include <boost/serialization/utility.hpp>
+#include <boost/serialization/vector.hpp>
+#include <boost/serialization/split_member.hpp>
+#include <boost/serialization/split_free.hpp>
+#include <boost/serialization/version.hpp>
+
+#include "boost/lambda/lambda.hpp"
+
+#include <boost/numeric/ublas/lu.hpp>
+#include <boost/numeric/ublas/io.hpp>
+
+#include <boost/iostreams/filtering_stream.hpp>
+#include <boost/iostreams/filter/gzip.hpp>
+#include <boost/iostreams/device/file.hpp>
+
+#include "boost/lexical_cast.hpp"
+#include "boost/pending/lowest_bit.hpp"
+
+#include "boost/regex.hpp"
+
+#endif
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-02-09 23:35:38
|
Revision: 5008
http://sourceforge.net/p/simupop/code/5008
Author: bpeng2000
Date: 2016-02-09 23:35:36 +0000 (Tue, 09 Feb 2016)
Log Message:
-----------
Use pre-compiled header and static library to suppress compiling warnings from boost and reduce build time
Modified Paths:
--------------
trunk/MANIFEST.in
trunk/setup.py
trunk/simuPOP_version.py
trunk/src/genoStru.cpp
trunk/src/genoStru.h
trunk/src/migrator.cpp
trunk/src/migrator.h
trunk/src/pedigree.h
trunk/src/penetrance.h
trunk/src/population.cpp
trunk/src/population.h
trunk/src/qtrait.h
trunk/src/selector.h
trunk/src/simuPOP_cfg.h
trunk/src/utility.cpp
Modified: trunk/MANIFEST.in
===================================================================
--- trunk/MANIFEST.in 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/MANIFEST.in 2016-02-09 23:35:36 UTC (rev 5008)
@@ -18,6 +18,7 @@
include config_linux.h
include config_mac.h
include src/simuPOP_cfg.h
+include src/boost_pch.h
include src/mutant_vector.h
include src/utility.h
include src/genoStru.h
Modified: trunk/setup.py
===================================================================
--- trunk/setup.py 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/setup.py 2016-02-09 23:35:36 UTC (rev 5008)
@@ -41,6 +41,7 @@
"""
import os, sys, platform, shutil, glob, re, tempfile, subprocess
import distutils.sysconfig
+from distutils.ccompiler import new_compiler
from distutils.errors import CompileError
from distutils.errors import DistutilsExecError
@@ -186,6 +187,9 @@
if distutils.sysconfig.get_config_var('CC') is not None:
USE_ICC = 'icc' in distutils.sysconfig.get_config_var('CC')
+if not USE_ICC and USE_OPENMP:
+ COMMON_MACROS.append(('_GLIBCXX_PARALLEL', None))
+
# simuPOP works with these boost versions. Newer versions will be used if these
# versions are not available, and will most likely work just fine.
boost_versions = ['1_35_0', '1_36_0', '1_37_0', '1_38_0', '1_39_0', '1_40_0',
@@ -385,7 +389,7 @@
# since it is troublesome to link to external gsl library,
# I embed some GSL files with simuPOP.
-LIB_FILES = [
+GSL_LIB_FILES = [
'gsl/sys/infnan.c',
'gsl/sys/coerce.c',
'gsl/sys/fdiv.c',
@@ -470,7 +474,9 @@
'gsl/cdf/gamma.c',
'gsl/cdf/poisson.c',
'gsl/error.c'
-] + [x for x in glob.glob(os.path.join(boost_serialization_dir, '*.cpp')) if 'xml' not in x and 'binary' not in x]\
+]
+
+LIB_FILES = [x for x in glob.glob(os.path.join(boost_serialization_dir, '*.cpp')) if 'xml' not in x and 'binary' not in x]\
+ [x for x in glob.glob(os.path.join(boost_iostreams_dir, '*.cpp')) if 'bzip' not in x]\
+ glob.glob(os.path.join(boost_regex_dir, '*.cpp'))
@@ -551,32 +557,39 @@
if not os.path.isdir('build'):
os.mkdir('build')
+COMMON_MACROS = [
+ ('BOOST_UBLAS_NDEBUG', None),
+ ('_HAS_ITERATOR_DEBUGGING', 0)
+ ]
+
+if os.name == 'nt':
+ COMMON_MACROS.extend([('BOOST_ALL_NO_LIB', None),
+ ('NO_ZLIB', 0), ('NO_BZIP' , 1),
+ # this one disables a lot of warnings about VC Checked iterators.
+ #('_SCL_SECURE_NO_WARNINGS', None)
+ ])
+
MACROS = {
'std': [('SIMUPOP_MODULE', 'simuPOP_std'),
- ('_SECURE_SCL', 1), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('_SECURE_SCL', 1)],
'op': [('SIMUPOP_MODULE', 'simuPOP_op'), ('OPTIMIZED', None),
- ('NDEBUG', None), ('BOOST_UBLAS_NDEBUG', None),
- ('_SECURE_SCL', 0), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('NDEBUG', None), ('_SECURE_SCL', 0)],
'la': [('SIMUPOP_MODULE', 'simuPOP_la'), ('LONGALLELE', None),
- ('_SECURE_SCL', 1), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('_SECURE_SCL', 1)],
'laop': [('SIMUPOP_MODULE', 'simuPOP_laop'), ('LONGALLELE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('BOOST_UBLAS_NDEBUG', None),
- ('_SECURE_SCL', 0), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
'ba': [('SIMUPOP_MODULE', 'simuPOP_ba'), ('BINARYALLELE', None),
- ('_SECURE_SCL', 1), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('_SECURE_SCL', 1)],
'baop': [('SIMUPOP_MODULE', 'simuPOP_baop'), ('BINARYALLELE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('BOOST_UBLAS_NDEBUG', None),
- ('_SECURE_SCL', 0), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
'mu': [('SIMUPOP_MODULE', 'simuPOP_mu'), ('MUTANTALLELE', None),
- ('_SECURE_SCL', 1), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('_SECURE_SCL', 1)],
'muop': [('SIMUPOP_MODULE', 'simuPOP_muop'), ('MUTANTALLELE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('BOOST_UBLAS_NDEBUG', None),
- ('_SECURE_SCL', 0), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
'lin': [('SIMUPOP_MODULE', 'simuPOP_lin'), ('LINEAGE', None),
- ('_SECURE_SCL', 1), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('_SECURE_SCL', 1)],
'linop': [('SIMUPOP_MODULE', 'simuPOP_linop'), ('LINEAGE', None),
- ('OPTIMIZED', None),('NDEBUG', None), ('BOOST_UBLAS_NDEBUG', None),
- ('_SECURE_SCL', 0), ('_HAS_ITERATOR_DEBUGGING', 0)],
+ ('OPTIMIZED', None),('NDEBUG', None), ('_SECURE_SCL', 0)],
}
@@ -640,7 +653,8 @@
if is_maverick():
common_extra_link_args = ['-stdlib=libstdc++']
common_extra_include_dirs = ['/usr/include/c++/4.2.1']
- common_extra_compile_args = ['-Wno-error=unused-command-line-argument-hard-error-in-future']
+ #common_extra_compile_args = ['-Wno-error=unused-command-line-argument-hard-error-in-future']
+ common_extra_compile_args = ['-Wno-error=unused-command-line-argument']
elif os.name == 'nt':
common_extra_link_args = []
common_extra_include_dirs = []
@@ -660,12 +674,18 @@
res['src'] = ['src/simuPOP_' + modu + '_wrap.cpp']
for src in SOURCE_FILES:
res['src'].append('build/%s/%s' % (modu, src))
- res['src'].extend(LIB_FILES)
+ #
+ # For some reason that I have not got a change to investigate, including GSL_LIB_FILES
+ # to LIB_FILES would lead to symbol not found and some error like that. I suspect
+ # that this is related to the order the library is linked, or some MACRO that is only
+ # used when the GSL files are compiled with simuPOP modules.
+ #
+ res['src'].extend(GSL_LIB_FILES)
# lib
if os.name == 'nt': # Windows, build zlib from source
- res['libraries'] = []
+ res['libraries'] = ['simuPOP_shared']
else:
- res['libraries'] = ['z']
+ res['libraries'] = ['z', 'simuPOP_shared']
if USE_OPENMP:
if USE_ICC:
res['libraries'].append('iomp5')
@@ -711,18 +731,8 @@
if USE_ICC:
res['extra_compile_args'].extend(['-wd981', '-wd191'])
# define_macros (deep copy)
- res['define_macros'] = [x for x in MACROS[modu]]
- res['define_macros'].extend([('SIMUPOP_VER', SIMUPOP_VER), ('SIMUPOP_REV', SIMUPOP_REV)])
- if os.name == 'nt':
- res['define_macros'].extend([('BOOST_ALL_NO_LIB', None),
- ('NO_ZLIB', 0), ('NO_BZIP' , 1),
- # this one disables a lot of warnings about VC Checked iterators.
- #('_SCL_SECURE_NO_WARNINGS', None)
- ])
- else:
- if not USE_ICC:
- if USE_OPENMP:
- res['define_macros'].append(('_GLIBCXX_PARALLEL', None))
+ res['define_macros'] = COMMON_MACROS + MACROS[modu]
+
res['undef_macros'] = []
return res
@@ -746,7 +756,40 @@
SIMUPOP_VER, SIMUPOP_REV = simuPOP_version()
# create source file for each module
MODULES = ['std', 'op', 'la', 'laop', 'ba', 'baop', 'mu', 'muop', 'lin', 'linop']
-
+ COMMON_MACROS.extend([
+ ('SIMUPOP_VER', SIMUPOP_VER),
+ ('SIMUPOP_REV', SIMUPOP_REV)
+ ])
+
+ try:
+ if os.name != 'nt' and (not os.path.isfile('src/boost_pch.h.pch') or \
+ os.path.getmtime('src/boost_pch.h.pch') > os.path.getmtime('src/boost_pch.h')):
+ c = new_compiler(verbose=1)
+ # allow compiling .h file
+ c.src_extensions.append('.hpp')
+ c.compile(['src/boost_pch.hpp'], output_dir='src',
+ include_dirs=['.', 'gsr', boost_include_dir] + common_extra_include_dirs,
+ extra_preargs = common_extra_compile_args,
+ macros=COMMON_MACROS)
+ except Exception as e:
+ # it is ok if boost_pch cannot be precompiled.
+ print('Failed to pre-compile boost headers: {}'.format(e))
+
+ try:
+ # try to get
+ print('Building a static library')
+ c = new_compiler(verbose=1)
+ # -w suppress all warnings caused by the use of boost libraries
+ objects = c.compile(LIB_FILES, extra_postargs=['-w', '-unknownnn'],
+ include_dirs=['gsl', 'gsl/specfunc', 'build', '.', boost_include_dir] + common_extra_include_dirs,
+ output_dir='build',
+ extra_preargs = common_extra_compile_args,
+ macros = COMMON_MACROS
+ )
+ c.create_static_lib(objects, "simuPOP_shared", output_dir='build')
+ except Exception as e:
+ sys.exit("Failed to build a shared supporting library: {}".format(e))
+
#
# Generate Wrapping files
#
@@ -803,7 +846,7 @@
Extension('simuPOP._gsl',
sources = GSL_FILES + ['src/gsl_wrap.c'],
include_dirs = ['gsl', 'gsl/specfunc', 'build', '.'] + common_extra_include_dirs,
- extra_compile_args = common_extra_compile_args,
+ extra_compile_args = common_extra_compile_args + ['-w'],
extra_link_args = common_extra_link_args,
)
]
Modified: trunk/simuPOP_version.py
===================================================================
--- trunk/simuPOP_version.py 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/simuPOP_version.py 2016-02-09 23:35:36 UTC (rev 5008)
@@ -1,2 +1,2 @@
-SIMUPOP_VER="1.1.8svn"
+SIMUPOP_VER="1.1.8"
SIMUPOP_REV="5003"
Modified: trunk/src/genoStru.cpp
===================================================================
--- trunk/src/genoStru.cpp 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/genoStru.cpp 2016-02-09 23:35:36 UTC (rev 5008)
@@ -25,7 +25,7 @@
#include "genoStru.h"
-#include "boost/lambda/lambda.hpp"
+#include "boost_pch.hpp"
using namespace boost::lambda;
Modified: trunk/src/genoStru.h
===================================================================
--- trunk/src/genoStru.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/genoStru.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -34,11 +34,7 @@
#include "utility.h"
#include "simuPOP_cfg.h"
-#include <boost/serialization/vector.hpp>
-#include <boost/serialization/version.hpp>
-#include <boost/serialization/split_member.hpp>
-#include <boost/serialization/split_free.hpp>
-
+#include "boost_pch.hpp"
#include <iterator>
using std::ostream;
using std::ostream_iterator;
Modified: trunk/src/migrator.cpp
===================================================================
--- trunk/src/migrator.cpp 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/migrator.cpp 2016-02-09 23:35:36 UTC (rev 5008)
@@ -26,9 +26,7 @@
#include "virtualSubPop.h"
#include "migrator.h"
-#include <boost/numeric/ublas/lu.hpp>
-#include <boost/numeric/ublas/io.hpp>
-
+#include "boost_pch.hpp"
namespace simuPOP {
Migrator::Migrator(const floatMatrix & rate, int mode, const uintList & toSubPops,
Modified: trunk/src/migrator.h
===================================================================
--- trunk/src/migrator.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/migrator.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -33,8 +33,9 @@
#include <list>
using std::list;
-#include <boost/numeric/ublas/matrix.hpp>
+#include "boost_pch.hpp"
+
#include <iostream>
using std::cout;
using std::endl;
Modified: trunk/src/pedigree.h
===================================================================
--- trunk/src/pedigree.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/pedigree.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -31,13 +31,7 @@
*/
#include "population.h"
-#if TR1_SUPPORT == 0
-# include <map>
-#elif TR1_SUPPORT == 1
-# include <unordered_map>
-#else
-# include <tr1/unordered_map>
-#endif
+#include "boost_pch.hpp"
namespace simuPOP {
Modified: trunk/src/penetrance.h
===================================================================
--- trunk/src/penetrance.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/penetrance.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -32,7 +32,7 @@
#include "utility.h"
#include "operator.h"
-#include "boost/tuple/tuple.hpp"
+#include "boost_pch.hpp"
#include <numeric>
using std::min;
Modified: trunk/src/population.cpp
===================================================================
--- trunk/src/population.cpp 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/population.cpp 2016-02-09 23:35:36 UTC (rev 5008)
@@ -27,9 +27,7 @@
#include "virtualSubPop.h"
// for file compression
-#include <boost/iostreams/filtering_stream.hpp>
-#include <boost/iostreams/filter/gzip.hpp>
-#include <boost/iostreams/device/file.hpp>
+#include "boost_pch.hpp"
#if PY_VERSION_HEX >= 0x03000000
# define PyInt_FromLong(x) PyLong_FromLong(x)
Modified: trunk/src/population.h
===================================================================
--- trunk/src/population.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/population.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -46,14 +46,7 @@
#include <deque>
using std::deque;
-#include <boost/archive/text_iarchive.hpp>
-#include <boost/archive/text_oarchive.hpp>
-#include <boost/serialization/utility.hpp>
-#include <boost/serialization/vector.hpp>
-#include <boost/serialization/split_member.hpp>
-#include <boost/serialization/split_free.hpp>
-#include <boost/serialization/version.hpp>
-
+#include "boost_pch.hpp"
#include "individual.h"
#include "virtualSubPop.h"
Modified: trunk/src/qtrait.h
===================================================================
--- trunk/src/qtrait.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/qtrait.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -32,7 +32,7 @@
#include "utility.h"
#include "operator.h"
-#include "boost/tuple/tuple.hpp"
+#include "boost_pch.hpp"
#include <numeric>
using std::min;
Modified: trunk/src/selector.h
===================================================================
--- trunk/src/selector.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/selector.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -32,7 +32,7 @@
#include "utility.h"
#include "operator.h"
-#include "boost/tuple/tuple.hpp"
+#include "boost_pch.hpp"
#include <numeric>
using std::min;
Modified: trunk/src/simuPOP_cfg.h
===================================================================
--- trunk/src/simuPOP_cfg.h 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/simuPOP_cfg.h 2016-02-09 23:35:36 UTC (rev 5008)
@@ -498,8 +498,9 @@
}
#pragma GCC system_header
-#include <boost/format.hpp>
+#include "boost_pch.hpp"
+
namespace simuPOP {
// standard library
Modified: trunk/src/utility.cpp
===================================================================
--- trunk/src/utility.cpp 2016-01-27 15:00:33 UTC (rev 5007)
+++ trunk/src/utility.cpp 2016-02-09 23:35:36 UTC (rev 5008)
@@ -52,7 +52,7 @@
using std::ifstream;
using std::ofstream;
-#include "boost/lexical_cast.hpp"
+#include "boost_pch.hpp"
// for data type lociList
#include "genoStru.h"
@@ -105,7 +105,7 @@
#include "boost/pending/lowest_bit.hpp"
using boost::lowest_bit;
-#include "boost/regex.hpp"
+
using boost::regex;
using boost::regex_match;
using boost::cmatch;
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-27 15:00:35
|
Revision: 5007
http://sourceforge.net/p/simupop/code/5007
Author: bpeng2000
Date: 2016-01-27 15:00:33 +0000 (Wed, 27 Jan 2016)
Log Message:
-----------
add a bit more test (using removeLoci with Fst)
Modified Paths:
--------------
trunk/test/test_08_stat.py
Modified: trunk/test/test_08_stat.py
===================================================================
--- trunk/test/test_08_stat.py 2016-01-27 14:56:20 UTC (rev 5006)
+++ trunk/test/test_08_stat.py 2016-01-27 15:00:33 UTC (rev 5007)
@@ -514,6 +514,10 @@
# change order, but the result should be the same
stat(pop, structure=[5, 3, 2])
self.assertAlmostEqual(pop.dvars().F_st, 0.0261665)
+ # or remove some loci
+ pop.removeLoci(keep=[5,2,3])
+ stat(pop, structure=ALL_AVAIL)
+ self.assertAlmostEqual(pop.dvars().F_st, 0.0261665)
def testHaploFreq(self):
'Testing calculation of haplotype frequency'
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-27 14:56:22
|
Revision: 5006
http://sourceforge.net/p/simupop/code/5006
Author: bpeng2000
Date: 2016-01-27 14:56:20 +0000 (Wed, 27 Jan 2016)
Log Message:
-----------
Add tests for Fst with a subset of loci
Modified Paths:
--------------
trunk/test/test_08_stat.py
Modified: trunk/test/test_08_stat.py
===================================================================
--- trunk/test/test_08_stat.py 2016-01-25 21:09:27 UTC (rev 5005)
+++ trunk/test/test_08_stat.py 2016-01-27 14:56:20 UTC (rev 5006)
@@ -494,6 +494,26 @@
self.assertAlmostEqual(pop.dvars().F_is, -0.2290202)
self.assertAlmostEqual(pop.dvars().F_it, -0.1034594)
self.assertAlmostEqual(pop.dvars().F_st, 0.1021633)
+ #
+ # test Fst with a subset of loci, and with unordered loci etc
+ pop = Population(size=[500,100,1000],
+ ploidy=2, loci = [10])
+ pop.setVirtualSplitter(RangeSplitter([[0,125], [125, 375], [375, 500],
+ [0, 50], [50, 80], [80, 100],
+ [0, 100],[100, 600], [600, 1000]]))
+ for genotype, subPop in zip(
+ [[0,0,1],[0,1,1],[1,1,0],[0,1,0],[0,0,1],[1,1,1],[0,1,1],[0,1,1],[1,0, 1]],
+ [(0, 0), (0, 1), (0, 2), (1, 3), (1, 4), (1, 5), (2, 6), (2, 7), (2, 8)]):
+ initGenotype(pop, genotype=genotype, subPops=[subPop])
+ stat(pop, structure=ALL_AVAIL)
+ self.assertAlmostEqual(pop.dvars().F_st, 0.0258261)
+ stat(pop, structure=[2])
+ self.assertAlmostEqual(pop.dvars().F_st, 0.0263055)
+ stat(pop, structure=[2, 3, 5])
+ self.assertAlmostEqual(pop.dvars().F_st, 0.0261665)
+ # change order, but the result should be the same
+ stat(pop, structure=[5, 3, 2])
+ self.assertAlmostEqual(pop.dvars().F_st, 0.0261665)
def testHaploFreq(self):
'Testing calculation of haplotype frequency'
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-25 21:09:29
|
Revision: 5005
http://sourceforge.net/p/simupop/code/5005
Author: bpeng2000
Date: 2016-01-25 21:09:27 +0000 (Mon, 25 Jan 2016)
Log Message:
-----------
Minor updates on users guide
Modified Paths:
--------------
trunk/doc/userGuide.lyx
Modified: trunk/doc/userGuide.lyx
===================================================================
--- trunk/doc/userGuide.lyx 2016-01-25 21:02:34 UTC (rev 5004)
+++ trunk/doc/userGuide.lyx 2016-01-25 21:09:27 UTC (rev 5005)
@@ -8803,12 +8803,12 @@
\begin_layout Standard
Backward migration matrices are widely used in theoretical population genetics
and coalescent based simulations.
- Instead of specifying the probability of migrating from another subpopulation
+ Instead of specifying the probability of migrating from one subpopulation
to another (namely how migration happens), such matrices specify the probabilit
y that individuals in a subpopulation originate from others (namely the
result of migration).
simuPOP simulates such models by converting backward migration matrices
- to foward ones using the theory describes below.
+ to foward ones using the theory described below.
Due to the limit of such models, simuPOP cannot simulate migration from/to
virtual subpopulatons, creation of new subpopulation, different source
and destination subpopulations, and will generate an error if the conversion
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-25 21:02:35
|
Revision: 5004
http://sourceforge.net/p/simupop/code/5004
Author: bpeng2000
Date: 2016-01-25 21:02:34 +0000 (Mon, 25 Jan 2016)
Log Message:
-----------
Starting the 1.1.8 release cycle
Modified Paths:
--------------
trunk/simuPOP_version.py
Modified: trunk/simuPOP_version.py
===================================================================
--- trunk/simuPOP_version.py 2016-01-21 05:26:44 UTC (rev 5003)
+++ trunk/simuPOP_version.py 2016-01-25 21:02:34 UTC (rev 5004)
@@ -1,2 +1,2 @@
-SIMUPOP_VER="1.1.7"
+SIMUPOP_VER="1.1.8svn"
SIMUPOP_REV="5003"
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-21 05:26:46
|
Revision: 5003
http://sourceforge.net/p/simupop/code/5003
Author: bpeng2000
Date: 2016-01-21 05:26:44 +0000 (Thu, 21 Jan 2016)
Log Message:
-----------
Minor fixes for windows and python3, prepare for simuPOP 1.1.7 release
Modified Paths:
--------------
trunk/MANIFEST.in
trunk/development/conda/meta.yaml
trunk/doc/refManual.lyx
trunk/doc/refManual.pdf
trunk/doc/userGuide.lyx
trunk/doc/userGuide.pdf
trunk/simuPOP_version.py
trunk/src/initializer.cpp
trunk/src/utility.h
Modified: trunk/MANIFEST.in
===================================================================
--- trunk/MANIFEST.in 2016-01-20 04:47:23 UTC (rev 5002)
+++ trunk/MANIFEST.in 2016-01-21 05:26:44 UTC (rev 5003)
@@ -322,7 +322,6 @@
prune boost_1_49_0/boost/spirit
prune boost_1_49_0/boost/test
prune boost_1_49_0/boost/thread
-prune boost_1_49_0/boost/typeof
prune boost_1_49_0/boost/units
prune boost_1_49_0/boost/variant
prune boost_1_49_0/boost/wave
Modified: trunk/development/conda/meta.yaml
===================================================================
--- trunk/development/conda/meta.yaml 2016-01-20 04:47:23 UTC (rev 5002)
+++ trunk/development/conda/meta.yaml 2016-01-21 05:26:44 UTC (rev 5003)
@@ -1,10 +1,10 @@
package:
name: simupop
- version: !!str 1.1.6
+ version: !!str 1.1.7
source:
- fn: simuPOP-1.1.6.tar.gz
- url: https://binstar.org/simupop/simuPOP/1.1.6/download/simuPOP-1.1.6-src.tar.gz
+ fn: simuPOP-1.1.7.tar.gz
+ url: https://binstar.org/bpeng/simuPOP/1.1.7/download/simuPOP-1.1.7-src.tar.gz
md5:
# patches:
# List any patch files here
Modified: trunk/doc/refManual.lyx
===================================================================
--- trunk/doc/refManual.lyx 2016-01-20 04:47:23 UTC (rev 5002)
+++ trunk/doc/refManual.lyx 2016-01-21 05:26:44 UTC (rev 5003)
@@ -5,7 +5,7 @@
\textclass manual
\begin_preamble
\renewcommand{\py@ptsize}{12pt}
-\setreleaseinfo{Release 1.1.6 (\mbox{Rev: 4972})}
+\setreleaseinfo{Release 1.1.7 (\mbox{Rev: 5000})}
% file revision $Rev: 3372$
\authoraddress{
{\bf Department of Epidemiology, U.T. M.D. Anderson Cancer Center}\\
Modified: trunk/doc/refManual.pdf
===================================================================
(Binary files differ)
Modified: trunk/doc/userGuide.lyx
===================================================================
--- trunk/doc/userGuide.lyx 2016-01-20 04:47:23 UTC (rev 5002)
+++ trunk/doc/userGuide.lyx 2016-01-21 05:26:44 UTC (rev 5003)
@@ -6,7 +6,7 @@
\begin_preamble
\renewcommand{\py@ptsize}{12pt}
-\setreleaseinfo{Release 1.1.6 (\mbox{Rev: 4972})}
+\setreleaseinfo{Release 1.1.7 (\mbox{Rev: 5000})}
% file revision $Rev: 3372$
\authoraddress{
{\bf Department of Epidemiology, U.T. M.D. Anderson Cancer Center}\\
Modified: trunk/doc/userGuide.pdf
===================================================================
(Binary files differ)
Modified: trunk/simuPOP_version.py
===================================================================
--- trunk/simuPOP_version.py 2016-01-20 04:47:23 UTC (rev 5002)
+++ trunk/simuPOP_version.py 2016-01-21 05:26:44 UTC (rev 5003)
@@ -1,2 +1,2 @@
-SIMUPOP_VER="1.1.6"
-SIMUPOP_REV="4972"
+SIMUPOP_VER="1.1.7"
+SIMUPOP_REV="5003"
Modified: trunk/src/initializer.cpp
===================================================================
--- trunk/src/initializer.cpp 2016-01-20 04:47:23 UTC (rev 5002)
+++ trunk/src/initializer.cpp 2016-01-21 05:26:44 UTC (rev 5003)
@@ -202,6 +202,10 @@
}
+#if PY_VERSION_HEX >= 0x03000000
+# define PyInt_FromLong(x) PyLong_FromLong(x)
+#endif
+
vectorf InitGenotype::getFreqOrProp(size_t locus, const vspID & vsp) const
{
pyFunc func = m_freq.func().isValid() ? m_freq.func() : m_prop.func();
Modified: trunk/src/utility.h
===================================================================
--- trunk/src/utility.h 2016-01-20 04:47:23 UTC (rev 5002)
+++ trunk/src/utility.h 2016-01-21 05:26:44 UTC (rev 5003)
@@ -92,9 +92,14 @@
// this is used to compare loci positions when loci are provided by
// (chr, pos) pair.
-#define PRECISION(f) ((round(f)*100000.)/100000.)
-//double(size_t(f)) == f ? f : (size_t((f) * 10000. + 0.5) / 10000.))
+#if defined (_WIN32) || defined (_WIN64)
+// windows VC does not have round function
+#define PRECISION(f) (double(size_t(f)) == f ? f : (size_t((f) * 100000. + 0.5) / 100000.))
+#else
+#define PRECISION(f) (round(f*100000.)/100000.)
+#endif
+
namespace simuPOP {
// ////////////////////////////////////////////////////////////
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-20 04:47:25
|
Revision: 5002
http://sourceforge.net/p/simupop/code/5002
Author: bpeng2000
Date: 2016-01-20 04:47:23 +0000 (Wed, 20 Jan 2016)
Log Message:
-----------
Minor updates of users guide
Modified Paths:
--------------
trunk/doc/userGuide.lyx
Modified: trunk/doc/userGuide.lyx
===================================================================
--- trunk/doc/userGuide.lyx 2016-01-20 04:40:18 UTC (rev 5001)
+++ trunk/doc/userGuide.lyx 2016-01-20 04:47:23 UTC (rev 5002)
@@ -7503,7 +7503,7 @@
\end_layout
\begin_layout Itemize
-By multiple allele frequency or proportion returned by a function passed
+By multiple allele frequencies or proportions returned by a function passed
to parameter
\family typewriter
freq
@@ -7512,7 +7512,7 @@
\family typewriter
prop
\family default
- (new in version ).
+ (new in version 1.1.7).
This function can accept parameters
\family typewriter
loc
@@ -7532,9 +7532,9 @@
\family typewriter
freq=lambda : random.random()
\family default
- so that a new frequency will be drawn from a uniform distribution for each
+ so that a new frequency is drawn from an uniform distribution for each
new locus.
- Note that simuPOP expects the return value from this function to be a list
+ Note that simuPOP expects the return value of this function to be a list
of frequencies for alleles 0, 1, ..., but treats a single return value
\emph on
x
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-20 04:40:21
|
Revision: 5001
http://sourceforge.net/p/simupop/code/5001
Author: bpeng2000
Date: 2016-01-20 04:40:18 +0000 (Wed, 20 Jan 2016)
Log Message:
-----------
Allow parameters freq and prop of operator InitGenotype to take functions to initialize genotype with locus or subpopulation specific frequencies.
Modified Paths:
--------------
trunk/doc/refManual.pdf
trunk/doc/userGuide.lyx
trunk/doc/userGuide.pdf
trunk/doc/userGuide.py
trunk/src/initializer.cpp
trunk/src/initializer.h
trunk/src/simuPOP_doc.i
trunk/src/utility.h
trunk/test/test_06_initialization.py
Modified: trunk/doc/refManual.pdf
===================================================================
(Binary files differ)
Modified: trunk/doc/userGuide.lyx
===================================================================
--- trunk/doc/userGuide.lyx 2016-01-18 22:42:45 UTC (rev 5000)
+++ trunk/doc/userGuide.lyx 2016-01-20 04:40:18 UTC (rev 5001)
@@ -7502,6 +7502,50 @@
to the second copy.
\end_layout
+\begin_layout Itemize
+By multiple allele frequency or proportion returned by a function passed
+ to parameter
+\family typewriter
+freq
+\family default
+ or
+\family typewriter
+prop
+\family default
+ (new in version ).
+ This function can accept parameters
+\family typewriter
+loc
+\family default
+,
+\family typewriter
+subPop
+\family default
+ or
+\family typewriter
+vsp
+\family default
+ and returns locus, subpopopulation or virtual subpopulation specific allele
+ frequencies.
+ For example, if you would like to initialize genotypes with random allele
+ frequency, you can set
+\family typewriter
+freq=lambda : random.random()
+\family default
+ so that a new frequency will be drawn from a uniform distribution for each
+ new locus.
+ Note that simuPOP expects the return value from this function to be a list
+ of frequencies for alleles 0, 1, ..., but treats a single return value
+\emph on
+x
+\emph default
+ as [
+\emph on
+x, 1-x
+\emph default
+] for simplicity.
+\end_layout
+
\begin_layout Standard
Parameter
\family typewriter
Modified: trunk/doc/userGuide.pdf
===================================================================
(Binary files differ)
Modified: trunk/doc/userGuide.py
===================================================================
--- trunk/doc/userGuide.py 2016-01-18 22:42:45 UTC (rev 5000)
+++ trunk/doc/userGuide.py 2016-01-20 04:40:18 UTC (rev 5001)
@@ -2082,42 +2082,58 @@
#begin_ignore
sim.setRNG(seed=12345)
#end_ignore
-pop = sim.Population(size=[20, 30], loci=[5, 7])
+pop = sim.Population(size=[2000, 3000], loci=[5, 7])
# by allele frequency
+def printFreq(pop, loci):
+ sim.stat(pop, alleleFreq=loci)
+ print(', '.join(['{:.3f}'.format(pop.dvars().alleleFreq[x][0]) for x in loci]))
+
sim.initGenotype(pop, freq=[.4, .6])
sim.dump(pop, max=6, structure=False)
-sim.stat(pop, alleleFreq=range(12))
-print(['%.2f' % pop.dvars().alleleFreq[x][0] for x in range(5)])
+printFreq(pop, range(5))
# by proportion
sim.initGenotype(pop, prop=[0.4, 0.6])
-sim.stat(pop, alleleFreq=range(12))
-print(['%.2f' % pop.dvars().alleleFreq[x][0] for x in range(5)])
+printFreq(pop, range(5))
# by haplotype frequency
sim.initGenotype(pop, freq=[.4, .6], haplotypes=[[1, 1, 0, 1], [0, 0, 1]])
sim.dump(pop, max=6, structure=False)
-sim.stat(pop, alleleFreq=range(12))
-print(['%.2f' % pop.dvars().alleleFreq[x][0] for x in range(5)])
+printFreq(pop, range(5))
# by haplotype proportion
sim.initGenotype(pop, prop=[0.4, 0.6], haplotypes=[[1, 1, 0], [0, 0, 1, 1]])
-sim.stat(pop, alleleFreq=range(12))
-print(['%.2f' % pop.dvars().alleleFreq[x][0] for x in range(5)])
+printFreq(pop, range(5))
# by genotype
pop = sim.Population(size=[2, 3], loci=[5, 7])
sim.initGenotype(pop, genotype=[1]*5 + [2]*7 + [3]*5 +[4]*7)
sim.dump(pop, structure=False)
#
# use virtual subpopulation
+pop = sim.Population(size=[2000, 3000], loci=[5, 7])
pop.setVirtualSplitter(sim.SexSplitter())
sim.initSex(pop)
sim.initGenotype(pop, genotype=range(10), loci=range(5))
# initialize all males
sim.initGenotype(pop, genotype=[2]*7, loci=range(5, 12),
subPops=[(0, 0), (1, 0)])
-sim.dump(pop, structure=False)
# assign genotype by proportions
pop.setVirtualSplitter(sim.ProportionSplitter([0.4, 0.6]))
sim.initGenotype(pop, freq=[0.2, 0.8], subPops=[(0,0)])
sim.initGenotype(pop, freq=[0.5, 0.5], subPops=[(0,1)])
+#
+# initialize by random allele frequency
+import random
+sim.initGenotype(pop, freq=lambda : random.random())
+printFreq(pop, range(5))
+# initialize with loci specific frequency. here
+# lambda loc: 0.01*loc is equivalent to
+# lambda loc: [0.01*loc, 1-0.01*loc]
+sim.initGenotype(pop,
+ freq=lambda loc: 0.01*loc)
+printFreq(pop, range(5))
+# initialize with VSP-specific frequency
+sim.initGenotype(pop,
+ freq=lambda vsp: [[0.2, 0.8], [0.5, 0.5]][vsp[1]],
+ subPops=[(0, 0), (0, 1)])
+
#end_file
#begin_file log/InitInfo.py
Modified: trunk/src/initializer.cpp
===================================================================
--- trunk/src/initializer.cpp 2016-01-18 22:42:45 UTC (rev 5000)
+++ trunk/src/initializer.cpp 2016-01-20 04:40:18 UTC (rev 5001)
@@ -156,8 +156,8 @@
}
-InitGenotype::InitGenotype(const vectorf & freq,
- const uintList & genotype, const vectorf & prop,
+InitGenotype::InitGenotype(const floatListFunc & freq,
+ const uintList & genotype, const floatListFunc & prop,
const intMatrix & haplotypes,
const lociList & loci,
const uintList & ploidy,
@@ -169,21 +169,21 @@
m_haplotypes(haplotypes.elems()), m_loci(loci),
m_ploidy(ploidy)
{
- for (size_t i = 0; i < m_freq.size(); ++i)
- PARAM_FAILIF(fcmp_lt(m_freq[i], 0) || fcmp_gt(m_freq[i], 1), ValueError,
+ for (size_t i = 0; i < m_freq.elems().size(); ++i)
+ PARAM_FAILIF(fcmp_lt(m_freq.elems()[i], 0) || fcmp_gt(m_freq.elems()[i], 1), ValueError,
"Allele frequency should be between 0 and 1");
- for (size_t i = 0; i < m_prop.size(); ++i)
- PARAM_FAILIF(fcmp_lt(m_prop[i], 0) || fcmp_gt(m_prop[i], 1), ValueError,
+ for (size_t i = 0; i < m_prop.elems().size(); ++i)
+ PARAM_FAILIF(fcmp_lt(m_prop.elems()[i], 0) || fcmp_gt(m_prop.elems()[i], 1), ValueError,
"Allele proportion should be between 0 and 1");
- PARAM_FAILIF(!m_haplotypes.empty() && !m_freq.empty() && m_haplotypes.size() != m_freq.size(),
+ PARAM_FAILIF(!m_haplotypes.empty() && !m_freq.elems().empty() && m_haplotypes.size() != m_freq.size(),
ValueError, "Haplotype frequency, if specified, should be specified for each haplotype.")
- PARAM_FAILIF(!m_haplotypes.empty() && !m_prop.empty() && m_haplotypes.size() != m_prop.size(),
+ PARAM_FAILIF(!m_haplotypes.empty() && !m_prop.elems().empty() && m_haplotypes.size() != m_prop.size(),
ValueError, "Haplotype proportion, if specified, should be specified for each haplotype.")
- PARAM_FAILIF(!m_freq.empty() && fcmp_ne(accumulate(m_freq.begin(), m_freq.end(), 0.), 1.0), ValueError,
+ PARAM_FAILIF(!m_freq.elems().empty() && fcmp_ne(accumulate(m_freq.elems().begin(), m_freq.elems().end(), 0.), 1.0), ValueError,
"Allele frequencies should add up to one.");
- PARAM_FAILIF(!m_prop.empty() && fcmp_ne(accumulate(m_prop.begin(), m_prop.end(), 0.), 1.0), ValueError,
+ PARAM_FAILIF(!m_prop.elems().empty() && fcmp_ne(accumulate(m_prop.elems().begin(), m_prop.elems().end(), 0.), 1.0), ValueError,
"Allele proportions should add up to one.");
- PARAM_ASSERT(!m_haplotypes.empty() || (!m_freq.empty() + !m_genotype.empty() + !m_prop.empty() == 1), ValueError,
+ PARAM_ASSERT(!m_haplotypes.empty() || (m_freq.valid() + !m_genotype.empty() + m_prop.valid() == 1), ValueError,
"Please specify one and only one of parameters freq, prop (or haplotypes) and genotype.");
}
@@ -192,9 +192,9 @@
{
string desc = "<simuPOP.InitGenotype> initialize individual genotype ";
- if (!m_freq.empty())
+ if (m_freq.valid())
desc += string("acccording to ") + (m_haplotypes.empty() ? " allele " : " haplotype ") + "frequencies.";
- else if (m_prop.empty())
+ else if (m_prop.valid())
desc += string("according to proportion of ") + (m_haplotypes.empty() ? " alleles." : " haplotypes.");
else
desc += "using specified genotype.";
@@ -202,6 +202,43 @@
}
+vectorf InitGenotype::getFreqOrProp(size_t locus, const vspID & vsp) const
+{
+ pyFunc func = m_freq.func().isValid() ? m_freq.func() : m_prop.func();
+
+ PyObject * args = PyTuple_New(func.numArgs());
+ for (size_t i = 0; i < func.numArgs(); ++i) {
+ const string & arg = func.arg(i);
+ if (arg == "loc")
+ PyTuple_SET_ITEM(args, i, PyInt_FromLong(locus));
+ else if (arg == "subPop")
+ PyTuple_SET_ITEM(args, i, PyInt_FromLong(vsp.subPop()));
+ else if (arg == "vsp") {
+ PyObject * v = PyTuple_New(2);
+ PyTuple_SET_ITEM(v, 0, PyInt_FromLong(vsp.subPop()));
+ if (vsp.isVirtual())
+ PyTuple_SET_ITEM(v, 1, PyInt_FromLong(vsp.virtualSubPop()));
+ else {
+ Py_INCREF(Py_None);
+ PyTuple_SET_ITEM(v, 1, Py_None);
+ }
+ // it steals reference count of v
+ PyTuple_SET_ITEM(args, i, v);
+ }
+ }
+ vectorf freq = func(PyObj_As_Array, args);
+ if (freq.size() == 1 && freq.back() != 1.0)
+ freq.push_back(1.0 - freq[0]);
+ for (size_t i = 0; i < freq.size(); ++i) {
+ if (fcmp_lt(freq[i], 0) || fcmp_gt(freq[i], 1))
+ throw RuntimeError("Returned allele frequency should be between 0 and 1.");
+ }
+ if (fcmp_ne(std::accumulate(freq.begin(), freq.end(), 0.0), 1.0))
+ throw RuntimeError("Returned allele frequency should add up to 1.");
+ Py_XDECREF(args);
+ return freq;
+}
+
bool InitGenotype::apply(Population & pop) const
{
const subPopList subPops = applicableSubPops(pop);
@@ -224,6 +261,19 @@
subPopList::const_iterator sp_end = subPops.end();
size_t sz = m_genotype.size();
+ // check function
+ int by_locus = false;
+ if (m_freq.func().isValid() || m_prop.func().isValid()) {
+ pyFunc func = m_freq.func().isValid() ? m_freq.func() : m_prop.func();
+ //
+ for (size_t fp = 0; fp < func.numArgs(); ++fp) {
+ if (func.arg(fp) == "loc")
+ by_locus = true;
+ else if (func.arg(fp) != "subPop" && func.arg(fp) != "vsp")
+ throw ValueError("Unacceptable parameter " + func.arg(fp));
+ }
+ }
+
DBG_WARNIF(sz >= 1 && sz != loci.size() && sz != loci.size() * pop.ploidy(),
(boost::format("Genotype [%1%] specified in operator InitGenotype has %2% alleles, which does not match"
" number of loci of individuals. This sequence will be truncated or repeated and may lead to erroneous results.")
@@ -233,59 +283,98 @@
pop.activateVirtualSubPop(*sp);
if (!m_genotype.empty()) {
// multi-thread write to compressed mutant is not allowed.
- {
- IndIterator it = pop.indIterator(sp->subPop());
- for (size_t ii = 0; it.valid(); ++it, ++ii)
- for (vectoru::iterator p = ploidy.begin(); p != ploidy.end(); ++p)
- for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc, ++idx)
- it->setAllele(TO_ALLELE(m_genotype[idx % sz]), *loc, static_cast<int>(*p));
-
- }
- } else if (!m_prop.empty()) {
+ IndIterator it = pop.indIterator(sp->subPop());
+ for (size_t ii = 0; it.valid(); ++it, ++ii)
+ for (vectoru::iterator p = ploidy.begin(); p != ploidy.end(); ++p)
+ for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc, ++idx)
+ it->setAllele(TO_ALLELE(m_genotype[idx % sz]), *loc, static_cast<int>(*p));
+ } else if (m_prop.valid()) {
WeightedSampler ws;
size_t sz = pop.subPopSize(*sp);
+ // different for each locus
if (m_haplotypes.empty()) {
// initialize by allele. Gurantee proportion at each locus.
for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc) {
- ws.set(m_prop.begin(), m_prop.end(), sz * ploidy.size());
+ if (m_prop.elems().empty()) {
+ vectorf prop = getFreqOrProp(*loc, *sp);
+ ws.set(prop.begin(), prop.end(), sz * ploidy.size());
+ } else
+ ws.set(m_prop.elems().begin(), m_prop.elems().end(), sz * ploidy.size());
IndIterator it = pop.indIterator(sp->subPop());
for (; it.valid(); ++it)
for (vectoru::iterator p = ploidy.begin(); p != ploidy.end(); ++p)
it->setAllele(TO_ALLELE(ws.draw()), *loc, static_cast<int>(*p));
}
+ // use haplotype
} else {
- ws.set(m_prop.begin(), m_prop.end(), sz * ploidy.size());
- {
+ IndIterator it = pop.indIterator(sp->subPop());
+ if (m_prop.elems().empty()) {
+ if (by_locus)
+ throw ValueError("Return value of user provided function cannot have parameter loc when haplotypes are provided.");
+ vectorf prop = getFreqOrProp(0, *sp);
+ ws.set(prop.begin(), prop.end(), sz * ploidy.size());
+ } else
+ ws.set(m_prop.elems().begin(), m_prop.elems().end(), sz * ploidy.size());
+ for (; it.valid(); ++it)
+ for (vectoru::iterator p = ploidy.begin(); p != ploidy.end(); ++p) {
+ const vectori & haplotype = m_haplotypes[ws.draw()];
+ size_t hapSz = haplotype.size();
+ size_t j = 0;
+ for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc, ++j)
+ it->setAllele(TO_ALLELE(haplotype[j % hapSz]), *loc, static_cast<int>(*p));
+ }
+ }
+ } else {
+ // use allele frequency
+ if (m_haplotypes.empty()) {
+ // if a funciton is provided.
+ if (m_freq.func().isValid()) {
+ for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc) {
+ WeightedSampler ws(getFreqOrProp(*loc, *sp));
+ //
+ IndIterator it = pop.indIterator(sp->subPop());
+ for (; it.valid(); ++it) {
+ for (vectoru::iterator p = ploidy.begin(); p != ploidy.end(); ++p) {
+ it->setAllele(TO_ALLELE(ws.draw()), *loc, static_cast<int>(*p));
+ }
+ }
+ }
+ } else {
+ // use specified frequencies
IndIterator it = pop.indIterator(sp->subPop());
- for (; it.valid(); ++it)
+ WeightedSampler ws(m_freq.elems());
+ for (; it.valid(); ++it) {
for (vectoru::iterator p = ploidy.begin(); p != ploidy.end(); ++p) {
- const vectori & haplotype = m_haplotypes[ws.draw()];
- size_t hapSz = haplotype.size();
- size_t j = 0;
- for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc, ++j)
- it->setAllele(TO_ALLELE(haplotype[j % hapSz]), *loc, static_cast<int>(*p));
+ for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc)
+ it->setAllele(TO_ALLELE(ws.draw()), *loc, static_cast<int>(*p));
}
+ }
}
- }
- } else {
- // openMP is disabled here because of unknown cause of initialization error for binary module
- {
+ // if user provides haplotypes
+ } else {
+ WeightedSampler ws;
IndIterator it = pop.indIterator(sp->subPop());
// m_freq can be empty if ....
- WeightedSampler ws(m_freq.empty() ? vectorf(m_haplotypes.size(), 1. / m_haplotypes.size()) : m_freq);
- for (; it.valid(); ++it)
+ if (!m_freq.valid()) {
+ vectorf freq = vectorf(m_haplotypes.size(), 1. / m_haplotypes.size());
+ ws.set(freq.begin(), freq.end());
+ } else if (m_freq.elems().empty()) {
+ if (by_locus)
+ throw ValueError("Return value of user provided function cannot have parameter loc when haplotypes are provided.");
+ vectorf freq = getFreqOrProp(0, *sp);
+ ws.set(freq.begin(), freq.end());
+ } else
+ ws.set(m_freq.elems().begin(), m_freq.elems().end());
+ //
+ for (; it.valid(); ++it) {
for (vectoru::iterator p = ploidy.begin(); p != ploidy.end(); ++p) {
- if (m_haplotypes.empty()) {
- for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc)
- it->setAllele(TO_ALLELE(ws.draw()), *loc, static_cast<int>(*p));
- } else {
- const vectori & haplotype = m_haplotypes[ws.draw()];
- size_t hapSz = haplotype.size();
- size_t j = 0;
- for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc, ++j)
- it->setAllele(TO_ALLELE(haplotype[j % hapSz]), *loc, static_cast<int>(*p));
- }
+ const vectori & haplotype = m_haplotypes[ws.draw()];
+ size_t hapSz = haplotype.size();
+ size_t j = 0;
+ for (vectoru::const_iterator loc = loci.begin(); loc != loci.end(); ++loc, ++j)
+ it->setAllele(TO_ALLELE(haplotype[j % hapSz]), *loc, static_cast<int>(*p));
}
+ }
}
}
pop.deactivateVirtualSubPop(sp->subPop());
Modified: trunk/src/initializer.h
===================================================================
--- trunk/src/initializer.h 2016-01-18 22:42:45 UTC (rev 5000)
+++ trunk/src/initializer.h 2016-01-20 04:40:18 UTC (rev 5001)
@@ -175,9 +175,16 @@
* you can use parameter \e prop to specified the exact proportions of
* alleles \c 0, \c 1, ..., although alleles with small proportions
* might not be assigned at all. Values of parameter \e prob or \e prop
- * should add up to 1. If parameter \e haplotypes is specified, it should
- * contain a list of haplotypes and parameter \e prob or \e prop specifies
- * frequencies or proportions of each haplotype. If \e loci, \e ploidy
+ * should add up to 1. In addition to a vector, parameter \e prob and
+ * \e prop can also be a function that accepts optional parameters
+ * \e loc, \e subPop or \e vsp and returns a list of requencies for
+ * alleles \c 0, \c 1, etc, or a number for frequency of allele \c 0
+ * as a speciail case for each locus, subpopulation (parameter \e subPop),
+ * or virtual subpopulations (parameter \e vsp, pass as a tuple). If
+ * parameter \e haplotypes is specified, it should contain a list of
+ * haplotypes and parameter \e prob or \e prop specifies frequencies or
+ * proportions of each haplotype (possibly diferently for each subpopulation
+ * but not each locus if the function form is used). If \e loci, \e ploidy
* and/or \e subPop are specified, only specified loci, ploidy, and
* individuals in these (virtual) subpopulations will be initialized.
* Parameter \e loci can be a list of loci indexes, names or \c ALL_AVAIL.
@@ -191,8 +198,9 @@
* similar to function \c Population.setGenotype() except that you can
* limit the initialization to certain \e loci and \e ploidy.
*/
- InitGenotype(const vectorf & freq = vectorf(),
- const uintList & genotype = vectoru(), const vectorf & prop = vectorf(),
+ InitGenotype(const floatListFunc & freq = vectorf(),
+ const uintList & genotype = vectoru(),
+ const floatListFunc & prop = vectorf(),
const intMatrix & haplotypes = intMatrix(),
const lociList & loci = lociList(),
const uintList & ploidy = uintList(),
@@ -205,7 +213,6 @@
{
}
-
/// HIDDEN Deep copy of the operator \c InitGenotype
virtual BaseOperator * clone() const
{
@@ -221,10 +228,13 @@
bool apply(Population & pop) const;
private:
+ vectorf getFreqOrProp(size_t loc, const vspID & vsp) const;
+
+private:
/// allele frequencies (assume all loci are the same for a subPop
- const vectorf m_freq;
+ const floatListFunc m_freq;
const vectoru m_genotype;
- const vectorf m_prop;
+ const floatListFunc m_prop;
const matrixi m_haplotypes;
//
const lociList m_loci;
Modified: trunk/src/simuPOP_doc.i
===================================================================
--- trunk/src/simuPOP_doc.i 2016-01-18 22:42:45 UTC (rev 5000)
+++ trunk/src/simuPOP_doc.i 2016-01-20 04:40:18 UTC (rev 5001)
@@ -1878,6 +1878,8 @@
%ignore simuPOP::floatListFunc::floatListFunc(double val);
+%ignore simuPOP::floatListFunc::valid() const;
+
%ignore simuPOP::floatListFunc::empty() const;
%ignore simuPOP::floatListFunc::size() const;
@@ -3724,20 +3726,26 @@
can use parameter prop to specified the exact proportions of
alleles 0, 1, ..., although alleles with small proportions might
not be assigned at all. Values of parameter prob or prop should
- add up to 1. If parameter haplotypes is specified, it should
- contain a list of haplotypes and parameter prob or prop specifies
- frequencies or proportions of each haplotype. If loci, ploidy
- and/or subPop are specified, only specified loci, ploidy, and
- individuals in these (virtual) subpopulations will be initialized.
- Parameter loci can be a list of loci indexes, names or ALL_AVAIL.
- If the length of a haplotype is not enough to fill all loci, the
- haplotype will be reused. If a list (or a single) haplotypes are
- specified without freq or prop, they are used with equal
- probability. In the last case, if a sequence of genotype is
- specified, it will be used repeatedly to initialize all alleles
- sequentially. This works similar to function
- Population.setGenotype() except that you can limit the
- initialization to certain loci and ploidy.
+ add up to 1. In addition to a vector, parameter prob and prop can
+ also be a function that accepts optional parameters loc, subPop or
+ vsp and returns a list of requencies for alleles 0, 1, etc, or a
+ number for frequency of allele 0 as a speciail case for each
+ locus, subpopulation (parameter subPop), or virtual subpopulations
+ (parameter vsp, pass as a tuple). If parameter haplotypes is
+ specified, it should contain a list of haplotypes and parameter
+ prob or prop specifies frequencies or proportions of each
+ haplotype (possibly diferently for each subpopulation but not each
+ locus if the function form is used). If loci, ploidy and/or subPop
+ are specified, only specified loci, ploidy, and individuals in
+ these (virtual) subpopulations will be initialized. Parameter loci
+ can be a list of loci indexes, names or ALL_AVAIL. If the length
+ of a haplotype is not enough to fill all loci, the haplotype will
+ be reused. If a list (or a single) haplotypes are specified
+ without freq or prop, they are used with equal probability. In
+ the last case, if a sequence of genotype is specified, it will be
+ used repeatedly to initialize all alleles sequentially. This works
+ similar to function Population.setGenotype() except that you can
+ limit the initialization to certain loci and ploidy.
";
@@ -7431,6 +7439,14 @@
";
+%feature("docstring") simuPOP::pyFunc::hasArg "
+
+Usage:
+
+ x.hasArg(arg)
+
+";
+
%feature("docstring") simuPOP::pyFunc::name "
Usage:
Modified: trunk/src/utility.h
===================================================================
--- trunk/src/utility.h 2016-01-18 22:42:45 UTC (rev 5000)
+++ trunk/src/utility.h 2016-01-20 04:40:18 UTC (rev 5001)
@@ -329,6 +329,11 @@
}
+ bool hasArg(const string & arg) const
+ {
+ return std::find(m_args.begin(), m_args.end(), arg) != m_args.end();
+ }
+
string name() const
{
return m_name;
@@ -1032,6 +1037,11 @@
return m_elems[i];
}
+ /// CPPONLY
+ bool valid() const
+ {
+ return !m_elems.empty() || m_func.func() != NULL;
+ }
/// CPPONLY
bool empty() const
Modified: trunk/test/test_06_initialization.py
===================================================================
--- trunk/test/test_06_initialization.py 2016-01-18 22:42:45 UTC (rev 5000)
+++ trunk/test/test_06_initialization.py 2016-01-20 04:40:18 UTC (rev 5001)
@@ -248,6 +248,95 @@
pop.subPopSize([2,1])), loci=[2,4,5], subPop=[[0,0], [2,1]])
self.assertGenotype(pop, 0, loci=[0,1,3,6,7])
+ def testInitGenotypeByFunc(self):
+ 'Test initialize genotype by function'
+ pop = Population(size=[500, 1000, 500], loci=[2,4,2])
+ # initialize all
+ initGenotype(pop, freq=lambda : [.2, .3, .5])
+ self.assertGenotypeFreq(pop, [.15, .25, .45],
+ [.25, .35, .55])
+ #
+ self.clearGenotype(pop)
+ #pop.dumpArray()
+ initGenotype(pop, freq=lambda : [.2, .3, .4, .1], loci=[2,4,6])
+ self.assertGenotypeFreq(pop, [.15, .25, .35, .05],
+ [.25, .35, .45, .15], loci=[2,4,6])
+ self.assertGenotype(pop, 0, loci=[0,1,3,5,7])
+ #
+ self.clearGenotype(pop)
+ #
+ initGenotype(pop, freq=lambda subPop: [.2,.8] if subPop == 0 else ([.8,.2] if subPop==1 else [.5,.5]))
+ self.assertGenotypeFreq(pop, [.15, .75], [.25, .85], subPop=[0])
+ self.assertGenotypeFreq(pop, [.75, .15], [.85, .25], subPop=[1])
+ self.assertGenotypeFreq(pop, [.45, .45], [.55, .55], subPop=[2])
+ #
+ # ploidy in InitByFreq'
+ pop = Population(size=[500, 1000, 500], loci=[2,4,2])
+ self.clearGenotype(pop)
+ initGenotype(pop, freq=lambda:[.2, .3, .5], loci=[2,4,6], ploidy=[0])
+ self.assertGenotypeFreq(pop, [.15, .25, .45], [.25, .35, .55],
+ loci=[2,4,6], atPloidy=0)
+ self.assertGenotype(pop, 0, loci=[0,3,5,7])
+ self.assertGenotype(pop, 0, atPloidy=1)
+ # virtual subPop
+ pop = Population(size=[5000, 10000, 5000], loci=[2,4,2], infoFields=['x'])
+ for ind in pop.individuals():
+ ind.setInfo(random.randint(10, 20), 'x')
+ pop.setVirtualSplitter(InfoSplitter('x', values=list(range(10, 15))))
+ initGenotype(pop, freq=
+ lambda subPop: [[0.2, 0.3, 0.5], [0.2,0.8], [0.5, 0.5]][subPop],
+ subPops=[[0,0],[1,1], [2,0]], loci=[2, 4, 6])
+ self.assertGenotypeFreq(pop, [.15, .25, .45], [.25, .35, .55],
+ subPop=[[0,0]], loci=[2,4,6])
+ self.assertGenotypeFreq(pop, [.15, .75], [.25, .85],
+ subPop=[[1,1]], loci=[2,4,6])
+ self.assertGenotypeFreq(pop, [.45, .45], [.55, .55],
+ subPop=[[2,0]], loci=[2,4,6])
+ # corner case
+ self.clearGenotype(pop)
+ initGenotype(pop,
+ freq=lambda subPop: [[0, 0, 1], [0, 1], [1]][subPop],
+ subPops=[[0,0],[1,1], [2,1]])
+ for ind in pop.individuals([0,0]):
+ if moduleInfo()['alleleType'] == 'binary':
+ for allele in ind.genotype():
+ self.assertEqual(allele, 1)
+ else:
+ for allele in ind.genotype():
+ self.assertEqual(allele, 2)
+ for ind in pop.individuals([1,1]):
+ for allele in ind.genotype():
+ self.assertEqual(allele, 1)
+ for ind in pop.individuals([2,1]):
+ for allele in ind.genotype():
+ self.assertEqual(allele, 0)
+ self.clearGenotype(pop)
+ self.assertRaises(ValueError,initGenotype, pop, freq=[-1,2])
+ # passing by VSP
+ self.clearGenotype(pop)
+ #def ff(vsp):
+ # print(vsp)
+ # return [[0, 0, 1], [0, 1], [1]][vsp[0]]
+ initGenotype(pop,
+ freq=lambda vsp: [[0, 0, 1], [0, 1], [1]][vsp[0]],
+ subPops=[[0,0],[1,1], [2,1]])
+ for ind in pop.individuals([0,0]):
+ if moduleInfo()['alleleType'] == 'binary':
+ for allele in ind.genotype():
+ self.assertEqual(allele, 1)
+ else:
+ for allele in ind.genotype():
+ self.assertEqual(allele, 2)
+ for ind in pop.individuals([1,1]):
+ for allele in ind.genotype():
+ self.assertEqual(allele, 1)
+ for ind in pop.individuals([2,1]):
+ for allele in ind.genotype():
+ self.assertEqual(allele, 0)
+ self.clearGenotype(pop)
+ self.assertRaises(ValueError,initGenotype, pop, freq=[-1,2])
+
+
def testInitByHaplotypes(self):
'Testing initialization by haplotypes (operator InitGenotype)'
pop = Population(size=[500, 1000, 500], loci=[2,4,2], infoFields=['x'])
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-18 22:42:48
|
Revision: 5000
http://sourceforge.net/p/simupop/code/5000
Author: bpeng2000
Date: 2016-01-18 22:42:45 +0000 (Mon, 18 Jan 2016)
Log Message:
-----------
Further improve Population.save so that it does not fail with unpickleable local variables.
Modified Paths:
--------------
trunk/src/utility.cpp
trunk/test/test_02_population.py
Modified: trunk/src/utility.cpp
===================================================================
--- trunk/src/utility.cpp 2016-01-16 17:58:23 UTC (rev 4999)
+++ trunk/src/utility.cpp 2016-01-18 22:42:45 UTC (rev 5000)
@@ -2406,21 +2406,40 @@
if (! pickle)
throw RuntimeError("Failed to import module pickle to serialize population variables.");
- // Some items in __builtins__ are not pickable so we will have to remove them
- // before picking.
- PyObject * builtins = PyString_FromString("__builtins__");
- if (PyDict_Contains(m_dict, builtins))
- PyDict_DelItem(m_dict, builtins);
- Py_DECREF(builtins);
-
// here we use version 2 because this is the latest version that supported by
// both python 2 and python 3, also because it is the one that handles simuPOP's
// defdict type using its __reduce__ interface.
PyObject * pres = PyObject_CallMethod(pickle, "dumps", "(Oi)", m_dict, 2);
if (pres == NULL) {
- PyErr_Print();
PyErr_Clear();
- throw RuntimeError("Failed to call pickle.dumps to save population variables.");
+ /* If the dictionary is not pickleable, we have to go a longer way and test if
+ each item of the dictionary is pickleable. */
+ PyObject * new_dict = PyDict_Copy(m_dict);
+ // Some items in __builtins__ are not pickable so we will have to remove them
+ // before picking.
+ PyObject *key, *value;
+ Py_ssize_t pos = 0;
+ while (PyDict_Next(m_dict, &pos, &key, &value)) {
+ // try to pickle each item
+ // if the key is not pickleable, remove the key from the copied dictionary
+ if (PyObject_CallMethod(pickle, "dumps", "(Oi)", key, 2) == NULL) {
+ PyErr_Clear();
+ PyDict_DelItem(new_dict, key);
+ }
+ if (PyObject_CallMethod(pickle, "dumps", "(Oi)", value, 2) == NULL) {
+ PyErr_Clear();
+ PyDict_DelItem(new_dict, key);
+ }
+ }
+ // now try to pickale the whole new object, in any case, the original
+ // dictionary is not touched.
+ pres = PyObject_CallMethod(pickle, "dumps", "(Oi)", new_dict, 2);
+ Py_DECREF(new_dict);
+ if (pres == NULL) {
+ PyErr_Print();
+ PyErr_Clear();
+ throw RuntimeError("Failed to call pickle.dumps to save population variables.");
+ }
}
Py_ssize_t sz = 0;
char * buf = NULL;
Modified: trunk/test/test_02_population.py
===================================================================
--- trunk/test/test_02_population.py 2016-01-16 17:58:23 UTC (rev 4999)
+++ trunk/test/test_02_population.py 2016-01-18 22:42:45 UTC (rev 5000)
@@ -1292,6 +1292,14 @@
pop1 = loadPopulation("popout")
self.assertEqual(a, pop1.dvars().alleleFreq[0][1])
self.assertEqual(pop, pop1)
+ #
+ # testing the save of a population with non-pickleable objects
+ pop.dvars().module_os = os
+ self.assertTrue('module_os' in pop.vars())
+ # module_os is not saved because it is a module
+ pop.save('popout')
+ pop1 = loadPopulation('popout')
+ self.assertFalse('module_os' in pop1.vars())
os.remove('popout')
def testCrossPlatformLoad(self):
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-16 17:58:25
|
Revision: 4999
http://sourceforge.net/p/simupop/code/4999
Author: bpeng2000
Date: 2016-01-16 17:58:23 +0000 (Sat, 16 Jan 2016)
Log Message:
-----------
Finalize documentation for the addition of operator BackwardMigrator
Modified Paths:
--------------
trunk/doc/refManual.lyx
trunk/doc/refManual.pdf
trunk/doc/userGuide.lyx
trunk/doc/userGuide.pdf
trunk/doc/userGuide.py
trunk/src/__init__.py
trunk/src/simuPOP_doc.i
trunk/src/utils.py
Modified: trunk/doc/refManual.lyx
===================================================================
--- trunk/doc/refManual.lyx 2016-01-16 17:31:00 UTC (rev 4998)
+++ trunk/doc/refManual.lyx 2016-01-16 17:58:23 UTC (rev 4999)
@@ -23,6 +23,7 @@
\chapter*{Front Matter\label{front}}
\fi
+\usepackage{listings}
\renewcommand{\lstlistlistingname}{List of Examples}
\renewcommand{\lstlistingname}{Example}
@@ -979,6 +980,13 @@
\backslash
+BackwardMigratorRef
+\end_layout
+
+\begin_layout Plain Layout
+
+
+\backslash
SplitSubPopsRef
\end_layout
@@ -1638,6 +1646,13 @@
\backslash
+simuPOPbackwardMigrateRef
+\end_layout
+
+\begin_layout Plain Layout
+
+
+\backslash
simuPOPmixedMutateRef
\end_layout
Modified: trunk/doc/refManual.pdf
===================================================================
(Binary files differ)
Modified: trunk/doc/userGuide.lyx
===================================================================
--- trunk/doc/userGuide.lyx 2016-01-16 17:31:00 UTC (rev 4998)
+++ trunk/doc/userGuide.lyx 2016-01-16 17:58:23 UTC (rev 4999)
@@ -24,6 +24,7 @@
\chapter*{Front Matter\label{front}}
\fi
+\usepackage{listings}
\renewcommand{\lstlistlistingname}{List of Examples}
\renewcommand{\lstlistingname}{Example}
@@ -8960,13 +8961,9 @@
\end_inset
(equal subpopulation size) so simuPOP will use
-\series bold
-
\begin_inset Formula $B$
\end_inset
-
-\series default
directly in this case.
Also note that
\begin_inset Formula $B$
@@ -8976,7 +8973,7 @@
\begin_inset Formula $S'$
\end_inset
-and
+ and
\begin_inset Formula $F$
\end_inset
@@ -27172,31 +27169,12 @@
\end_layout
\begin_layout Standard
-This module uses either rpy or matplotlib as the underlying plotting library.
- rpy is used by defaut because it is easier to use and support more plotting
- options (
-\family typewriter
-InfoPlotter
-\family default
-,
-\family typewriter
-BoxPlotter
-\family default
-,
-\family typewriter
-HistPlotter
-\family default
- and
-\family typewriter
-QQPlotter
-\family default
- are only available for rpy).
- matplotlib will be used if rpy is not available, or if matplotlib is set
- to be the system plotter (
-\family typewriter
-simuOpt.setOptions(plotter='matplotlib')
-\family default
-).
+simuPOP version 1.1.6 and earlier supports both rpy and matplotlib as the
+ underlying plotting library.
+ However, because of bugs in rpy2 and difficulties in supporting rpy, rpy2
+ and matplotlib, rpy/rpy2 support is removed in simuPOP 1.1.7.
+ Please use simuPOP 1.1.6 if you are interested in using rpy2.
+
\end_layout
\begin_layout Standard
@@ -27958,270 +27936,6 @@
\end_layout
-\begin_layout Subsection
-Histograms, QQ plots and boxplots (operator
-\family typewriter
-plotter.InfoPlotter
-\family default
-,
-\family typewriter
-plotter.HistPlotter
-\family default
-,
-\family typewriter
-plotter.QQPlotter
-\family default
-).
-\end_layout
-
-\begin_layout Standard
-Class
-\family typewriter
-plotter.InfoPlotter
-\family default
- can be used to draw figures for information fields of individuals in one
- or more subpopulations.
-
-\family typewriter
-plotter.HistPlotter
-\family default
- and
-\family typewriter
-plotter.QQPlotter
-\family default
- are two special cases of this class that uses functions
-\family typewriter
-hist
-\family default
- and
-\family typewriter
-qqnorm
-\family default
- respectively.
- Although an
-\family typewriter
-InfoPlotter
-\family default
- using a
-\family typewriter
-boxplot
-\family default
- function could be used, a specialized
-\family typewriter
-plotter.BoxPlotter
-\family default
- is defined so that mutliple boxplot whiskers could be drawn in the same
- plot.
-\end_layout
-
-\begin_layout Standard
-Using the same evolutionary process as Example
-\begin_inset CommandInset ref
-LatexCommand ref
-reference "ScatterPlotter"
-
-\end_inset
-
-, Example
-\begin_inset CommandInset ref
-LatexCommand ref
-reference "HistPlotter"
-
-\end_inset
-
- uses a
-\family typewriter
-HistPlotter
-\family default
- to plot the histograms (Figure
-\begin_inset CommandInset ref
-LatexCommand ref
-reference "fig:HistPlotter"
-
-\end_inset
-
-), a
-\family typewriter
-QQPlotter
-\family default
- to plot QQ plot (Figure
-\begin_inset CommandInset ref
-LatexCommand ref
-reference "fig:QQPlotter"
-
-\end_inset
-
-), and a
-\family typewriter
-BoxPlotter
-\family default
- to plot the boxplots (Figure
-\begin_inset CommandInset ref
-LatexCommand ref
-reference "fig:BoxPlotter"
-
-\end_inset
-
-) of individual ancestry values.
- By defining two virtual subpopulations by sex, the
-\family typewriter
-HistPlotter
-\family default
- and
-\family typewriter
-QQPlotter
-\family default
- plots two histograms and two QQ plots, one for males and one for females.
- Different colors are used for these figures.
- Note that these plots use the special expression value for parameter
-\family typewriter
-main
-\family default
- so that generation number can appear in the titles.
- The same technique is used in the
-\family typewriter
-dev_print_file
-\family default
- parameter of the
-\family typewriter
-BoxPlotter
-\family default
-, which overrides the default filename derived from parameter
-\family typewriter
-saveAs
-\family default
-.
-
-\begin_inset CommandInset include
-LatexCommand lstinputlisting
-filename "log/HistPlotter.log"
-lstparams "caption={Use HistPlotter to plot the histogram of individual ancestries.},label={HistPlotter}"
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Standard
-\begin_inset Float figure
-wide false
-sideways false
-status open
-
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:HistPlotter"
-
-\end_inset
-
-Histogram of individual ancestry values
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/HistPlotter_2_0.png
- lyxscale 40
- width 80text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Standard
-\begin_inset Float figure
-wide false
-sideways false
-status open
-
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:QQPlotter"
-
-\end_inset
-
-QQ plot of individual ancestry values
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/QQPlotter_2_0.png
- lyxscale 40
- width 80text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\begin_inset Float figure
-wide false
-sideways false
-status open
-
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:BoxPlotter"
-
-\end_inset
-
-Boxplot of individual ancestry values
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/Gen2.png
- lyxscale 40
- width 80text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
\begin_layout Section
Module
\family typewriter
@@ -28555,14 +28269,6 @@
demonstrates how to use function sampling.indexToID to prepare a pedigree
and how to use sampling.DrawPedigree to plot it.
-\begin_inset CommandInset include
-LatexCommand lstinputlisting
-filename "log/plotPedigree.log"
-lstparams "caption={Prepare and plot Pedigrees},label={plotPedigree}"
-
-\end_inset
-
-
\end_layout
\begin_layout Standard
@@ -28579,47 +28285,6 @@
it is common that one parent would have multiple spouses.
\end_layout
-\begin_layout Standard
-\begin_inset Float figure
-wide false
-sideways false
-status open
-
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:Pedigree"
-
-\end_inset
-
-Structure of a small complete pedigree
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/pedigree.png
- lyxscale 40
- width 60text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
\begin_layout Subsection
Sampling affected sibpairs (class
\family typewriter
@@ -28684,47 +28349,6 @@
\end_layout
-\begin_layout Standard
-\begin_inset Float figure
-wide false
-sideways false
-status open
-
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:affectedSibpair"
-
-\end_inset
-
-Affected sibpairs drawn from a small pedigree
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/affectedSibpair.png
- lyxscale 40
- width 60text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
\begin_layout Subsection
Sampling nuclear families (class
\family typewriter
@@ -28791,47 +28415,6 @@
\end_layout
-\begin_layout Standard
-\begin_inset Float figure
-wide false
-sideways false
-status open
-
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:nuclearFamily"
-
-\end_inset
-
-Nuclear families drawn from a small pedigree
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/nuclerFamily.png
- lyxscale 40
- width 60text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
\begin_layout Subsection
Sampling three-generation families (class
\family typewriter
@@ -28910,47 +28493,6 @@
\end_layout
-\begin_layout Standard
-\begin_inset Float figure
-wide false
-sideways false
-status open
-
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:threeGenFamily"
-
-\end_inset
-
-Three generation families drawn from a small pedigree
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/threeGenFamily.png
- lyxscale 40
- width 60text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
\begin_layout Subsection
Sampling different types of samples (class
\family typewriter
@@ -29022,46 +28564,9 @@
\end_layout
\begin_layout Standard
-\begin_inset Float figure
-wide false
-sideways false
-status open
-\begin_layout Plain Layout
-\begin_inset Caption Standard
-
-\begin_layout Plain Layout
-\begin_inset CommandInset label
-LatexCommand label
-name "fig:combinedSampling"
-
-\end_inset
-
-An affected sibpair and a nuclear family drawn from a small pedigree
\end_layout
-\end_inset
-
-
-\end_layout
-
-\begin_layout Plain Layout
-\align center
-\begin_inset Graphics
- filename log/combinedSampling.png
- lyxscale 40
- width 60text%
-
-\end_inset
-
-
-\end_layout
-
-\end_inset
-
-
-\end_layout
-
\begin_layout Subsection
Sampling from subpopulations and virtual subpopulations *
\end_layout
@@ -30377,23 +29882,7 @@
simuCDCV.py -h
\family default
.
- The help message is displayed in Example
-\begin_inset CommandInset ref
-LatexCommand ref
-reference "reichHelp"
-
-\end_inset
-
-.
-\begin_inset CommandInset include
-LatexCommand lstinputlisting
-filename "log/simuCDCV.hlp"
-lstparams "caption={Help information for the simuCDCV script},label={reichHelp}"
-
-\end_inset
-
-
\end_layout
\begin_layout Enumerate
Modified: trunk/doc/userGuide.pdf
===================================================================
(Binary files differ)
Modified: trunk/doc/userGuide.py
===================================================================
--- trunk/doc/userGuide.py 2016-01-16 17:31:00 UTC (rev 4998)
+++ trunk/doc/userGuide.py 2016-01-16 17:58:23 UTC (rev 4999)
@@ -4840,10 +4840,7 @@
traj = simulateForwardTrajectory(N=[2000, 4000], fitness=[1, 0.99, 0.98],
beginGen=0, endGen=100, beginFreq=[0.2, 0.3],
endFreq=[[0.1, 0.11], [0.2, 0.21]])
-# rpy syntax
-#traj.plot('log/forwardTrajectory.png', plot_ylim=[0, 0.5], col_sp=['red', 'blue'],
-# plot_main='Simulated Trajectory (forward-time)')
-# matplotlib syntax
+#
traj.plot('log/forwardTrajectory.png', set_ylim_top=0.5,
plot_c_sp=['r', 'b'], set_title_label='Simulated Trajectory (forward-time)')
pop = sim.Population(size=[2000, 4000], loci=10, infoFields='fitness')
@@ -5076,11 +5073,11 @@
ExponentialGrowthModel(r=0.01, NT=[2000, 4000]),
AdmixtureModel(model=('HI', 0, 1, 0.8, 'admixed'), T=10)
]).plot('log/MultiStage.png')
-#OutOfAfricaModel(10000).plot('log/OutOfAfrica.png')
+OutOfAfricaModel(10000).plot('log/OutOfAfrica.png')
#OutOfAfricaModel(10000, scale=10).plot('log/ScaledOutOfAfrica.png')
-#SettlementOfNewWorldModel(10000).plot('log/SettlementOfNewWorld.png')
+SettlementOfNewWorldModel(10000).plot('log/SettlementOfNewWorld.png')
#SettlementOfNewWorldModel(10000, scale=10).plot('log/ScaledSettlementOfNewWorld.png')
-#CosiModel(20000).plot('log/Cosi.png')
+CosiModel(20000).plot('log/Cosi.png')
#CosiModel(20000, scale=10).plot('log/ScaledCosi.png')
#end_ignore
#end_file
@@ -5209,6 +5206,39 @@
)
#end_file
+#begin_file log/varPlotter.py
+#begin_ignore
+import simuOpt
+simuOpt.setOptions(quiet=True)
+#end_ignore
+import simuPOP as sim
+#begin_ignore
+sim.setRNG(seed=12345)
+#end_ignore
+from simuPOP.plotter import VarPlotter
+pop = sim.Population(size=1000, loci=2)
+simu = sim.Simulator(pop, rep=3)
+simu.evolve(
+ initOps=[
+ sim.InitSex(),
+ sim.InitGenotype(genotype=[1, 2, 2, 1])
+ ],
+ matingScheme=sim.RandomMating(ops=sim.Recombinator(rates=0.01)),
+ postOps=[
+ sim.Stat(LD=[0, 1]),
+ #
+ VarPlotter('LD[0][1]', step=5, update=40, saveAs='log/varplot.png',
+ legend=['Replicate %d' % x for x in range(3)],
+ set_ylabel_ylabel='LD between marker 1 and 2',
+ set_title_label='LD decay',
+ set_ylim_bottom=0, set_ylim_top=0.25,
+ plot_linestyle_rep=['-', ':', '-.'],
+ ),
+ ],
+ gen=100
+)
+#end_file
+
#begin_file log/varPlotByDim.py
#begin_ignore
import simuOpt
Modified: trunk/src/__init__.py
===================================================================
--- trunk/src/__init__.py 2016-01-16 17:31:00 UTC (rev 4998)
+++ trunk/src/__init__.py 2016-01-16 17:58:23 UTC (rev 4999)
@@ -358,7 +358,7 @@
if not simuOptions['Quiet']:
info = moduleInfo()
- print("simuPOP Version %s : Copyright (c) 2004-2011 Bo Peng" % (__version__))
+ print("simuPOP Version %s : Copyright (c) 2004-2016 Bo Peng" % (__version__))
# compile date, compiler etc are macros that are replaced during compile time.
print("Revision %d (%s) for Python %s (%dbit, %d%s)" % \
(info['revision'], info['date'], info['python'], info['wordsize'], info['threads'],
Modified: trunk/src/simuPOP_doc.i
===================================================================
--- trunk/src/simuPOP_doc.i 2016-01-16 17:31:00 UTC (rev 4998)
+++ trunk/src/simuPOP_doc.i 2016-01-16 17:58:23 UTC (rev 4999)
@@ -55,6 +55,93 @@
";
+%feature("docstring") simuPOP::BackwardMigrator "
+
+Details:
+
+ This operator migrates individuals between all available or
+ specified subpopulations, according to a backward migration
+ matrix. It differs from Migrator in how migration matrixes are
+ interpreted. Due to the limit of this model, this operator does
+ not support migration by information field, migration by count
+ (mode = BY_COUNT), migration from virtual subpopulations,
+ migration between different number of subpopulations, and the
+ creation of new subpopulation, as operator Migrator provides. In
+ contrast to a forward migration matrix where $m_{ij}$ is
+ considered the probability (proportion or count) of individuals
+ migrating from subpopulation i to j, elements in a reverse
+ migration matrix $m_{ij}$ is considered the probability
+ (proportion or count) of individuals migrating from subpopulation
+ j to i, namely the probability (proportion or count) of
+ individuals originats from subpopulation j. If migration is
+ applied by probability, the row of the migration matrix
+ corresponding to a destination subpopulation is intepreted as
+ probabilities to orignate from each source subpopulation. Each
+ individual's source subpopulation is assigned randomly according
+ to these probabilities. Note that the probability of originating
+ from the present subpopulation is automatically calculated so the
+ corresponding matrix elements are ignored. If migration is
+ applied by proportion, the row of the migration matrix
+ corresponding to a destination subpopulation is intepreted as
+ proportions to originate from each source subpopulation. The
+ number of migrants from each source subpopulation is determined
+ before random indidividuals are chosen to migrate. Unlike the
+ forward migration matrix that describes how migration should be
+ performed, the backward migration matrix describes the result of
+ migration. The underlying forward migration matrix is calculated
+ at each generation and is in theory not the same across
+ generations. This operator calculates the corresponding forward
+ migration matrix from backward matrix and current population size.
+ This process is not always feasible so an error will raise if no
+ valid ending population size or forward migration matrix could be
+ determined. Please refer to the simuPOP user's guide for an
+ explanation of the theory behind forward and backward migration
+ matrices.
+
+";
+
+%feature("docstring") simuPOP::BackwardMigrator::BackwardMigrator "
+
+Usage:
+
+ BackwardMigrator(rate=[], mode=BY_PROBABILITY, begin=0, end=-1,
+ step=1, at=[], reps=ALL_AVAIL, subPops=ALL_AVAIL,
+ infoFields=\"migrate_to\")
+
+Details:
+
+ Create a BackwardMigrator that moves individuals between subPop
+ subpopulations randomly according to a backward migration matrix
+ rate. The size of the matrix should match the number of
+ subpopulations. Depending on the value of parameter mode,
+ elements in the migration matrix (rate) are interpreted as either
+ the probabilities to originate from source subpopulations (mode =
+ BY_PROBABILITY) or proportions of individuals originate from the
+ source (virtual) subpopulations (mode = BY_PROPORTION). Migration
+ by count is not supported by this operator. Please refer to
+ operator BaseOperator for a detailed explanation for all
+ parameters.
+
+";
+
+%feature("docstring") simuPOP::BackwardMigrator::~BackwardMigrator "
+
+Description:
+
+ destructor
+
+Usage:
+
+ x.~BackwardMigrator()
+
+";
+
+%feature("docstring") simuPOP::BackwardMigrator::clone "Obsolete or undocumented function."
+
+%feature("docstring") simuPOP::BackwardMigrator::apply "Obsolete or undocumented function."
+
+%feature("docstring") simuPOP::BackwardMigrator::describe "Obsolete or undocumented function."
+
%feature("docstring") simuPOP::BaseMutator "
Details:
@@ -756,7 +843,7 @@
distrubution to find the next true event. Also, for the cases of
p=0.5, random bits are generated. This class maintain a two
dimensional table: a vector of probabilities cross expected number
- of trials p1 p2 p3 p4 p5 trial 1 trial 2 ... trial N We expect
+ of trials p1 p2 p3 p4 p5 trial 1 trial 2 ... trial N We expect
that N is big (usually populaiton size) and p_i are small using
fast bernulliTrial method for fix p, we can fill up this table
very quickly column by column This class will provide easy access
@@ -3320,9 +3407,9 @@
this class implements a C++ iterator class that iterate through
individuals in a (sub)population. If allInds are true, the
visiblility of individuals will not be checked. Note that
- IndividualIterator *will* iterate through only visible
- individuals, and allInds is only provided when we know in advance
- that all individuals are visible. This is a way to obtain better
+ IndividualIteratorwill iterate through only visible individuals,
+ and allInds is only provided when we know in advance that all
+ individuals are visible. This is a way to obtain better
performance in simple cases.
";
@@ -4608,8 +4695,7 @@
them and assign affection status accordingly. ADDITIVE,
multiplicative, and a heterogeneour multi-locus model are
supported. Please refer to Neil Rish (1989) \"Linkage Strategies
- for Genetically Complex Traits\" for some analysis of these
- models.
+ for Genetically Complex Traits\" for some analysis of these models.
";
@@ -5265,11 +5351,11 @@
it is applied to this population. If stopOnKeyStroke is False
(default), it will always pause a population when it is applied,
if this parameter is set to True, the operator will pause a
- population if *any* key has been pressed. If a specific character
- is set, the operator will stop when this key has been pressed.
- This allows, for example, the use of several pause operators to
- pause different populations. After a population has been paused,
- a message will be displayed (unless prompt is set to False) and
+ population if any key has been pressed. If a specific character is
+ set, the operator will stop when this key has been pressed. This
+ allows, for example, the use of several pause operators to pause
+ different populations. After a population has been paused, a
+ message will be displayed (unless prompt is set to False) and
tells you how to proceed. You can press 's' to stop the evolution
of this population, 'S' to stop the evolution of all populations,
or 'p' to enter a Python shell. The current population will be
@@ -5771,9 +5857,8 @@
parameter ploidy. This operator is by default applied to
individuals in the first subpopulation but you can apply it to a
different or more than one (virtual) subpopulations using
- parameter *subPops* (``AllAvail`` is also accepted). Please refer
- to class BaseOperator for detailed descriptions of other
- parameters.
+ parameter subPops (AllAvail is also accepted). Please refer to
+ class BaseOperator for detailed descriptions of other parameters.
";
@@ -7035,9 +7120,9 @@
%ignore simuPOP::Population::setDict(PyObject *dict);
-%ignore simuPOP::Population::varsAsString() const;
+%ignore simuPOP::Population::varsAsString(bool use_pickle=false) const;
-%ignore simuPOP::Population::varsFromString(const string &vars);
+%ignore simuPOP::Population::varsFromString(const string &vars, bool use_pickle=false);
%feature("docstring") simuPOP::Population::evaluate "Obsolete or undocumented function."
@@ -8512,14 +8597,11 @@
Note:
- conversion tract length is usually short, and is estimated to be
- between 337 and 456 bp, with overall range between maybe 50 - 2500
- bp. This is usually not enough to convert, for example, two
- adjacent markers from the HapMap dataset.There is no recombination
- between sex chromosomes (Chromosomes X and Y), although
- recombination is possible between pesudoautosomal regions on these
- chromosomes. If such a feature is required, you will have to
- simulate the pesudoautosomal regions as separate chromosomes.
+ There is no recombination between sex chromosomes (Chromosomes X
+ and Y), although recombination is possible between pesudoautosomal
+ regions on these chromosomes. If such a feature is required, you
+ will have to simulate the pesudoautosomal regions as separate
+ chromosomes.
";
@@ -9335,14 +9417,12 @@
%ignore simuPOP::SharedVariables::asString() const;
-%feature("docstring") simuPOP::SharedVariables::fromString "
+%ignore simuPOP::SharedVariables::fromString(const string &vars);
-Usage:
+%ignore simuPOP::SharedVariables::topickle() const;
- x.fromString(vars)
+%ignore simuPOP::SharedVariables::frompickle(const string &vars);
-";
-
%ignore simuPOP::simpleStmt;
%feature("docstring") simuPOP::simpleStmt::simpleStmt "
@@ -11431,10 +11511,10 @@
Individuals without parents are assumed to be in the top-most
ancestral generation. This is the case for individuals in the top-
most ancestral generation if the file is saved by function
- ``Pedigree.save()``, and for individuals who only appear as
- another individual's parent, if the file is saved by operator
- ``PedigreeTagger``. The order at which offsprng is specified is
- not important because this function essentially creates a top-most
+ Pedigree.save(), and for individuals who only appear as another
+ individual's parent, if the file is saved by operator
+ PedigreeTagger. The order at which offsprng is specified is not
+ important because this function essentially creates a top-most
ancestral generation using IDs without parents, and creates the
next generation using offspring of these parents, and so on until
all generations are recreated. That is to say, if you have a
Modified: trunk/src/utils.py
===================================================================
--- trunk/src/utils.py 2016-01-16 17:31:00 UTC (rev 4998)
+++ trunk/src/utils.py 2016-01-16 17:58:23 UTC (rev 4999)
@@ -831,9 +831,8 @@
self.traj[gen].append([x for x in spFreq])
def plot(self, filename=None, **kwargs):
- '''Plot simulated Trajectory using ``R`` through a Python module ``rpy``,
- or ``matplotlib`` if ``rpy`` is not available. The function will return
- silently if module ``rpy`` and ``matplotlib`` cannot be imported.
+ '''Plot simulated Trajectory using ``matplotlib``. The function will
+ return silently if module ``matplotlib`` cannot be imported.
This function will use different colors to plot trajectories at
different loci. The trajectories are plotted from generation 0 to
@@ -845,98 +844,16 @@
your version of R or functions of matplotlib.
This function makes use of the derived keyword parameter feature of
- module ``plotter``. The ``rpy`` version allows prefixes ``par``,
- ``plot``, ``lines`` and ``dev_print``, and suffixes ``loc`` and ``sp``.
- For example, you could use parameter ``plot_ylim`` to reset the default
- value of ``ylim`` in R function ``plot``. The ``matplotlib`` version allows
- prefix ``figure``, ``plot``, ``set_title``, ``set_ylim``, ``set_xlabel``,
- and ``set_ylabel``, and suffixes ``loc`` and ``sp``.
- '''
+ module ``plotter``, with prefix ``figure``, ``plot``, ``set_title``,
+ ``set_ylim``, ``set_xlabel``, and ``set_ylabel``, and suffixes ``loc``
+ and ``sp``. '''
plotter = simuOptions['Plotter']
- if plotter is None:
- try:
- import rpy
- self._rpy_plot(filename, rpy2=False, **kwargs)
- except ImportError:
- try:
- import rpy2
- self._rpy_plot(filename, rpy2=True, **kwargs)
- except ImportError:
- try:
- import matplotlib
- self._mat_plot(filename, **kwargs)
- except ImportError:
- sys.stderr.write('Failed to draw figure: rpy or matplotlib is not available.')
- elif plotter == 'rpy':
- self._rpy_plot(filename, rpy2=False, **kwargs)
- elif plotter == 'rpy2':
- self._rpy_plot(filename, rpy2=True, **kwargs)
- else:
+ try:
+ import matplotlib
self._mat_plot(filename, **kwargs)
+ except ImportError:
+ sys.stderr.write('Failed to draw figure: rpy or matplotlib is not available.')
- def _rpy_plot(self, filename, rpy2, **kwargs):
- import plotter
- #
- args = plotter.DerivedArgs(
- defaultFuncs=['plot', 'lines'],
- allFuncs = ['par', 'plot', 'lines', 'dev_print'],
- suffixes = ['loc', 'sp'],
- defaultParams = {
- 'plot_ylim': [0, 1],
- 'plot_ylab': 'Allele frequency',
- 'plot_xlab': 'Generation'
- },
- **kwargs
- )
- # device
- plotter.newDevice()
- #
- gens = self.traj.keys()
- gens.sort()
- plotter.r.par(**args.getArgs('par', None))
- plotter.r.plot(gens[0], 0,
- **args.getArgs('plot', None, xlim=(min(gens), max(gens)), type='n'))
- allLines = []
- for gen in gens:
- for sp in range(len(self.traj[gen])):
- if len(allLines) <= sp:
- allLines.append([])
- lines = allLines[sp]
- if len(lines) != len(self.traj[gen][0]):
- if len(lines) != 0:
- for loc, line in enumerate(lines):
- # remove leading zero's
- for start in range(len(line)):
- if line[start] != 0:
- break
- if start == len(line) - 1:
- continue
- r.lines(x=range(beginGen - len(line) + start, gen), y=line[start:],
- **args.getArgs('lines', None, sp=sp, loc=loc))
- # the first point
- allLines[sp] = [[x] for x in self.traj[gen][sp]]
- else:
- for loc, x in enumerate(self.traj[gen][sp]):
- lines[loc].append(x)
- # plot the lines again
- for sp, lines in enumerate(allLines):
- for loc, line in enumerate(lines):
- # remove leading zero's
- for start in range(len(line)):
- if line[start] != 0:
- break
- if start == len(line) - 1:
- continue
- beginGen = gen - len(line) + start
- plotter.r.lines(x=range(beginGen, gen), y=line[start:],
- **args.getArgs('lines', None, sp=sp, loc=loc))
- #
- plotter.saveFigure(**args.getArgs('dev_print', None, file=filename))
- if rpy2:
- plotter.r['dev.off']()
- else:
- plotter.r.dev_off()
-
def _mat_plot(self, filename, **kwargs):
import matplotlib.pylab as plt
import plotter
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-16 17:31:03
|
Revision: 4998
http://sourceforge.net/p/simupop/code/4998
Author: bpeng2000
Date: 2016-01-16 17:31:00 +0000 (Sat, 16 Jan 2016)
Log Message:
-----------
Document and test BackwardMigrator
Modified Paths:
--------------
trunk/MANIFEST.in
trunk/doc/userGuide.lyx
trunk/doc/userGuide.py
trunk/src/migrator.cpp
trunk/test/test_12_migration.py
Modified: trunk/MANIFEST.in
===================================================================
--- trunk/MANIFEST.in 2016-01-16 07:13:35 UTC (rev 4997)
+++ trunk/MANIFEST.in 2016-01-16 17:31:00 UTC (rev 4998)
@@ -10,7 +10,6 @@
include simuPOP_version.py
include MANIFEST.in
include setup.py
-include SConstruct
# source header files
include config.h
Modified: trunk/doc/userGuide.lyx
===================================================================
--- trunk/doc/userGuide.lyx 2016-01-16 07:13:35 UTC (rev 4997)
+++ trunk/doc/userGuide.lyx 2016-01-16 17:31:00 UTC (rev 4998)
@@ -8748,6 +8748,291 @@
\end_layout
\begin_layout Subsection
+Migration using backward migration matrix (operator
+\family typewriter
+BackwardMigrator
+\family default
+)
+\end_layout
+
+\begin_layout Standard
+Backward migration matrices are widely used in theoretical population genetics
+ and coalescent based simulations.
+ Instead of specifying the probability of migrating from another subpopulation
+ to another (namely how migration happens), such matrices specify the probabilit
+y that individuals in a subpopulation originate from others (namely the
+ result of migration).
+ simuPOP simulates such models by converting backward migration matrices
+ to foward ones using the theory describes below.
+ Due to the limit of such models, simuPOP cannot simulate migration from/to
+ virtual subpopulatons, creation of new subpopulation, different source
+ and destination subpopulations, and will generate an error if the conversion
+ process fails.
+\end_layout
+
+\begin_layout Standard
+To explain the differences between forward and backward migration matrices,
+ let us assume that there are
+\begin_inset Formula $d$
+\end_inset
+
+ subpopulations with population sizes
+\begin_inset Formula $S=\left[S_{1},S_{2},...,S_{d}\right]$
+\end_inset
+
+, and a forward migration matrix
+\begin_inset Formula
+\[
+F=\left[\begin{array}{cccc}
+f_{11} & f_{12} & \cdots & f_{1d}\\
+f_{21} & f_{22} & \cdots & f_{2d}\\
+\vdots & & & \vdots\\
+f_{d1} & f_{d2} & \cdots & f_{dd}
+\end{array}\right]
+\]
+
+\end_inset
+
+where
+\begin_inset Formula $f_{ij}$
+\end_inset
+
+ is the probability that an individual will migration from subpopulation
+
+\begin_inset Formula $i$
+\end_inset
+
+ to
+\begin_inset Formula $j$
+\end_inset
+
+.
+ After migration happens, subppulation sizes are changed to
+\begin_inset Formula $S'=\left[S'_{1},S'_{2},...,S'_{d}\right]$
+\end_inset
+
+, and the origin of individuals in each subpopulation can be described by
+ the backward migration matrix
+\begin_inset Formula
+\[
+B=\left[\begin{array}{cccc}
+b_{11} & b_{12} & \cdots & b_{1d}\\
+b_{21} & b_{22} & \cdots & b_{2d}\\
+\vdots & & & \vdots\\
+b_{d1} & b_{d2} & \cdots & b_{dd}
+\end{array}\right]
+\]
+
+\end_inset
+
+where
+\begin_inset Formula $b_{ij}$
+\end_inset
+
+ is the probability that an individual in subpopulation
+\begin_inset Formula $i$
+\end_inset
+
+ originates from subpopulation
+\begin_inset Formula $j$
+\end_inset
+
+.
+
+\end_layout
+
+\begin_layout Standard
+These qualities can be derived from original population sizes and the forward
+ migration matrix.
+ That is to say, the size of new subpopulation
+\begin_inset Formula $k$
+\end_inset
+
+ is the sum of all migrants to this subpopulation
+\begin_inset Formula
+\[
+S'_{k}={\displaystyle \sum_{i=1}^{d}}S_{i}f_{ik}
+\]
+
+\end_inset
+
+and the size of the original population
+\begin_inset Formula $k$
+\end_inset
+
+ is the sum of all migrants from this subpopulation
+\begin_inset Formula
+\[
+S{}_{k}={\displaystyle \sum_{i=1}^{d}}S'_{i}b_{ik}
+\]
+
+\end_inset
+
+and the composition of subpopulation
+\begin_inset Formula $k$
+\end_inset
+
+ (e.g.
+ individuals originate from subpopulation
+\begin_inset Formula $j$
+\end_inset
+
+) is
+\begin_inset Formula
+\[
+b_{kj}=\dfrac{S_{j}f_{jk}}{S'_{k}}
+\]
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+In matrix form, these formulas can be written as
+\begin_inset Formula
+\[
+S'=F^{T}S
+\]
+
+\end_inset
+
+
+\begin_inset Formula
+\[
+S=B^{T}S'
+\]
+
+\end_inset
+
+and
+\end_layout
+
+\begin_layout Standard
+\begin_inset Formula
+\[
+B=diag(S')^{-1}F^{T}diag(S)
+\]
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+Therefore, given a backward migration matrix
+\begin_inset Formula $B$
+\end_inset
+
+ and current population size
+\begin_inset Formula $S$
+\end_inset
+
+, we can derive a forward migration matrix using
+\begin_inset Formula
+\[
+S'=\left(B^{T}\right)^{-1}S
+\]
+
+\end_inset
+
+and
+\end_layout
+
+\begin_layout Standard
+\begin_inset Formula
+\[
+F=diag(S)^{-1}B^{T}diag(S')
+\]
+
+\end_inset
+
+Note that
+\begin_inset Formula $F=B$
+\end_inset
+
+ is always true if
+\begin_inset Formula $B$
+\end_inset
+
+ is symmetric and
+\begin_inset Formula $S_{i}=S_{j}$
+\end_inset
+
+ (equal subpopulation size) so simuPOP will use
+\series bold
+
+\begin_inset Formula $B$
+\end_inset
+
+
+\series default
+ directly in this case.
+ Also note that
+\begin_inset Formula $B$
+\end_inset
+
+ might not be inversable and
+\begin_inset Formula $S'$
+\end_inset
+
+and
+\begin_inset Formula $F$
+\end_inset
+
+ might be invalid (e.g.
+ negative population size or forward migration rate) for given
+\begin_inset Formula $B$
+\end_inset
+
+ and
+\begin_inset Formula $S$
+\end_inset
+
+.
+ simuPOP will terminate with an error message in these cases.
+
+\end_layout
+
+\begin_layout Standard
+The following example
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "backwardMigration"
+
+\end_inset
+
+ demonstrates how to use a backward migration matrix to perform migration.
+ It initializes all individuals with indexes of subpopulations they belong
+ to before migration and calculates the percent of individuals from each
+ source population using a PyOperator with function originOfInds.
+ The so-called overseved backward migration matrix is similar to specified
+ migration matrix despite of stochastic effects.
+ This example also uses turnOnDebug function to let the operator print the
+ expected subpopulation size (
+\begin_inset Formula $S'$
+\end_inset
+
+) and calculate forward migration matrix (
+\begin_inset Formula $F$
+\end_inset
+
+) at each generation, which, as expected, vary from generation to generation.
+
+\end_layout
+
+\begin_layout Standard
+\begin_inset CommandInset include
+LatexCommand lstinputlisting
+filename "log/backwardMigrate.log"
+lstparams "caption={Migration using a backward migration matrix},label={backwardMigration}"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsection
Split subpopulations (operators
\family typewriter
SplitSubPops
Modified: trunk/doc/userGuide.py
===================================================================
--- trunk/doc/userGuide.py 2016-01-16 07:13:35 UTC (rev 4997)
+++ trunk/doc/userGuide.py 2016-01-16 17:31:00 UTC (rev 4998)
@@ -2698,6 +2698,53 @@
)
#end_file
+
+
+#begin_file log/backwardMigrate.py
+#begin_ignore
+import simuOpt
+simuOpt.setOptions(quiet=True)
+#end_ignore
+import simuPOP as sim
+#begin_ignore
+sim.setRNG(seed=12345)
+#end_ignore
+sim.turnOnDebug('DBG_MIGRATOR')
+pop = sim.Population(size=[10000, 5000, 8000], infoFields=['migrate_to', 'migrate_from'])
+def originOfInds(pop):
+ print('Observed backward migration matrix at generation {}'.format(pop.dvars().gen))
+ for sp in range(pop.numSubPop()):
+ # get source subpop for all individuals in subpopulation i
+ origins = pop.indInfo('migrate_from', sp)
+ spSize = pop.subPopSize(sp)
+ B_sp = [origins.count(j) * 1.0 /spSize for j in range(pop.numSubPop())]
+ print(' ' + ', '.join(['{:.3f}'.format(x) for x in B_sp]))
+ return True
+
+pop.evolve(
+ initOps=sim.InitSex(),
+ preOps=
+ # mark the source subpopulation of each individual
+ [sim.InitInfo(i, subPops=i, infoFields='migrate_from') for i in range(3)] + [
+ # perform migration
+ sim.BackwardMigrator(rate=[
+ [0, 0.04, 0.02],
+ [0.05, 0, 0.02],
+ [0.02, 0.01, 0]
+ ]),
+ # calculate and print observed backward migration matrix
+ sim.PyOperator(func=originOfInds),
+ # calculate population size
+ sim.Stat(popSize=True),
+ # and print it
+ sim.PyEval(r'"Pop size after migration: {}\n".format(", ".join([str(x) for x in subPopSize]))'),
+ ],
+ matingScheme=sim.RandomMating(),
+ gen = 5
+)
+#end_file
+
+
#begin_file log/splitByInfo.py
#begin_ignore
import simuOpt
Modified: trunk/src/migrator.cpp
===================================================================
--- trunk/src/migrator.cpp 2016-01-16 07:13:35 UTC (rev 4997)
+++ trunk/src/migrator.cpp 2016-01-16 17:31:00 UTC (rev 4998)
@@ -300,7 +300,6 @@
}
}
// inverse
- DBG_DO(DBG_MIGRATOR, cerr << "Transpose of backward migration matrix B is " << Bt << endl);
boost::numeric::ublas::permutation_matrix<std::size_t> pm(sz);
int res = lu_factorize(Bt, pm);
if (res != 0)
@@ -309,7 +308,9 @@
// backsubstite to get the inverse
lu_substitute(Bt, pm, m_inverse_rate);
+ /*
DBG_DO(DBG_MIGRATOR, cerr << "Inverse of B' is " << m_inverse_rate << endl);
+ */
}
Modified: trunk/test/test_12_migration.py
===================================================================
--- trunk/test/test_12_migration.py 2016-01-16 07:13:35 UTC (rev 4997)
+++ trunk/test/test_12_migration.py 2016-01-16 17:31:00 UTC (rev 4998)
@@ -23,6 +23,80 @@
sys.argv=new_argv
from simuPOP import *
+import numpy as np
+from numpy.linalg import inv
+
+def Sp_From_Forward_and_S(FM, S):
+ # calculate predicted subpopulation size from forward migration matrix (FM)
+ # and current population size (S) S'=F^T S
+ #
+ # We need to normalize FM(i,i) so that row sum is 1
+ return np.array(
+ [[1-sum(x)+k if j==i else k for j,k in enumerate(x)]
+ for i,x in enumerate(FM)]).transpose().dot(np.array(S))
+
+# for example: the following equals 2000, 4000, 4000
+# Sp_From_Forward_and_S([ [0, .05, .05],
+# [0.025, 0, 0],
+# [0.025, 0, 0] ],
+# [2000, 4000, 4000])
+
+def BM_From_Forward_and_S(FM, S):
+ # calculate backward migration matrix from forward migration matrix
+ # and current population size
+ #
+ # B = diag(Sp)^-1 F^T diag(S)
+ Sp = Sp_From_Forward_and_S(FM, S)
+ Ft = np.array([[1-sum(x)+k if j==i else k for j,k in enumerate(x)]
+ for i,x in enumerate(FM)]).transpose()
+ return np.diag(1. / Sp).dot(Ft).dot(np.diag(S))
+
+# the following get
+# [[ 0.9 , 0.05 , 0.05 ],
+# [ 0.025, 0.975, 0. ],
+# [ 0.025, 0. , 0.975]]
+#
+#BM_From_Forward_and_S([ [0, .05, .05],
+# [0.025, 0, 0],
+# [0.025, 0, 0] ],
+# [2000, 4000, 4000])
+
+def Sp_From_Backward_and_S(BM, S):
+ # Calculate predicted subpopulation size from backward migration matrix
+ # and current population size: S'=B^T^-1 * S
+ #
+ # normalize
+ B = np.array([[1-sum(x)+k if j==i else k for j,k in enumerate(x)]
+ for i,x in enumerate(BM)])
+ return inv(np.array(B.transpose())).dot(np.array(S))
+
+# the following yield [2000, 4000, 4000]
+#Sp_From_Backward_and_S([[ 0.9 , 0.05 , 0.05 ],
+# [ 0.025, 0.975, 0. ],
+# [ 0.025, 0. , 0.975]],
+# [2000, 4000, 4000])
+
+def FM_From_Backward_and_S(BM, S):
+ #
+ # F = diag(S)^-1 B^T diag(S')
+ #
+ Sp = Sp_From_Backward_and_S(BM, S)
+ Bt = np.array([[1-sum(x)+k if j==i else k for j,k in enumerate(x)]
+ for i,x in enumerate(BM)]).transpose()
+ return np.diag(1. / np.array(S)).dot(Bt).dot(np.diag(Sp))
+
+# the following yields
+# [[ 0.9 , 0.05 , 0.05 ],
+# [ 0.025, 0.975, 0. ],
+# [ 0.025, 0. , 0.975]]
+
+#FM_From_Backward_and_S([[ 0.9 , 0.05 , 0.05 ],
+# [ 0.025, 0.975, 0. ],
+# [ 0.025, 0. , 0.975]],
+# [2000, 4000, 4000])
+
+
+
class TestMigrator(unittest.TestCase):
def testmigrateByCounts(self):
@@ -44,7 +118,7 @@
[50, 0, 0],
[50, 0, 0] ], BY_COUNTS)
- def testmigrateByProportion(self):
+ def testMigrateByProportion(self):
'Testing migrate by proportion'
pop = Population(size=[2000,4000,4000], loci=[2], infoFields=['migrate_to'])
# now if we want to inject a mutation whenever fixation happens
@@ -60,6 +134,8 @@
[0, 0.25, 0] ])
self.assertEqual(pop.subPopSizes(), (2000, 4500, 3500))
+
+
def testmigrateByIndInfo(self):
'Testing migrate by indinfo'
pop = Population(size=[2000, 4000, 4000], loci=[2], infoFields=['migrate_to'])
@@ -374,7 +450,57 @@
for ind in (10, 11, 12, 13, 18, 19, 20):
self.assertEqual(pop.individual(ind).genotype(), [0]*16)
+ def testMigrateByBackwardProportion(self):
+ 'Testing migrate by proportion'
+ pop = Population(size=[2000,4000,4000], loci=[2], infoFields=['migrate_to'])
+ # now if we want to inject a mutation whenever fixation happens
+ self.assertRaises(ValueError, migrate, pop, rate=[1, 0.1])
+ backwardMigrate(pop, mode=BY_PROPORTION,
+ rate =[[ 0.9 , 0.05 , 0.05 ],
+ [ 0.025, 0.975, 0. ],
+ [ 0.025, 0. , 0.975]]
+ )
+ # is not exactly 2000, 4000, 4000 due to numeric issue.
+ self.assertEqual(pop.subPopSizes(), (1998, 4001, 4001))
+ #
+ pop = Population(size=[2000,4000,4000], loci=[2], infoFields=['migrate_to'])
+ backwardMigrate(pop, mode=BY_PROPORTION,
+ rate = [[ 0.5 , 0.5 , 0. ],
+ [ 0.11111111, 0.66666667, 0.22222222],
+ [ 0.14285714, 0. , 0.85714286]])
+ self.assertEqual(pop.subPopSizes(), (2002, 4498, 3500))
+ def testMigrateByBackwardProbability(self):
+ 'Testing migrate by probability'
+ pop = Population(size=[2000,4000,4000], loci=[2], infoFields=['migrate_to'])
+ # now if we want to inject a mutation whenever fixation happens
+ backwardMigrate(pop, mode=BY_PROBABILITY,
+ rate = [[ 0.9 , 0.05 , 0.05 ],
+ [ 0.025, 0.975, 0. ],
+ [ 0.025, 0. , 0.975]])
+ # print pop.subPopSizes()
+ self.assertTrue(abs(pop.subPopSize(0) - 2000) < 100,
+ "Expression abs(pop.subPopSize(0) - 2000) (test value %f) be less than 100. This test may occasionally fail due to the randomness of outcome." % (abs(pop.subPopSize(0) - 2000)))
+ self.assertTrue(abs(pop.subPopSize(1) - 4000) < 100,
+ "Expression abs(pop.subPopSize(1) - 4000) (test value %f) be less than 100. This test may occasionally fail due to the randomness of outcome." % (abs(pop.subPopSize(1) - 4000)))
+ self.assertTrue(abs(pop.subPopSize(2) - 4000) < 100,
+ "Expression abs(pop.subPopSize(2) - 4000) (test value %f) be less than 100. This test may occasionally fail due to the randomness of outcome." % (abs(pop.subPopSize(2) - 4000)))
+
+ pop = Population(size=[2000,4000,4000], loci=[2], infoFields=['migrate_to'])
+ backwardMigrate(pop, mode=BY_PROBABILITY,
+ rate = [[ 0.5 , 0.5 , 0. ],
+ [ 0.11111111, 0.66666667, 0.22222222],
+ [ 0.14285714, 0. , 0.85714286]])
+ # print pop.subPopSizes()
+ self.assertTrue(abs(pop.subPopSize(0) - 2000) < 100,
+ "Expression abs(pop.subPopSize(0) - 2000) (test value %f) be less than 100. This test may occasionally fail due to the randomness of outcome." % (abs(pop.subPopSize(0) - 2000)))
+ self.assertTrue(abs(pop.subPopSize(1) - 4500) < 100,
+ "Expression abs(pop.subPopSize(1) - 4500) (test value %f) be less than 100. This test may occasionally fail due to the randomness of outcome." % (abs(pop.subPopSize(1) - 4500)))
+ self.assertTrue(abs(pop.subPopSize(2) - 3500) < 100,
+ "Expression abs(pop.subPopSize(2) - 3500) (test value %f) be less than 100. This test may occasionally fail due to the randomness of outcome." % (abs(pop.subPopSize(2) - 3500)))
+
+
+
if __name__ == '__main__':
unittest.main()
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2016-01-16 07:13:38
|
Revision: 4997
http://sourceforge.net/p/simupop/code/4997
Author: bpeng2000
Date: 2016-01-16 07:13:35 +0000 (Sat, 16 Jan 2016)
Log Message:
-----------
Initial support for migration with backward migration matrix
Modified Paths:
--------------
trunk/src/__init__.py
trunk/src/migrator.cpp
trunk/src/migrator.h
Modified: trunk/src/__init__.py
===================================================================
--- trunk/src/__init__.py 2015-11-01 03:09:25 UTC (rev 4996)
+++ trunk/src/__init__.py 2016-01-16 07:13:35 UTC (rev 4997)
@@ -166,6 +166,7 @@
'InfoExec',
#
'Migrator',
+ 'BackwardMigrator',
'MergeSubPops',
'SplitSubPops',
'ResizeSubPops',
@@ -245,6 +246,7 @@
'tagID',
#
'migrate',
+ 'backwardMigrate',
'resizeSubPops',
'splitSubPops',
'mergeSubPops',
@@ -1132,6 +1134,10 @@
'Function form of operator ``Migrator``.'
Migrator(*args, **kwargs).apply(pop)
+def backwardMigrate(pop, *args, **kwargs):
+ 'Function form of operator ``BackwardMigrator``.'
+ BackwardMigrator(*args, **kwargs).apply(pop)
+
def splitSubPops(pop, *args, **kwargs):
'''Split subpopulations (*subPops*) of population *pop* according to either
*sizes* or *proportions* of the resulting subpopulations, or an information
Modified: trunk/src/migrator.cpp
===================================================================
--- trunk/src/migrator.cpp 2015-11-01 03:09:25 UTC (rev 4996)
+++ trunk/src/migrator.cpp 2016-01-16 07:13:35 UTC (rev 4997)
@@ -26,6 +26,9 @@
#include "virtualSubPop.h"
#include "migrator.h"
+#include <boost/numeric/ublas/lu.hpp>
+#include <boost/numeric/ublas/io.hpp>
+
namespace simuPOP {
Migrator::Migrator(const floatMatrix & rate, int mode, const uintList & toSubPops,
@@ -145,7 +148,6 @@
}
}
-
for (size_t from = 0, fromEnd = fromSubPops.size(); from < fromEnd; ++from) {
size_t spFrom = fromSubPops[from].subPop();
// rateSize might be toSize + 1, the last one is from->from
@@ -264,6 +266,243 @@
}
+
+BackwardMigrator::BackwardMigrator(const floatMatrix & rate, int mode,
+ int begin, int end, int step, const intList & at,
+ const intList & reps, const subPopList & subPops, const stringList & infoFields)
+ : BaseOperator("", begin, end, step, at, reps, subPops, infoFields),
+ m_rate(rate.elems()), m_inverse_rate(rate.elems().size(), rate.elems().size()), m_mode(mode), m_symmetric_matrix(true)
+{
+ DBG_FAILIF(!subPops.empty() && subPops.size() != m_rate.size(),
+ ValueError, "Length of param subPop must match rows of rate matrix.");
+ DBG_FAILIF(m_mode != BY_PROBABILITY && m_mode != BY_PROPORTION,
+ ValueError, "Only BY_PROBABILITY and BY_PROPORTION is supported by BackMigrator");
+
+ // find the inverse of B'
+ size_t sz = m_rate.size();
+ boost::numeric::ublas::matrix<double> Bt(sz, sz);
+ for (size_t i = 0; i < sz; ++i) {
+ if (m_rate[i].size() != sz)
+ throw ValueError("A m by m matrix is required for backward migration matrix.");
+ //
+ for (size_t j = 0; j < sz; ++j) {
+ if (i == j) {
+ Bt(i, i) = accumulate(m_rate[i].begin(), m_rate[i].end(), 0.0) - m_rate[i][i];
+ Bt(i, i) = 1. - Bt(i, i);
+ } else
+ Bt(i, j) = m_rate[j][i];
+ if (fcmp_lt(Bt(i, j), 0.))
+ throw ValueError((boost::format("Backward migration rate should be positive. %d provided.") % Bt(i,j)).str());
+ if (fcmp_gt(Bt(i, j), 1.))
+ throw ValueError((boost::format("Backward migration rate should be in the range of [0,1], %d provided.") % Bt(i,j)).str());
+ if (j > i && m_rate[i][j] != m_rate[j][i])
+ m_symmetric_matrix = false;
+ }
+ }
+ // inverse
+ DBG_DO(DBG_MIGRATOR, cerr << "Transpose of backward migration matrix B is " << Bt << endl);
+ boost::numeric::ublas::permutation_matrix<std::size_t> pm(sz);
+ int res = lu_factorize(Bt, pm);
+ if (res != 0)
+ throw RuntimeError("Failed to convert backward matrix to forward migration matrix. (Matrix is not inversable).");
+ m_inverse_rate.assign(boost::numeric::ublas::identity_matrix<double>(sz));
+ // backsubstite to get the inverse
+ lu_substitute(Bt, pm, m_inverse_rate);
+
+ DBG_DO(DBG_MIGRATOR, cerr << "Inverse of B' is " << m_inverse_rate << endl);
+}
+
+
+string BackwardMigrator::describe(bool /* format */) const
+{
+ return "<simuPOP.BackwardMigrator>";
+}
+
+
+bool BackwardMigrator::apply(Population & pop) const
+{
+ // set info of individual
+ size_t info = pop.infoIdx(infoField(0));
+
+ subPopList VSPs = applicableSubPops(pop);
+ if (VSPs.size() <= 1)
+ return true;
+
+ DBG_FAILIF(VSPs.size() != m_rate.size(),
+ ValueError, "Number of 'from' subpopulations should match number of rows of migration rate matrix.");
+
+ vectoru subPops;
+ for (size_t i = 0; i < VSPs.size(); ++i) {
+ DBG_FAILIF(VSPs[i].isVirtual(), ValueError,
+ "BackwardMigrator does not support virtual subpupulations.")
+ DBG_FAILIF(m_rate[i].size() != VSPs.size(), ValueError,
+ "A square matrix is required for BackwardMigrator")
+ subPops.push_back(VSPs[i].subPop());
+ }
+
+ // assign individuals their own subpopulation ID
+ for (size_t sp = 0; sp < pop.numSubPop(); ++sp) {
+ RawIndIterator it = pop.rawIndBegin(sp);
+ RawIndIterator it_end = pop.rawIndEnd(sp);
+ if (numThreads() > 1) {
+#ifdef _OPENMP
+ size_t popSize = it_end - it;
+# pragma omp parallel firstprivate(it, it_end)
+ {
+ size_t id = omp_get_thread_num();
+ it = it + id * (popSize / numThreads());
+ it_end = id == numThreads() - 1 ? it_end : it + popSize / numThreads();
+ for (; it != it_end; ++it)
+ it->setInfo(static_cast<double>(sp), info);
+ }
+#endif
+ } else {
+ for (; it != it_end; ++it)
+ it->setInfo(static_cast<double>(sp), info);
+ }
+ }
+
+ DBG_FAILIF(pop.hasActivatedVirtualSubPop(), RuntimeError,
+ "Migration can not be applied to activated virtual subpopulations");
+
+ // subpopulation size S (before migration)
+ vectoru S;
+ for (size_t i = 0; i < subPops.size(); ++i)
+ S.push_back(pop.subPopSize(subPops[i]));
+
+ // now, we need to calculate a forward migration matrix from the backward one
+ // the formula is
+
+ // S' = B'^-1 * S
+ // F = diag(S)^-1 B' diag(S')
+ //
+ // but for the special case of equal population size and symmetric matrix,
+ // two matrices are the same
+ bool simple_case = m_symmetric_matrix;
+ // symmtrix matrix, check if equal population size
+ if (simple_case) {
+ for (size_t i = 1; i < subPops.size(); ++i)
+ if (S[i] != S[i-1])
+ simple_case = false;
+ }
+
+ size_t sz = m_rate.size();
+ // if not the simple case, we need to calculate rate...
+ matrixf migrationRate;
+ if (simple_case)
+ migrationRate = m_rate;
+ else {
+ // with Bt^-1, we can calculate expected population size
+ vectorf Sp(sz);
+ for (size_t i = 0; i < sz; ++i) {
+ Sp[i] = 0;
+ for (size_t j = 0; j < sz; ++j)
+ Sp[i] += m_inverse_rate(i, j) * S[j];
+ }
+ DBG_DO(DBG_MIGRATOR, cerr << "Expected next population size is " << Sp << endl);
+ for (size_t i = 0; i < sz; ++i) {
+ if (Sp[i] <= 0)
+ throw RuntimeError((boost::format("Failed to calculate forward migration matrix: negative expected population size %f for subpopulation %d")
+ % Sp[i] % i).str());
+ }
+ // now, F = diag(S)^-1 * BT * diag (S')
+ migrationRate = m_rate;
+ for (size_t i = 0; i < sz; ++i)
+ for (size_t j = 0; j < sz; ++j)
+ // m_rate[i][i] might not be defined.
+ migrationRate[i][j] = m_rate[j][i] * Sp[j] / S[i];
+ }
+
+ // check parameters
+ for (size_t i = 0; i < sz; ++i) {
+ for (size_t j = 0; j < sz; ++j) {
+ if (fcmp_lt(migrationRate[i][j], 0.))
+ throw ValueError("Converted forward migration rate should be positive.");
+ if (fcmp_gt(migrationRate[i][j], 1.))
+ throw ValueError("Converted forward migration rate should be in the range of [0,1]");
+ }
+ // look for from=to cell.
+ double sum = accumulate(migrationRate[i].begin(), migrationRate[i].end(), 0.0);
+ //
+ double & self = migrationRate[i][i];
+ sum -= self;
+ if (fcmp_gt(sum, 1.0))
+ throw ValueError("Sum of migrate rate from one subPop should <= 1");
+ // reset to-my-self probability/proportion
+ self = 1.0 - sum;
+
+ }
+ if (! simple_case) {
+ DBG_DO(DBG_MIGRATOR, cerr << "Forward migration matrix is " << migrationRate << endl);
+ }
+
+ for (size_t from = 0, fromEnd = subPops.size(); from < fromEnd; ++from) {
+ size_t spFrom = subPops[from];
+ // rateSize might be toSize + 1, the last one is from->from
+ size_t toIndex;
+
+ size_t spSize = pop.subPopSize(spFrom);
+
+ if (m_mode == BY_PROBABILITY) {
+ WeightedSampler ws(migrationRate[from]);
+
+ // for each individual, migrate according to migration probability
+ for (IndIterator ind = pop.indIterator(spFrom); ind.valid(); ++ind) {
+ //toIndex = getRNG().randIntByFreq( rateSize, &migrationRate[from][0] ) ;
+ toIndex = ws.draw();
+
+ DBG_ASSERT(toIndex < migrationRate[from].size(), ValueError,
+ "Return index out of range.");
+
+ if (toIndex < sz && subPops[toIndex] != spFrom)
+ ind->setInfo(static_cast<double>(subPops[toIndex]), info);
+ }
+ } else {
+ // 2nd, or 3rd method
+ // first find out how many people will move to other subPop
+ // then randomly assign individuals to move
+ vectoru toNum(sz);
+ // in case that to sub is not in from sub, the last added
+ // element is not used. sum of toNum is not spSize.
+ for (size_t i = 0; i < sz; ++i)
+ toNum[i] = static_cast<ULONG>(spSize * migrationRate[from][i]);
+ // create a vector and assign indexes, then random shuffle
+ // and assign info
+ vectoru toIndices(spSize);
+ size_t k = 0;
+ for (size_t i = 0; i < sz && k < spSize; ++i)
+ for (size_t j = 0; j < toNum[i] && k < spSize; ++j)
+ toIndices[k++] = subPops[i];
+
+ // the rest of individuals will stay in their original subpopulation.
+ while (k < spSize)
+ toIndices[k++] = spFrom;
+
+ getRNG().randomShuffle(toIndices.begin(), toIndices.end());
+ IndIterator ind = pop.indIterator(spFrom);
+ // set info
+ for (size_t i = 0; ind.valid(); ++i, ++ind)
+ // The previous migration_to value, if set by a previous vsp, will be overridden.
+ ind->setInfo(static_cast<double>(toIndices[i]), info);
+ }
+ } // for all subPop.
+
+ // do migration.
+ size_t oldNumSubPop = pop.numSubPop();
+ pop.setSubPopByIndInfo(infoField(0));
+ // try to keep number of subpopulations.
+ if (pop.numSubPop() < oldNumSubPop) {
+ vectorf split(oldNumSubPop - pop.numSubPop() + 1, 0);
+ split[0] = static_cast<double>(pop.subPopSize(pop.numSubPop() - 1));
+ pop.splitSubPop(pop.numSubPop() - 1, split);
+ }
+ DBG_ASSERT(pop.numSubPop() >= oldNumSubPop, RuntimeError,
+ "Migrator should not decrease number of subpopulations.");
+
+ return true;
+}
+
+
struct compareVSP
{
int operator()(const vspID & v1, const vspID & v2)
Modified: trunk/src/migrator.h
===================================================================
--- trunk/src/migrator.h 2015-11-01 03:09:25 UTC (rev 4996)
+++ trunk/src/migrator.h 2016-01-16 07:13:35 UTC (rev 4997)
@@ -33,6 +33,8 @@
#include <list>
using std::list;
+#include <boost/numeric/ublas/matrix.hpp>
+
#include <iostream>
using std::cout;
using std::endl;
@@ -164,6 +166,100 @@
};
+/** This operator migrates individuals between all available or specified
+ * subpopulations, according to a backward migration matrix. It differs from
+ * \c Migrator in how migration matrixes are interpreted. Due to the limit
+ * of this model, this operator does not support migration by information
+ * field, migration by count (\e mode = \c BY_COUNT), migration from virtual
+ * subpopulations, migration between different number of subpopulations,
+ * and the creation of new subpopulation, as operator \c Migrator provides.
+ *
+ * In contrast to a forward migration matrix where $m_{ij}$ is considered
+ * the probability (proportion or count) of individuals migrating from subpopulation
+ * \c i to \c j, elements in a reverse migration matrix $m_{ij}$ is considered
+ * the probability (proportion or count) of individuals migrating from subpopulation
+ * \c j to \c i, namely the probability (proportion or count) of individuals
+ * originats from subpopulation \c j.
+ *
+ * If migration is applied by probability, the row of the migration matrix
+ * corresponding to a destination subpopulation is intepreted as probabilities to
+ * orignate from each source subpopulation. Each individual's source
+ * subpopulation is assigned randomly according to these probabilities. Note
+ * that the probability of originating from the present subpopulation is
+ * automatically calculated so the corresponding matrix elements are ignored.
+ *
+ * If migration is applied by proportion, the row of the migration matrix
+ * corresponding to a destination subpopulation is intepreted as proportions
+ * to originate from each source subpopulation. The number of migrants from each
+ * source subpopulation is determined before random indidividuals are
+ * chosen to migrate.
+ *
+ * Unlike the forward migration matrix that describes how migration should
+ * be performed, the backward migration matrix describes the result of
+ * migration. The underlying forward migration matrix is calculated at
+ * each generation and is in theory not the same across generations.
+ *
+ * This operator calculates the corresponding forward migration matrix
+ * from backward matrix and current population size. This process is not
+ * always feasible so an error will raise if no valid ending population
+ * size or forward migration matrix could be determined. Please refer to
+ * the simuPOP user's guide for an explanation of the theory behind forward
+ * and backward migration matrices.
+ */
+class BackwardMigrator : public BaseOperator
+{
+public:
+ /** Create a BackwardMigrator that moves individuals between \e subPop
+ * subpopulations randomly according to a backward migration matrix \e rate.
+ * The size of the matrix should match the number of subpopulations.
+ *
+ * Depending on the value of parameter \e mode, elements in the migration
+ * matrix (\e rate) are interpreted as either the probabilities to originate
+ * from source subpopulations (\e mode = \c BY_PROBABILITY) or proportions of
+ * individuals originate from the source (virtual) subpopulations (\e mode
+ * = \c BY_PROPORTION). Migration by count is not supported by this operator.
+ *
+ * Please refer to operator \c BaseOperator for a detailed explanation for
+ * all parameters.
+ */
+ BackwardMigrator(const floatMatrix & rate = floatMatrix(), int mode = BY_PROBABILITY,
+ int begin = 0, int end = -1, int step = 1, const intList & at = vectori(),
+ const intList & reps = intList(), const subPopList & subPops = subPopList(),
+ const stringList & infoFields = vectorstr(1, "migrate_to"));
+
+ /// destructor
+ virtual ~BackwardMigrator()
+ {
+ };
+
+ /// HIDDEN Deep copy of a Migrator
+ virtual BaseOperator * clone() const
+ {
+ return new BackwardMigrator(*this);
+ }
+
+
+ /// HIDDEN apply the Migrator to populaiton \e pop.
+ virtual bool apply(Population & pop) const;
+
+ /// HIDDEN
+ string describe(bool format = true) const;
+
+protected:
+ /// migration rate. its meaning is controled by m_mode
+ const matrixf m_rate;
+
+ boost::numeric::ublas::matrix<double> m_inverse_rate;
+
+ bool m_symmetric_matrix;
+
+ /// asProbability (1), asProportion (2),
+ const int m_mode;
+};
+
+
+
+
/** Split a given list of subpopulations according to either sizes of the
* resulting subpopulations, proportion of individuals, or an information
* field. The resulting subpopulations will have the same name as the
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2015-11-01 03:09:28
|
Revision: 4996
http://sourceforge.net/p/simupop/code/4996
Author: bpeng2000
Date: 2015-11-01 03:09:25 +0000 (Sun, 01 Nov 2015)
Log Message:
-----------
Remove rpy/rpy2 support from simuPOP
Modified Paths:
--------------
trunk/doc/refManual.lyx
trunk/doc/userGuide.py
trunk/simuOpt.py
trunk/src/plotter.py
trunk/src/sampling.py
Modified: trunk/doc/refManual.lyx
===================================================================
--- trunk/doc/refManual.lyx 2015-10-31 00:32:59 UTC (rev 4995)
+++ trunk/doc/refManual.lyx 2015-11-01 03:09:25 UTC (rev 4996)
@@ -1,5 +1,5 @@
-#LyX 2.0 created this file. For more info see http://www.lyx.org/
-\lyxformat 413
+#LyX 2.1 created this file. For more info see http://www.lyx.org/
+\lyxformat 474
\begin_document
\begin_header
\textclass manual
@@ -50,13 +50,13 @@
\font_roman default
\font_sans default
\font_typewriter beramono
+\font_math auto
\font_default_family default
\use_non_tex_fonts false
\font_sc false
\font_osf false
\font_sf_scale 100
\font_tt_scale 70
-
\graphics default
\default_output_format default
\output_sync 0
@@ -81,15 +81,24 @@
\pdf_quoted_options "linkcolor=TitleColor,urlcolor=LinkColor"
\papersize default
\use_geometry false
-\use_amsmath 1
-\use_esint 0
-\use_mhchem 1
-\use_mathdots 1
-\cite_engine natbib_authoryear
+\use_package amsmath 1
+\use_package amssymb 1
+\use_package cancel 1
+\use_package esint 0
+\use_package mathdots 1
+\use_package mathtools 1
+\use_package mhchem 1
+\use_package stackrel 1
+\use_package stmaryrd 1
+\use_package undertilde 1
+\cite_engine natbib
+\cite_engine_type authoryear
+\biblio_style plainnat
\use_bibtopic false
\use_indices false
\paperorientation portrait
\suppress_date false
+\justification true
\use_refstyle 0
\index Index
\shortcut idx
@@ -2176,20 +2185,6 @@
\backslash
-simuPOPplotternewDeviceRef
-\end_layout
-
-\begin_layout Plain Layout
-
-
-\backslash
-simuPOPplottersaveFigureRef
-\end_layout
-
-\begin_layout Plain Layout
-
-
-\backslash
simuPOPplotterVarPlotterRef
\end_layout
@@ -2200,34 +2195,6 @@
simuPOPplotterScatterPlotterRef
\end_layout
-\begin_layout Plain Layout
-
-
-\backslash
-simuPOPplotterInfoPlotterRef
-\end_layout
-
-\begin_layout Plain Layout
-
-
-\backslash
-simuPOPplotterHistPlotterRef
-\end_layout
-
-\begin_layout Plain Layout
-
-
-\backslash
-simuPOPplotterQQPlotterRef
-\end_layout
-
-\begin_layout Plain Layout
-
-
-\backslash
-simuPOPplotterBoxPlotterRef
-\end_layout
-
\end_inset
@@ -2390,13 +2357,6 @@
simuPOPsamplingdrawCombinedSamplesRef
\end_layout
-\begin_layout Plain Layout
-
-
-\backslash
-simuPOPsamplingplotPedigreeRef
-\end_layout
-
\end_inset
Modified: trunk/doc/userGuide.py
===================================================================
--- trunk/doc/userGuide.py 2015-10-31 00:32:59 UTC (rev 4995)
+++ trunk/doc/userGuide.py 2015-11-01 03:09:25 UTC (rev 4996)
@@ -4783,7 +4783,7 @@
#begin_file log/forwardTrajectory.py
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='matplotlib')
+simuOpt.setOptions(quiet=True)
import simuPOP as sim
#begin_ignore
sim.setRNG(seed=12345)
@@ -4823,7 +4823,7 @@
#begin_file log/backTrajectory.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='matplotlib')
+simuOpt.setOptions(quiet=True)
#end_ignore
import simuPOP as sim
#begin_ignore
@@ -4977,7 +4977,7 @@
#begin_file log/demoModel.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='matplotlib')
+simuOpt.setOptions(quiet=True)
#end_ignore
import simuPOP as sim
from simuPOP.demography import *
@@ -5044,7 +5044,7 @@
#begin_file log/demoEventModel.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='matplotlib')
+simuOpt.setOptions(quiet=True)
#end_ignore
import simuPOP as sim
from simuPOP.demography import *
@@ -5133,7 +5133,7 @@
#begin_file log/varPlotByRep.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='matplotlib')
+simuOpt.setOptions(quiet=True)
#end_ignore
import simuPOP as sim
#begin_ignore
@@ -5149,15 +5149,6 @@
matingScheme=sim.RandomMating(),
postOps=[
sim.Stat(alleleFreq=range(4)),
- # rpy syntax
- #VarPlotter('[alleleFreq[x][0] for x in range(4)]', byRep=True,
- # update=10, saveAs='log/rpy_byRep.png',
- # legend=['Locus %d' % x for x in range(4)],
- # ylab='Allele frequency',
- # ylim=[0, 1],
- # main_rep=['Genetic drift, replicate %d' % x for x in range(3)],
- #),
- # matplotlib syntax
VarPlotter('[alleleFreq[x][0] for x in range(4)]', byRep=True,
update=10, saveAs='log/varplot_byRep.png',
figure_figsize=(10, 8),
@@ -5174,7 +5165,7 @@
#begin_file log/varPlotByDim.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='matplotlib')
+simuOpt.setOptions(quiet=True)
#end_ignore
import simuPOP as sim
#begin_ignore
@@ -5227,7 +5218,7 @@
#begin_file log/ScatterPlotter.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='matplotlib')
+simuOpt.setOptions(quiet=True)
#end_ignore
import simuPOP as sim
#begin_ignore
@@ -5273,62 +5264,7 @@
)
#end_file
-#begin_file log/HistPlotter.py
-#begin_ignore
-import simuOpt
-simuOpt.setOptions(quiet=True, plotter='rpy2')
-#end_ignore
-import simuPOP as sim
-#begin_ignore
-sim.setRNG(seed=12345)
-#end_ignore
-import simuPOP as sim
-from simuPOP.plotter import HistPlotter, QQPlotter, BoxPlotter
-import random
-pop = sim.Population([500], infoFields=['x', 'y', 'anc'])
-# Defines VSP 0, 1, 2, 3, 4 by anc.
-pop.setVirtualSplitter(sim.SexSplitter())
-def passInfo(x, y, anc):
- 'Parental fields will be passed as tuples'
- off_anc = (anc[0] + anc[1])/2.
- off_x = (x[0] + x[1])/2 + random.normalvariate(off_anc - 0.5, 0.1)
- off_y = (y[0] + y[1])/2 + random.normalvariate(0, 0.1)
- return off_x, off_y, off_anc
-
-pop.evolve(
- initOps=[
- sim.InitSex(),
- # random geographic location
- sim.InitInfo(random.random, infoFields=['x', 'y']),
- # anc is 0 or 1
- sim.InitInfo(lambda : random.randint(0, 1), infoFields='anc')
- ],
- matingScheme=sim.RandomMating(ops=[
- sim.MendelianGenoTransmitter(),
- sim.PyTagger(passInfo)]),
- postOps=[
- HistPlotter(infoFields='anc',
- subPops=[(0,0), (0,1)], col_sp=['blue', 'red'],
- saveAs='log/HistPlotter.png',
- main="!'Histogram of ancestry values at generation %d' % gen",
- ),
- QQPlotter(infoFields='anc',
- subPops=[(0,0), (0,1)], col_sp=['blue', 'red'],
- saveAs='log/QQPlotter.png',
- main="!'QQ plot of ancestry values at generation %d' % gen",
- ),
- BoxPlotter(infoFields='anc',
- subPops=[(0,0), (0,1)],
- saveAs='whatever',
- dev_print_file='!"log/Gen%d.png" % gen',
- main="!'Boxplots of ancestry values at generation %d' % gen",
- ),
- ],
- gen = 5,
-)
-#end_file
-
#begin_file log/getParam.py
#begin_ignore
import simuOpt
@@ -5476,16 +5412,17 @@
#end_file
-#begin_file log/plotPedigree.py
+
+#begin_file log/sampleAffectedSibpair.py
#begin_ignore
import simuOpt
simuOpt.setOptions(quiet=True)
#end_ignore
import simuPOP as sim
#begin_ignore
-sim.setRNG(seed=12311)
+sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import indexToID #, plotPedigree
+from simuPOP.sampling import indexToID
pop = sim.Population(size=15, loci=5, infoFields=['father_idx', 'mother_idx'], ancGen=2)
pop.evolve(
preOps=[
@@ -5502,25 +5439,10 @@
print(pop.infoFields())
# save this population for future use
pop.save('log/pedigree.pop')
-# draw pedigree, this function is now deprecated because of the removal of
-# rpy/rpy2 support
-#plotPedigree(pop, filename='log/pedigree.png')
-#end_file
-
-#begin_file log/sampleAffectedSibpair.py
-#begin_ignore
-import simuOpt
-simuOpt.setOptions(quiet=True)
-#end_ignore
-import simuPOP as sim
-#begin_ignore
-sim.setRNG(seed=12347)
-#end_ignore
-from simuPOP.sampling import drawAffectedSibpairSample #, plotPedigree
+from simuPOP.sampling import drawAffectedSibpairSample
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawAffectedSibpairSample(pop, families=2)
-#plotPedigree(sample, filename='log/affectedSibpair.png')
#end_file
@@ -5533,7 +5455,7 @@
#begin_ignore
sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import drawNuclearFamilySample #, plotPedigree
+from simuPOP.sampling import drawNuclearFamilySample
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawNuclearFamilySample(pop, families=2, numOffspring=(2,4),
affectedParents=(1,2), affectedOffspring=(1, 3))
@@ -5547,7 +5469,6 @@
ped1 = sample.extractIndividuals(IDs=0, idField='ped_id')
# print the ID of all individuals in the first pedigree
print([ind.ind_id for ind in ped1.allIndividuals()])
-#plotPedigree(sample, filename='log/nuclerFamily.png')
#end_file
@@ -5561,11 +5482,10 @@
#begin_ignore
sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import drawThreeGenFamilySample #, plotPedigree
+from simuPOP.sampling import drawThreeGenFamilySample
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawThreeGenFamilySample(pop, families=2, numOffspring=(1, 3),
pedSize=(8, 15), numOfAffected=(2, 5))
-#plotPedigree(sample, filename='log/threeGenFamily.png')
#end_file
@@ -5578,13 +5498,12 @@
#begin_ignore
sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import drawCombinedSample, AffectedSibpairSampler, NuclearFamilySampler #, plotPedigree
+from simuPOP.sampling import drawCombinedSample, AffectedSibpairSampler, NuclearFamilySampler
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawCombinedSample(pop, samplers = [
AffectedSibpairSampler(families=1),
NuclearFamilySampler(families=1, numOffspring=(2,4), affectedParents=(1,2), affectedOffspring=(1,3))
])
-#plotPedigree(sample, filename='log/combinedSampling.png')
#end_file
Modified: trunk/simuOpt.py
===================================================================
--- trunk/simuOpt.py 2015-10-31 00:32:59 UTC (rev 4995)
+++ trunk/simuOpt.py 2015-11-01 03:09:25 UTC (rev 4996)
@@ -195,10 +195,7 @@
commandline option ``--gui``.
plotter
- The plotting library used to draw figures. simuPOP will try to
- use ``rpy``, ``rpy2``, and ``matplotlib`` if the module is
- available. You can set this parameter to one of them to force the
- use a particular library.
+ (Deprecated)
quiet
If set to ``True``, suppress the banner message when a simuPOP module
@@ -296,10 +293,8 @@
simuOptions['NumThreads'] = numThreads
elif numThreads is not None:
raise TypeError('An integer number is expected for parameter numThreads.')
- if plotter in ['rpy', 'rpy2', 'matplotlib']:
- simuOptions['Plotter'] = plotter
- elif plotter is not None:
- raise TypeError('Only rpy and matplotlib are allowed for parameter plotter.')
+ if plotter is not None:
+ sys.stderr.write('WARNING: plotter option is deprecated because of the removal of rpy/rpy2 support\n')
#
# define some validataion functions
Modified: trunk/src/plotter.py
===================================================================
--- trunk/src/plotter.py 2015-10-31 00:32:59 UTC (rev 4995)
+++ trunk/src/plotter.py 2015-11-01 03:09:25 UTC (rev 4996)
@@ -26,38 +26,20 @@
#
'''
-This module defines several utility functions and Python operators that make
-use of the Python rpy or rpy2 module (http://rpy.sourceforge.net) to plot
-expressions and information fields of evolving populations using a popular
-statistical analysis language R (http://www.r-project.org) or matplotlib.
-You can specify the plotting library using simuOpt.setOptions(plotter='matplotlib').
+This module defines several utility functions and Python operators that plot
+expressions and information fields of evolving populations using matplotlib.
+Plotting through rpy and rpy2 was supported in version 1.1.6 and earlier but
+was removed in version 1.1.7 due to stability problems of rpy2.
-Each operator calls a sequence of R or matplotlib functions to draw and save
+Each operator calls a sequence of matplotlib functions to draw and save
figures. A special parameter passing mechanism is used so that you can specify
-arbitrary parameters to these functions. For example, you can use parameter
-``par_mfrow=[2,2]`` to pass ``mfrow=[2,2]`` to function ``par``, and use
-``lty_rep=[1,2]`` to pass ``lty=1`` and ``lty=2`` to specify different line
-types for different replicates if ``rpy`` is used. The help message of each class
-will describe which and in what sequence these R or matplotlib functions are
-called to help you figure out which parameters are allowed. The syntax
-will be different if different plotting module is used.
+arbitrary parameters to these functions.
'''
__all__ = [
- 'newDevice',
- 'saveFigure',
'DerivedArgs',
'VarPlotter',
'ScatterPlotter',
- 'InfoPlotter',
- 'HistPlotter',
- 'QQPlotter',
- 'BoxPlotter',
- # export essential piece of rpy or rpy2 so that other modules could use them
- 'r',
- # exported with rpy. Not meaningful if rpy2 is used.
- 'with_mode',
- 'NO_CONVERSION'
]
from math import ceil, sqrt
@@ -65,177 +47,14 @@
from simuOpt import simuOptions
-# used by rpy2
-def my_py2ri(o):
- '''Convert tuple to R object. This is defined because the default
- py2ri function does not support tuple type.'''
- if isinstance(o, (tuple, list)):
- if all([isinstance(x, int) for x in o]):
- res = ro.vectors.IntVector(o)
- elif all([isinstance(x, (int, float)) for x in o]):
- res = ro.vectors.FloatVector(o)
- elif all([isinstance(x, str) for x in o]):
- res = ro.vectors.StrVector(o)
- elif isinstance(o, tuple):
- raise RuntimeError('Failed to convert tuple {}'.format(o))
- else:
- # default_py2ri can handle list of
- # heterogeneous type, I guess
- res = ro.default_py2ri(o)
- else:
- res = ro.default_py2ri(o)
- return res
+import matplotlib.pylab as plt
-if simuOptions['Plotter'] is None:
- try:
- import rpy_options
- rpy_options.set_options(VERBOSE = False)
- from rpy import r, with_mode, NO_CONVERSION
- plotter = 'rpy'
- except ImportError, e:
- try:
- import rpy2.robjects as ro
- from rpy2.robjects import r, conversion
- conversion.py2ri = my_py2ri
- with_mode = None
- NO_CONVERSION = None
- plotter = 'rpy2'
- except ImportError, e:
- try:
- import matplotlib.pylab as plt
- plotter = 'matplotlib'
- except:
- print('Neither rpy nor matplotlib is available.')
- raise e
-elif simuOptions['Plotter'] == 'rpy':
- import rpy_options
- rpy_options.set_options(VERBOSE = False)
- from rpy import r, with_mode, NO_CONVERSION
- plotter = 'rpy'
-elif simuOptions['Plotter'] == 'rpy2':
- import rpy2.robjects as ro
- from rpy2.robjects import r, conversion
- conversion.py2ri = my_py2ri
- with_mode = None
- NO_CONVERSION = None
- plotter = 'rpy2'
-else:
- import matplotlib.pylab as plt
- plotter = 'matplotlib'
-
-
-# if under windows, fix a bug with rpy which uses blocking i/o so
-# R figure will not be refreshed in time. See
-# https://Stat.ethz.ch/pipermail/r-devel/2006-January/036049.html
-# for details.
-if os.name == 'nt':
- # In addition to options(windowsBuffered=False), I find that I also need to
- # call windows.options(buffered=False) to make functions such as hist work.
- #
- # This function is only available for R 2.9.0 (rev 48333)
- if plotter == 'rpy':
- r.options(windowsBuffered=False)
- if int(r.R_Version()['svn rev']) >= 48333:
- r.windows_options(buffered=False)
- elif plotter == 'rpy2':
- if int(r('R.Version()$"svn rev"')[0]) >= 48333:
- r("windows.options(buffered=False)")
-
-
from simuPOP import PyOperator, ALL_AVAIL
-def newDevice():
- '''Create a new graphics window and return its device number in R. This
- function essentially calls ``getOption('device')()`` in R.
- '''
- if plotter == 'rpy':
- # open a new window
- try:
- r('getOption("device")()')
- except Exception as e:
- raise RuntimeError("Failed to get R version to start a graphical device: {}".format(e));
- # get device number
- device = r.dev_cur()
- if device == 0:
- raise RuntimeError('Can not open new device')
- return device
- elif plotter == 'rpy2':
- # open a new window
- try:
- # 46754 is the revision number for R 2.8.0
- # For R >= 2.8.0, getOption('device') returns a function
- r('getOption("device")()')
- except Exception as e:
- raise RuntimeError("Failed to get R version to start a graphical device: {}".format(e));
- # get device number
- device = r['dev.cur']()
- if device == 0:
- raise RuntimeError('Can not open new device')
- return device
-
- else:
- return plt.figure()
-
-def saveFigure(file=None, **kwargs):
- '''
- Save current figure into ``file``. File format and graphics device are
- determined by file extension. Supported file formats include ``pdf``,
- ``png``, ``bmp``, ``jpg (jpeg)``, ``tif (tiff)``, and ``eps``, which
- correspond to R devices ``pdf``, ``png``, ``bmp``, ``jpeg``, ``tiff``
- and ``postscript``. A postscript device will be used if there is no file
- extension or the file extension is not recognizable. Additional keyword
- parameters will be passed to the underlying ``dev.print`` function.
- '''
- if file is None:
- return
- if 'rpy' in plotter:
- filename, ext = os.path.splitext(file)
- dirname = os.path.dirname(file)
- if dirname != '' and not os.path.isdir(dirname):
- # this might fail and raise an error
- os.makedirs(dirname)
- # default extension and format
- if ext == '':
- ext = '.eps'
- #
- params = {}
- # these two parameters have to be specified for raster formats
- try:
- # I need to use this more lengthy form because some
- # functions are not available in, for example, R 2.6.2
- if ext.lower() == '.pdf':
- device = r.pdf
- elif ext.lower() == '.png':
- device = r.png
- params = {'width': 800, 'height': 600}
- elif ext.lower() == '.bmp':
- device = r.bmp
- params = {'width': 800, 'height': 600}
- elif ext.lower() in ['.jpg', '.jpeg']:
- device = r.jpeg
- params = {'width': 800, 'height': 600}
- elif ext.lower() in ['.tif', '.tiff']:
- device = r.tiff
- params = {'width': 800, 'height': 600}
- elif ext.lower() == '.eps':
- device = r.postscript
- except Exception, e:
- print e
- print 'Can not determine which device to use to save file %s. A postscript driver is used.' % name
- device = r.postscript
- params.update(kwargs)
- if plotter == 'rpy':
- r.dev_print(file=file, device=device, **params)
- else:
- r['dev.print'](file=file, device=device, **params)
- else:
- plt.savefig(file)
-
-
class DerivedArgs:
'''This class implements the derived keyword argument handling mechanism that
is used by all classes defined in this module. It is provided for users who
- would like to use this mechanism for their own rpy or matplotlib-related
+ would like to use this mechanism for their own matplotlib-related
operators.
An derived keyword argument is an argument that is prefixed with a function
@@ -417,30 +236,6 @@
various formats if parameter ``saveAs`` is specified. File format is
determined by file extension. After the evolution, the graphic device could
be left open (``leaveOpen``).
-
- If the rpy module is used (default if rpy if available), arbitrary keyword
- parameters could be specified and be passed to the underlying R drawing functions
- ``plot`` and ``lines``. These parameters could be used to specify line type (``lty``),
- color (``col``), title (``main``), limit of x and y axes (``xlim`` and
- ``ylim``) and many other options (see R manual for details). As a special
- case, multiple values can be passed to each replicate and/or dimension if
- the name of a parameter ends with ``_rep``, ``_dim``, or ``_repdim``
- For example, ``lty_rep=range(1, 5)`` will pass parameters ``lty=1``, ...
- ``lty=4`` to four replicates. You can also pass parameters to specific
- R functions such as ``par``, ``plot``, ``lines``, ``legend``, ``dev_print``
- by prefixing parameter names with a function name. For example,
- ``dev_print_width=300`` will pass ``width=300`` to function ``dev.print()``
- when you save your figures using this function. In addition, if the value
- of a parameter is a string starting with ``!``, the evaluated result of
- the remaining string will be used as parameter value. Further customization
- of your figures could be achieved by writing your own hook functions that
- will be called before and after a figure is drawn, and after each ``plot``
- call.
-
- This opertor calls R functions ``par``, ``plot``, ``lines``, ``legend``,
- and ``dev.print``. Functions ``plot`` and ``lines`` are the default
- destination for keyword arguments and the ones that accept list parameters
- to customize lines by replicate and/or dimension.
'''
def __init__(self, expr, win=0, update=1, byRep=False, byDim=False,
saveAs="", leaveOpen=False, legend=[], preHook=None, postHook=None,
@@ -505,18 +300,15 @@
plotHook
A function that, if given, will be called after each ``plot``
- function. The ``r`` object from the ``rpy`` module or ``Figure``
- object from the ``matplotlib`` module , generation
- list, data being plotted, replicate number (if applicable) and
- dimension index (if applicable) will be passed as keyword arguments
- ``r``, ``gen``, ``data``, ``rep`` (optional) and ``dim``
- (optional).
+ function. The ``Figure`` object from the ``matplotlib`` module ,
+ generation list, data being plotted, replicate number (if applicable)
+ and dimension index (if applicable) will be passed as keyword arguments
+ ``gen``, ``data``, ``rep`` (optional) and ``dim`` (optional).
kwargs
Additional keyword arguments that will be interpreted and sent to
- underlying R or matplotlib functions. These arguments could have prefixes
- (destination function names) ``plot_``, ``lines_``, ``par_``,
- ``legend_`` and ``dev_print_``, and suffixes (list parameters)
+ underlying matplotlib functions. These arguments could have prefixes
+ (destination function names) and suffixes (list parameters)
``_rep``, ``_dim``, and ``_repdim`` for the ``rpy`` option. Arguments
without prefixes are sent to functions ``plot`` and ``lines``. String
values with a leading ``!`` will be replaced by its evaluated result
@@ -534,35 +326,19 @@
self.preHook = preHook
self.postHook = postHook
self.plotHook = plotHook
- if 'rpy' in plotter:
- self.args = DerivedArgs(
- defaultFuncs = ['plot', 'lines'],
- allFuncs = ['par', 'plot', 'lines', 'dev_print', 'legend'],
- suffixes = ['rep', 'dim', 'repdim'],
- defaultParams = {
- 'plot_type': 'l',
- 'plot_xlab': 'Generation',
- 'plot_ylab': '',
- 'lines_lty': 1,
- 'legend_bty': 'n',
- 'legend_x': 'topright',
- },
- **kwargs
- )
- else:
- self.args = DerivedArgs(
- defaultFuncs = ['figure', 'plot'],
- allFuncs = ['figure', 'plot', 'set_title', 'set_xlabel',
- 'set_ylabel', 'set_ylim', 'legend'],
- suffixes = ['rep', 'dim', 'repdim'],
- defaultParams = {
- 'plot_linestyle': '-',
- 'set_xlabel_xlabel': 'Generation',
- 'set_ylabel_ylabel': '',
- 'set_title_label': '',
- },
- **kwargs
- )
+ self.args = DerivedArgs(
+ defaultFuncs = ['figure', 'plot'],
+ allFuncs = ['figure', 'plot', 'set_title', 'set_xlabel',
+ 'set_ylabel', 'set_ylim', 'legend'],
+ suffixes = ['rep', 'dim', 'repdim'],
+ defaultParams = {
+ 'plot_linestyle': '-',
+ 'set_xlabel_xlabel': 'Generation',
+ 'set_ylabel_ylabel': '',
+ 'set_title_label': '',
+ },
+ **kwargs
+ )
# internal flags
self.nRep = 0
self.reps = [] # allows specification of selected replicates
@@ -574,21 +350,13 @@
self.gen = []
self.data = []
# when apply is called, self._rpy_plot is called.
- PyOperator.__init__(self, func=self._rpy_plot if 'rpy' in plotter else self._mat_plot,
+ PyOperator.__init__(self, func=self._mat_plot,
begin=begin, end=end, step=step, at=at, reps=reps,
subPops=ALL_AVAIL, infoFields=[])
def __del__(self):
- # Close the device if needed.
- if plotter == 'rpy':
- if not self.leaveOpen and hasattr(self, 'device'):
- r.dev_off(self.device)
- elif plotter == 'rpy2':
- if not self.leaveOpen and hasattr(self, 'device'):
- r['dev.off'](self.device)
- else:
- if not self.leaveOpen:
- plt.close()
+ if not self.leaveOpen:
+ plt.close()
def _pushData(self, gen, rep, data):
'''Push history data to self.data for later retrieval. Set self.min and
@@ -639,133 +407,6 @@
else:
return self.data[rep]
- def _rpy_plot(self, pop):
- "Evaluate expression in pop and save result. Plot all data if needed"
- gen = pop.dvars().gen
- rep = pop.dvars().rep
- # push data
- self._pushData(gen, rep, pop.evaluate(self.expr))
- # Draw a plot only when
- # 1. There are at least two obervations.
- # 2. rep is the last recorded replicate.
- # 3. we are self.update away from last plot.
- if len(self.gen) <= 1 or rep + 1 != len(self.data) or \
- (self.update >= 1 and gen < self.lastPlot + self.update):
- # do not plot
- return True
- else:
- self.lastPlot = gen
- # create a new graphical device if needed
- if not hasattr(self, 'device'):
- self.device = newDevice()
- else: # if there are multiple devices, set it back
- if plotter == 'rpy':
- r.dev_set(self.device)
- else:
- r['dev.set'](self.device)
- # call the preHook function if given
- if self.preHook is not None:
- self.preHook(r)
- # figure out the dimension of data
- if self.nRep == 1:
- self.byRep = False
- if self.nDim == 1:
- self.byDim = False
- # needs subplots?
- nPlots = 1
- if self.byDim:
- nPlots *= self.nDim
- if self.byRep:
- nPlots *= self.nRep
- # try to use colors
- if self.nDim > 1 and not self.byDim:
- self.args.addDefault(col_dim=r.rainbow(self.nDim))
- if self.nRep > 1 and not self.byRep:
- self.args.addDefault(col_rep=r.rainbow(self.nRep))
- # suggest how to arrange subplots
- if nPlots > 1:
- nrow = int(ceil(sqrt(nPlots)))
- ncol = int(ceil(nPlots/float(nrow)))
- if nrow > ncol:
- nrow, ncol = ncol, nrow
- self.args.addDefault(par_mfrow=[nrow, ncol])
- # users might set additional parameters and override calculated mfrow.
- r.par(**self.args.getArgs('par', pop))
- # now plot.
- if self.byRep:
- # handle each replicate separately
- for rep_idx,rep in enumerate(self.reps):
- if self.byDim:
- # separate plot for each dim
- for dim in range(self.nDim):
- data = self._getData(rep, dim)
- r.plot(self.gen, data,
- **self.args.getArgs('plot', pop, rep=rep_idx, dim=dim,
- repdim=self.nDim*rep_idx + dim, ylim=[self.min, self.max]))
- if self.plotHook is not None:
- self.plotHook(r=r, gen=self.gen, data=data, rep=rep, dim=dim)
- else:
- # all var in one subplot
- data = self._getData(rep, 0)
- r.plot(self.gen, data,
- **self.args.getArgs('plot', pop, rep=rep_idx, dim=0,
- repdim=self.nDim * rep_idx, ylim=[self.min, self.max]))
- if self.plotHook is not None:
- self.plotHook(r=r, gen=self.gen, data=data, rep=rep)
- for dim in range(1, self.nDim):
- r.lines(self.gen, self._getData(rep, dim),
- **self.args.getArgs('lines', pop, rep=rep_idx, dim=dim,
- repdim=self.nDim * rep_idx + dim))
- if len(self.legend) > 0:
- args = self.args.getLegendArgs('lines', pop, ['lty', 'col', 'lwd'],
- 'rep', range(self.nRep))
- args.update(self.args.getArgs('legend', pop))
- r.legend(legend=self.legend, **args)
- else:
- # all replicate in one figure
- if self.byDim:
- for dim in range(self.nDim):
- data = self._getData(self.reps[0], dim)
- r.plot(self.gen, data,
- **self.args.getArgs('plot', pop, rep=self.reps[0], dim=dim, repdim=dim,
- ylim=[self.min, self.max]))
- if self.plotHook is not None:
- self.plotHook(r=r, gen=self.gen, data=data, dim=dim)
- for rep_idx,rep in enumerate(self.reps[1:]):
- r.lines(self.gen, self._getData(rep, dim),
- **self.args.getArgs('lines', pop, rep=rep_idx+1, dim=dim,
- repdim=self.nDim * (rep_idx + 1) + dim))
- if len(self.legend) > 0:
- args = self.args.getLegendArgs('lines', pop, ['lty', 'col', 'lwd'],
- 'rep', range(self.nRep))
- args.update(self.args.getArgs('legend', pop))
- r.legend(legend=self.legend, **args)
- else:
- data = self._getData(0, 0)
- r.plot(self.gen, data,
- **self.args.getArgs('plot', pop, rep=0, dim=0, repdim=0,
- ylim=[self.min, self.max]))
- if self.plotHook is not None:
- self.plotHook(r=r, gen=self.gen, data=data)
- for rep_idx,rep in enumerate(self.reps):
- for dim in range(self.nDim):
- r.lines(self.gen, self._getData(rep, dim),
- **self.args.getArgs('lines', pop, rep=rep_idx, dim=dim,
- repdim=self.nDim * rep_idx + dim))
- if len(self.legend) > 0:
- args = self.args.getLegendArgs('lines', pop, ['lty', 'col', 'lwd'],
- ['rep', 'dim'], [(x,y) for x in range(self.nRep) for y in range(self.nDim)])
- args.update(self.args.getArgs('legend', pop))
- r.legend(legend=self.legend, **args)
- # call the postHook function if given
- if self.postHook is not None:
- self.postHook(r)
- if self.saveAs != '':
- file, ext = os.path.splitext(self.saveAs)
- filename = '%s_%d%s' % (file, gen, ext)
- saveFigure(**self.args.getArgs('dev_print', pop, file=filename))
- return True
-
def _mat_plot(self, pop):
"Evaluate expression in pop and save result. Plot all data if needed"
gen = pop.dvars().gen
@@ -917,18 +558,7 @@
Population, using values at two information fields as their x- and y-axis.
Arbitrary keyword parameters could be specified and be passed to the
- underlying R drawing functions ``plot`` and ``points``. These parameters
- could be used to specify point type (``pch``), color (``col``),
- title (``main``), limit of x and y axes (``xlim`` and ``ylim``) and many
- other options (see R manual for details). You can also pass parameters
- to specific R functions such as ``par``, ``plot``, ``points``, ``legend``,
- ``pdf`` by prefixing parameter names with a function name. For example,
- ``par_mar=[1]*4`` will pass ``par=[1]*4`` to function ``par()`` which is
- called before a figure is drawn. (Note that the function to save a figure
- is ``dev.print`` so parameters such as ``dev_print_width`` should be
- used.) Further customization of your figures could be achieved by writing
- your own hook functions that will be called before and after a figure is
- drawn.
+ underlying matplotlib drawing functions ``plot`` and ``scatter``.
The power of this operator lies in its ability to differentiate individuals
from different (virtual) subpopulations. If you specify IDs of (virtual)
@@ -941,11 +571,6 @@
individuals with red. In addition, if the value of a parameter is a string
starting with ``!``, the evaluated result of the remaining string will be
used as parameter value.
-
- This opertor calls R functions ``par``, ``plot``, ``points``, ``legend``,
- and ``dev.print``. Functions ``plot`` and ``points`` are the default
- destination for keyword arguments and the ones that accept list parameters
- to customize lines by (virtual) subpopulation.
'''
def __init__(self, infoFields=[], saveAs="", leaveOpen=False, legend=[],
preHook=None, postHook=None, begin=0, end=-1, step=1,
@@ -980,22 +605,18 @@
preHook
A function that, if given, will be called before the figure is
- draw. The ``r`` object from the ``plotter`` module or the ``Figure``
+ draw. The ``plotter`` module or the ``Figure``
object from matplotlib will be passed to this function.
postHook
A function that, if given, will be called after the figure is
- drawn. The ``r`` object from the ``plotter`` module or the ``Figure``
+ drawn. The ``plotter`` module or the ``Figure``
object from matplotlib will be passed to this function.
kwargs
Additional keyword arguments that will be interpreted and sent to
- underlying R functions. These arguments could have prefixes
- (destination function names) ``plot_``, ``points_``, ``par_``,
- ``legend_`` and ``dev_print_``, and suffixes (list parameters)
- ``_sp``. Arguments without prefixes are sent to functions
- ``plot`` and ``points``. String values with a leading ``!`` will be
- replaced by its evaluated result against the current population.
+ underlying matplotlib functions. String values with a leading ``!``
+ will be replaced by its evaluated result against the current population.
'''
# parameters
self.infoFields = infoFields
@@ -1007,106 +628,29 @@
self.preHook = preHook
self.postHook = postHook
self.subPops = subPops
- if 'rpy' in plotter:
- self.args = DerivedArgs(
- defaultFuncs = ['plot', 'points'],
- allFuncs = ['par', 'plot', 'points', 'dev_print', 'legend'],
- suffixes = ['sp'],
- defaultParams = {
- 'legend_bty': 'n',
- 'plot_xlab': self.infoFields[0],
- 'plot_ylab': self.infoFields[1],
- 'legend_x': 'topright',
- },
- **kwargs)
- else:
- self.args = DerivedArgs(
- defaultFuncs = ['figure', 'plot'],
- allFuncs = ['figure', 'plot', 'scatter', 'set_title', 'set_xlabel',
- 'set_ylabel', 'set_ylim', 'legend'],
- suffixes = ['sp'],
- defaultParams = {
- 'plot_linestyle': '-',
- 'set_xlabel_xlabel': self.infoFields[0],
- 'set_ylabel_ylabel': self.infoFields[1],
- 'set_title_label': '',
- },
- **kwargs
- )
+ self.args = DerivedArgs(
+ defaultFuncs = ['figure', 'plot'],
+ allFuncs = ['figure', 'plot', 'scatter', 'set_title', 'set_xlabel',
+ 'set_ylabel', 'set_ylim', 'legend'],
+ suffixes = ['sp'],
+ defaultParams = {
+ 'plot_linestyle': '-',
+ 'set_xlabel_xlabel': self.infoFields[0],
+ 'set_ylabel_ylabel': self.infoFields[1],
+ 'set_title_label': '',
+ },
+ **kwargs
+ )
if len(self.subPops) > 1:
- if 'rpy' in plotter:
- self.args.addDefault(
- pch_sp = range(1, len(self.subPops) + 1),
- col_sp = r.rainbow(len(self.subPops)))
- else:
- cm = plt.get_cmap('gist_rainbow')
- self.args.addDefault(
- scatter_c_sp=[cm(i*1.0/len(self.subPops)) for i in range(len(self.subPops))])
+ cm = plt.get_cmap('gist_rainbow')
+ self.args.addDefault(
+ scatter_c_sp=[cm(i*1.0/len(self.subPops)) for i in range(len(self.subPops))])
# when apply is called, self._rpy_plot is called.
- PyOperator.__init__(self, func=self._rpy_plot if 'rpy' in plotter else self._mat_plot,
+ PyOperator.__init__(self, func=self._mat_plot,
begin=begin, end=end, step=step, at=at, reps=reps)
- def __del__(self):
- # Close the device if needed.
- if not self.leaveOpen and hasattr(self, 'device'):
- if plotter == 'rpy':
- r.dev_off(self.device)
- elif plotter == 'rpy2':
- r['dev.off'](self.device)
-
- def _rpy_plot(self, pop):
- "Evaluate expression in pop and save result. Plot all data if needed"
- gen = pop.dvars().gen
- rep = pop.dvars().rep
- # create a new graphical device if needed
- if not hasattr(self, 'device'):
- self.device = newDevice()
- else: # if there are multiple devices, set it back
- if plotter == 'rpy':
- r.dev_set(self.device)
- else:
- r['dev.set'](self.device)
- # call the preHook function if given
- if self.preHook is not None:
- self.preHook(r)
- # call par in case some parameter is provided
- parParam = self.args.getArgs('par', pop)
- if len(parParam) > 0:
- r.par(**parParam)
- #
- x = pop.indInfo(self.infoFields[0])
- y = pop.indInfo(self.infoFields[1])
- xlim = [min(x), max(x)]
- ylim = [min(y), max(y)]
- # if there is no subpopulation, easy
- if len(self.subPops) == 0:
- r.plot(x, y,
- **self.args.getArgs('plot', pop, type='p', xlim=xlim, ylim=ylim))
- else:
- parPlot = self.args.getArgs('plot', pop, type='n', xlim=xlim, ylim=ylim)
- parPlot['type'] = 'n'
- r.plot(x[0], y[0], **parPlot)
- for idx,sp in enumerate(self.subPops):
- x = pop.indInfo(self.infoFields[0], sp)
- y = pop.indInfo(self.infoFields[1], sp)
- r.points(x, y, **self.args.getArgs('points', pop, sp=idx))
- # legend
- if len(self.legend) > 0:
- args = self.args.getLegendArgs('points', pop, ['col', 'pch', 'lwd', 'cex'],
- 'sp', range(len(self.subPops)))
- args.update(self.args.getArgs('legend', pop))
- r.legend(legend=self.legend, **args)
- # call the postHook function if given
- if self.postHook is not None:
- self.postHook(r)
- if self.saveAs != '':
- file, ext = os.path.splitext(self.saveAs)
- filename = '%s_%d_%d%s' % (file, gen, rep, ext)
- saveFigure(**self.args.getArgs('dev_print', pop, file=filename))
- return True
-
def _mat_plot(self, pop):
"Evaluate expression in pop and save result. Plot all data if needed"
gen = pop.dvars().gen
@@ -1157,464 +701,3 @@
self.device.savefig(filename)
return True
-class InfoPlotter(PyOperator):
- '''
- This operator uses a R function such as ``hist`` and ``qqplot`` to plot
- properties of one or more information fields of individuals in one or more
- (virtual) subpopulations. Separate subplots are used for different
- information fields and subpopulations.
-
- This operator essentially gets values of information fields and sends them
- to a R function such as ``hist``. The resulting figures could be customized
- by additional keyword parameters and various hook functions. For example,
- a ``qqline`` function could be called in a ``plotHook`` function to add a
- QQ line to a ``qqnorm`` plot. The ``plotHook`` can be used to draw the
- whole (sub)plot if no R function is specified for parameter ``func``.
-
- Besides regular keyword parameters, keyword parameters ending in ``_sp``,
- ``_fld`` or ``_spfld`` are expected to have multiple values which will be
- used for differnt subpopulations, information fields, and their
- combinations. You can also specify which function the keyword should be
- sent by prefixing a function name to the parameter name. For example,
- ``pch_fld=[1, 2]`` will use different symbols for different information
- fields, and ``par_mar=[1]*4`` will send parameter ``mar=[1]*4`` to function
- ``par``. In addition, if the value of a parameter is a string starting with
- ``!``, the evaluated result of the remaining string will be used as
- parameter value.
-
- This opertor calls R functions ``par``, ``dev.print``, and a user-specified
- function. Additional keyword arguments without function prefix will be sent
- to this function.
- '''
- def __init__(self, func=None, infoFields=[], saveAs="", leaveOpen=False,
- preHook=None, postHook=None, plotHook = None, begin=0,
- end=-1, step=1, at=[], reps=ALL_AVAIL, subPops=[], **kwargs):
- '''
- func
- Name of the R or matplotlib function that will be called to draw figures
- from values of given information fields. No R function will be called
- if it is not specified. In this case, a ``plotHook`` can be used
- to plot passed values.
-
- infoFields
- Information fields whose values will be sent to the specified
- plotting function.
-
- subPops
- A list of subpopulations and virtual subpopulations. Each
- subpopulation will be plotted in a separate subplot.
-
- saveAs
- Save figures in files saveAs_gen_rep.ext (e.g. ``figure_10_0.eps``
- if ``saveAs='figure.eps'``). If ext is given, a corresponding
- device will be used. Otherwise, a default postscript driver will be
- used. Currently supported formats include ``.pdf``, ``.png``,
- ``.bmp``, ``.jpg``, and ``.tif``. The default filename could be
- overriden by derived argument ``dev_print_file``.
-
- leaveOpen
- Whether or not leave the plot open when plotting is done. Default
- to ``False`` functions. If this option is set to ``True``, you will
- have to close the graphic device explicitly using function
- ``r.dev_off()``. Note that leaving the device open allows
- further manipuation of the figures outside of this operator.
-
- preHook
- A function that, if given, will be called before the figure is
- draw. The ``r`` object for ``rpy`` will be passed to this function.
-
- postHook
- A function that, if given, will be called after the figure is
- drawn. The ``r`` object for ``rpy`` will be passed to this function.
-
- plotHook
- A function that, if given, will be called after each ``plot``
- function. The ``r`` object from the ``rpy`` module , data being
- plotted, name of the information field and index of subpopulation
- (in parameter ``subPops``, if applicable) will be passed with
- keywords ``r``, ``data``, ``field`` and ``subPop`` (optional)
- respectively.
-
- kwargs
- Additional keyword arguments that will be interpreted and sent to
- underlying R functions. These arguments could have prefixes
- (destination function names) ``par_``, ``dev_print_`` and the
- function you specify (parameter ``func``), and suffixes (list
- parameters) ``_sp``, ``_fld``, and ``_spfld``. Arguments without
- prefixes are sent to the user specified function. String values
- with a leading ``!`` will be replaced by its evaluated result
- against the current population.
- '''
- # parameters
- if type(infoFields) == type(''):
- self.infoFields = [infoFields]
- else:
- self.infoFields = infoFields
- self.func = func
- if len(self.infoFields) == 0:
- raise RuntimeError('At least one information field should be given')
- self.saveAs = saveAs
- self.leaveOpen = leaveOpen
- self.preHook = preHook
- self.postHook = postHook
- self.plotHook = plotHook
- self.subPops = subPops
- if 'rpy' in plotter:
- if self.func is not None:
- if plotter == 'rpy':
- self.rfunc = r(self.func)
- else:
- self.rfunc = r[self.func]
- self.args = DerivedArgs(
- defaultFuncs = [] if self.func is None else [self.func],
- allFuncs = ['par', 'dev_print', 'legend'] + ([] if self.func is None else [self.func]),
- suffixes = ['sp', 'fld', 'spfld'],
- defaultParams = {'legend_x': 'topright'},
- **kwargs)
- else:
- raise RuntimeError('InfoPlotter does not yet support matplotlib plotter')
- # when apply is called, self._rpy_plot is called.
- PyOperator.__init__(self, func=self._rpy_plot,
- begin=begin, end=end, step=step, at=at, reps=reps)
-
-
- def __del__(self):
- # Close the device if needed.
- if plotter == 'rpy':
- if not self.leaveOpen and hasattr(self, 'device'):
- r.dev_off(self.device)
- elif plotter == 'rpy2':
- if not self.leaveOpen and hasattr(self, 'device'):
- r['dev.off'](self.device)
- else:
- if not self.leaveOpen:
- plt.close()
-
- def _rpy_plot(self, pop):
- "Evaluate expression in pop and save result. Plot all data if needed"
- gen = pop.dvars().gen
- rep = pop.dvars().rep
- # create a new graphical device if needed
- if not hasattr(self, 'device'):
- self.device = newDevice()
- else: # if there are multiple devices, set it back
- if plotter == 'rpy':
- r.dev_set(self.device)
- else:
- r['dev.set'](self.device)
- # call the preHook function if given
- if self.preHook is not None:
- self.preHook(r)
- # subplots?
- nPlots = len(self.infoFields)
- if len(self.subPops) > 1:
- nPlots *= len(self.subPops)
- # call par in case some parameter is provided
- if nPlots > 1:
- nrow = int(ceil(sqrt(nPlots)))
- ncol = int(ceil(nPlots/float(nrow)))
- if nrow > ncol:
- nrow, ncol = ncol, nrow
- self.args.addDefault(par_mfrow=[nrow, ncol])
- #
- r.par(**self.args.getArgs('par', pop))
- #
- for fldIdx,fld in enumerate(self.infoFields):
- # if there is no subpopulation, easy
- if len(self.subPops) == 0:
- val = pop.indInfo(fld)
- if self.func is not None:
- self.rfunc(val, **self.args.getArgs(self.func, pop, fld=fldIdx, sp=0,
- spfld=fldIdx, main='%s at gen %d' % (fld, gen), xlab=fld, ylab=self.func))
- if self.plotHook is not None:
- self.plotHook(r=r, data=val, field=fld)
- else:
- for spIdx,sp in enumerate(self.subPops):
- val = pop.indInfo(fld, sp)
- if self.func is not None:
- self.rfunc(val, **self.args.getArgs(self.func, pop,
- fld=fldIdx, sp=spIdx, spfld=len(self.infoFields)*spIdx + fldIdx,
- main='%s in %s at gen %d' % (fld, pop.subPopName(sp), gen),
- xlab=fld, ylab=self.func))
- if self.plotHook is not None:
- self.plotHook(r=r, data=val, field=fld, subPop=sp)
- # call the postHook function if given
- if self.postHook is not None:
- self.postHook(r)
- if self.saveAs != '':
- file, ext = os.path.splitext(self.saveAs)
- filename = '%s_%d_%d%s' % (file, gen, rep, ext)
- saveFigure(**self.args.getArgs('dev_print', pop, file=filename))
- return True
-
-class HistPlotter(InfoPlotter):
- '''An ``InfoPlotter`` that uses function ``hist``.'''
- def __init__(self, *args, **kwargs):
- '''Returns an ``InfoPlotter`` that uses R function ``hist`` to draw
- histogram of individual information fields of specified (virtual)
- subpopulations. Please see ``InfoPlotter`` for details.
- '''
- InfoPlotter.__init__(self, 'hist', *args, **kwargs)
-
-
-class QQPlotter(InfoPlotter):
- '''An ``InfoPlotter`` that uses function ``qqnorm``.'''
- def __init__(self, *args, **kwargs):
- '''Returns an ``InfoPlotter`` that uses R function ``qqnorm`` to draw
- qq plot of individual information fields of specified (virtual)
- subpopulations. Please see ``InfoPlotter`` for details.
- '''
- InfoPlotter.__init__(self, 'qqnorm', *args, **kwargs)
-
-
-class BoxPlotter(PyOperator):
- '''
- This operator draws boxplots of one or more information fields of
- individuals in one or more (virtual) subpopulations of a population.
- Although a ``InfoPlotter`` with ``func=boxplot`` could be used to plot
- boxplots for each information field and/or subpopulation, this class allows
- multiple whiskers to share one plot. How the whiskers are oraganized is
- controlled by parameters ``byField`` and ``bySubPop``.
-
- This operator essentially gets values of information fields and sends them
- to boxplots. Individual ownerships (subpopulation or field) are also passed
- so that multiple whiskers could be drawn in the same plot. The resulting
- figures could be customized by additional keyword parameters and various
- hook functions.
-
- Besides regular keyword parameters, keyword parameters ending in ``_sp``,
- ``_fld`` or ``_spfld`` are expected to have multiple values which will be
- used for differnt subpopulations, information fields, and their
- combinations. You can also specify which function the keyword should be
- sent by prefixing a function name to the parameter name. For example,
- ``pch_fld=[1, 2]`` will use different symbols for different information
- fields, and ``par_mar=[1]*4`` will send parameter ``mar=[1]*4`` to function
- ``par``. In addition, if the value of a parameter is a string starting with
- ``!``, the evaluated result of the remaining string will be used as
- parameter value.
-
- This opertor calls R functions ``par``, ``boxplot`` and ``dev.print``.
- Keyword parameters without function prefix will be passed to ``boxplot``.
- '''
- def __init__(self, infoFields=[], byField=False, bySubPop=False, saveAs="",
- leaveOpen=False, preHook=None, postHook=None, plotHook = None,
- begin=0, end=-1, step=1, at=[], reps=ALL_AVAIL, subPops=[], **kwargs):
- '''
- infoFields
- Information fields whose values will be sent to R function
- ``boxplot``.
-
- subPops
- A list of subpopulations and virtual subpopulations. Separate
- whiskers will be drawn for individuals in these subpopulations.
-
- byField
- If multiple information fields are specified, separate the whiskers
- different subplots if this parameter is ``True``.
-
- bySubPop
- If multiple (virtual) subpopulations are specified, separate the
- whiskers to different subplots if this parameter is ``True``.
-
- saveAs
- Save figures in files saveAs_gen_rep.ext (e.g. ``figure_10_0.eps``
- if ``saveAs='figure.eps'``). If ext is given, a corresponding
- device will be used. Otherwise, a default postscript driver will be
- used. Currently supported formats include ``.pdf``, ``.png``,
- ``.bmp``, ``.jpg``, and ``.tif``. The default filename could be
- overriden by derived argument ``dev_print_file``.
-
- leaveOpen
- Whether or not leave the plot open when plotting is done. Default
- to ``False`` functions. If this option is set to ``True``, you will
- have to close the graphic device explicitly using function
- ``r.dev_off()``. Note that leaving the device open allows
- further manipuation of the figures outside of this operator.
-
- preHook
- A function that, if given, will be called before the figure is
- draw. The ``r`` object for ``rpy`` will be passed to this function.
-
- postHook
- A function that, if given, will be called after the figure is
- drawn. The ``r`` object for ``rpy`` will be passed to this function.
-
- plotHook
- A function that, if given, will be called after each ``plot``
- function. The ``r`` object from the ``rpy`` module, current field
- and subpopulation will be passed with keywords ``r``, ``field`` and
- ``subPop`` if applicable.
-
- kwargs
- Additional keyword arguments that will be interpreted and sent to
- underlying R functions. These arguments could have prefixes
- (destination function names) ``plot_``, ``boxplot_``, ``par_``,
- and ``dev_print_``, and suffixes (list parameters) ``_sp``,
- ``_fld`` and ``_spfld``. Arguments without prefixes are sent to
- function ``boxplot``. String values with a leading ``!`` will be
- replaced by its evaluated result against the current population.
- '''
- # parameters
- if type(infoFields) == type(''):
- self.infoFields = [infoFields]
- else:
- self.infoFields = infoFields
- if len(self.infoFields) == 0:
- raise RuntimeError('At least one information field should be given')
- if 'rpy' not in plotter:
- raise RuntimeError('BoxPlotter function does not support matplotlib plotter')
- self.saveAs = saveAs
- self.leaveOpen = leaveOpen
- self.preHook = preHook
- self.postHook = postHook
- self.plotHook = plotHook
- self.subPops = subPops
- self.byField = byField
- if len(self.infoFields) == 1:
- self.byField = False
- self.bySubPop = bySubPop
- if len(self.subPops) <= 1:
- self.bySubPop = False
- self.args = DerivedArgs(
- defaultFuncs = ['boxplot'],
- allFuncs = ['par', 'boxplot', 'dev_print', 'legend'],
- suffixes = ['sp', 'fld', 'spfld'],
- defaultParams = {'legend_x': 'topright'},
- **kwargs)
- # when apply is called, self.plot is called, additional keyword
- # parameters are passed by kwargs.
- PyOperator.__init__(self, func=self._rpy_plot, begin=begin, end=end,
- step=step, at=at, reps=reps)
-
- def __del__(self):
- # Close the device if needed.
- if not self.leaveOpen and hasattr(self, 'device'):
- if plotter == 'rpy':
- r.dev_off(self.device)
- elif plotter == 'rpy2':
- r['dev.off'](self.device)
-
-
- def _rpy_plot(self, pop):
- "Evaluate expression in pop and save result. Plot all data if needed"
- gen = pop.dvars().gen
- rep = pop.dvars().rep
- # create a new graphical device if needed
- if not hasattr(self, 'device'):
- self.device = newDevice()
- else: # if there are multiple devices, set it back
- if plotter == 'rpy':
- r.dev_set(self.device)
- else:
- r['dev.set'](self.device)
- # call the preHook function if given
- if self.preHook is not None:
- self.preHook(r)
- # subplots?
- nPlots = 1
- if len(self.infoFields) > 1 and self.byField:
- nPlots *= len(self.infoFields)
- if len(self.subPops) > 1 and self.bySubPop:
- nPlots *= len(self.subPops)
- # call par in case some parameter is provided
- if nPlots > 1:
- nrow = int(ceil(sqrt(nPlots)))
- ncol = int(ceil(nPlots/float(nrow)))
- if nrow > ncol:
- nrow, ncol = ncol, nrow
- self.args.addDefault(par_mfrow=[nrow, ncol])
- #
- r.par(**self.args.getArgs('par', pop))
- #
- if self.byField:
- for fldIdx,fld in enumerate(self.infoFields):
- if self.bySubPop:
- # multiple Field and subpop, each has its own subplot
- for spIdx,sp in enumerate(self.subPops):
- val = pop.indInfo(fld, sp)
- r.boxplot(val, **self.args.getArgs('boxplot', pop,
- fld=fld, sp=sp, spfld=len(self.infoFields)*spIdx + fldIdx,
- main='%s in %s at gen %d' % (fld, pop.subPopName(sp), gen)))
- if self.plotHook is not None:
- self.plotHook(r=r, field=fld, subPop=spIdx)
- else:
- # combine data
- data = []
- owner = []
- if len(self.subPops) == 0:
- data = pop.indInfo(fld)
- owner = [fld]*len(data)
- ...
[truncated message content] |
|
From: <bpe...@us...> - 2015-10-31 00:33:02
|
Revision: 4995
http://sourceforge.net/p/simupop/code/4995
Author: bpeng2000
Date: 2015-10-31 00:32:59 +0000 (Sat, 31 Oct 2015)
Log Message:
-----------
Try to fix rpy/rpy2 related problems
Modified Paths:
--------------
trunk/doc/userGuide.py
trunk/src/plotter.py
trunk/src/sampling.py
Modified: trunk/doc/userGuide.py
===================================================================
--- trunk/doc/userGuide.py 2015-10-30 17:42:45 UTC (rev 4994)
+++ trunk/doc/userGuide.py 2015-10-31 00:32:59 UTC (rev 4995)
@@ -4823,7 +4823,7 @@
#begin_file log/backTrajectory.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='rpy2')
+simuOpt.setOptions(quiet=True, plotter='matplotlib')
#end_ignore
import simuPOP as sim
#begin_ignore
@@ -4843,11 +4843,8 @@
traj = simulateBackwardTrajectory(N=Nt, fitness=fitness, nLoci=2,
endGen=1000, endFreq=[0.1, 0.2])
# matplotlib syntax
-#traj.plot('log/backTrajectory.png', set_ylim_top=0.3, set_ylim_bottom=0,
-# plot_c_loc=['r', 'b'], set_title_label='Simulated Trajectory (backward-time)')
-# rpy syntax
-traj.plot('log/backTrajectory.png', plot_ylim=[0, 0.3], plot_xlim=[0, 1000],
- col_loc=['red', 'blue'], plot_main='Simulated Trajectory (backward-time)')
+traj.plot('log/backTrajectory.png', set_ylim_top=0.3, set_ylim_bottom=0,
+ plot_c_loc=['r', 'b'], set_title_label='Simulated Trajectory (backward-time)')
print('Trajectory simulated with length %s ' % len(traj.traj))
pop = sim.Population(size=Nt(0), loci=[1]*2)
@@ -5177,7 +5174,7 @@
#begin_file log/varPlotByDim.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='rpy2')
+simuOpt.setOptions(quiet=True, plotter='matplotlib')
#end_ignore
import simuPOP as sim
#begin_ignore
@@ -5211,33 +5208,17 @@
matingScheme=sim.RandomMating(),
postOps=[
sim.Stat(alleleFreq=range(4)),
- # rpy and rpy2 syntax
VarPlotter('[alleleFreq[x][0] for x in range(4)]', byDim=True,
- update=10, saveAs='log/rpy_byDim.png',
+ update=10, saveAs='log/varplot_byDim.png',
legend=['Replicate %d' % x for x in range(3)],
- ylab='Allele frequency',
- ylim=[0, 1],
- main_dim=['Genetic drift, freq=%.1f' % ((x+1)*0.10) for x in range(4)],
- col_rep=['red', 'blue', 'black'],
- lty_rep=[1, 2, 3],
- # the default png dimension is 800x600
- dev_print_width=600, dev_print_height=500,
- # do not draw axes in r.plot, leaving the job to drawFrame
- plot_axes=False,
- plotHook = rpy_drawFrame,
+ set_ylabel_ylabel='Allele frequency',
+ set_ylim_bottom=0, set_ylim_top=1,
+ set_title_label_dim=['Genetic drift, freq=%.1f' % ((x+1)*0.10) for x in range(4)],
+ plot_c_rep=['red', 'blue', 'black'],
+ plot_linestyle_rep=['-', '-.', ':'],
+ figure_figsize=(10,8),
+ plotHook = mat_drawFrame,
),
- # matplot lib syntax
- #VarPlotter('[alleleFreq[x][0] for x in range(4)]', byDim=True,
- # update=10, saveAs='log/varplot_byDim.png',
- # legend=['Replicate %d' % x for x in range(3)],
- # set_ylabel_ylabel='Allele frequency',
- # set_ylim_bottom=0, set_ylim_top=1,
- # set_title_label_dim=['Genetic drift, freq=%.1f' % ((x+1)*0.10) for x in range(4)],
- # plot_c_rep=['red', 'blue', 'black'],
- # plot_linestyle_rep=['-', '-.', ':'],
- # figure_figsize=(10,8),
- # plotHook = mat_drawFrame,
- #),
],
gen=100
)
@@ -5246,7 +5227,7 @@
#begin_file log/ScatterPlotter.py
#begin_ignore
import simuOpt
-simuOpt.setOptions(quiet=True, plotter='rpy2')
+simuOpt.setOptions(quiet=True, plotter='matplotlib')
#end_ignore
import simuPOP as sim
#begin_ignore
@@ -5278,26 +5259,14 @@
sim.MendelianGenoTransmitter(),
sim.PyTagger(passInfo)]),
postOps=[
- # rpy syntax
ScatterPlotter(['x', 'y'],
saveAs = 'log/ScatterPlotter.png',
subPops = [(0, 0), (0, 1), (0, 2), (0, 3), (0, 4)],
- ylim = [0, 1.2],
- main = "!'Ancestry distribution of individuals at generation %d' % gen",
+ set_ylim_bottom = 0, set_ylim_top=1.2,
+ set_title_label = "!'Ancestry distribution of individuals at generation %d' % gen",
legend = ['anc < 0.2', '0.2 <= anc < 0.4', '0.4 <= anc < 0.6',
'0.6 <= anc < 0.8', '0.8 <= anc'],
- plot_axes = False,
- par_mar = [0, 0, 2, 0],
),
- # matplotlib syntax
- #ScatterPlotter(['x', 'y'],
- # saveAs = 'log/ScatterPlotter.png',
- # subPops = [(0, 0), (0, 1), (0, 2), (0, 3), (0, 4)],
- # set_ylim_bottom = 0, set_ylim_top=1.2,
- # set_title_label = "!'Ancestry distribution of individuals at generation %d' % gen",
- # legend = ['anc < 0.2', '0.2 <= anc < 0.4', '0.4 <= anc < 0.6',
- # '0.6 <= anc < 0.8', '0.8 <= anc'],
- #),
],
gen = 5,
@@ -5516,7 +5485,7 @@
#begin_ignore
sim.setRNG(seed=12311)
#end_ignore
-from simuPOP.sampling import indexToID, plotPedigree
+from simuPOP.sampling import indexToID #, plotPedigree
pop = sim.Population(size=15, loci=5, infoFields=['father_idx', 'mother_idx'], ancGen=2)
pop.evolve(
preOps=[
@@ -5533,8 +5502,9 @@
print(pop.infoFields())
# save this population for future use
pop.save('log/pedigree.pop')
-# draw pedigree
-plotPedigree(pop, filename='log/pedigree.png')
+# draw pedigree, this function is now deprecated because of the removal of
+# rpy/rpy2 support
+#plotPedigree(pop, filename='log/pedigree.png')
#end_file
@@ -5547,10 +5517,10 @@
#begin_ignore
sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import drawAffectedSibpairSample, plotPedigree
+from simuPOP.sampling import drawAffectedSibpairSample #, plotPedigree
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawAffectedSibpairSample(pop, families=2)
-plotPedigree(sample, filename='log/affectedSibpair.png')
+#plotPedigree(sample, filename='log/affectedSibpair.png')
#end_file
@@ -5563,7 +5533,7 @@
#begin_ignore
sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import drawNuclearFamilySample, plotPedigree
+from simuPOP.sampling import drawNuclearFamilySample #, plotPedigree
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawNuclearFamilySample(pop, families=2, numOffspring=(2,4),
affectedParents=(1,2), affectedOffspring=(1, 3))
@@ -5577,7 +5547,7 @@
ped1 = sample.extractIndividuals(IDs=0, idField='ped_id')
# print the ID of all individuals in the first pedigree
print([ind.ind_id for ind in ped1.allIndividuals()])
-plotPedigree(sample, filename='log/nuclerFamily.png')
+#plotPedigree(sample, filename='log/nuclerFamily.png')
#end_file
@@ -5591,11 +5561,11 @@
#begin_ignore
sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import drawThreeGenFamilySample, plotPedigree
+from simuPOP.sampling import drawThreeGenFamilySample #, plotPedigree
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawThreeGenFamilySample(pop, families=2, numOffspring=(1, 3),
pedSize=(8, 15), numOfAffected=(2, 5))
-plotPedigree(sample, filename='log/threeGenFamily.png')
+#plotPedigree(sample, filename='log/threeGenFamily.png')
#end_file
@@ -5608,13 +5578,13 @@
#begin_ignore
sim.setRNG(seed=12347)
#end_ignore
-from simuPOP.sampling import drawCombinedSample, AffectedSibpairSampler, NuclearFamilySampler, plotPedigree
+from simuPOP.sampling import drawCombinedSample, AffectedSibpairSampler, NuclearFamilySampler #, plotPedigree
pop = sim.loadPopulation('log/pedigree.pop')
sample = drawCombinedSample(pop, samplers = [
AffectedSibpairSampler(families=1),
NuclearFamilySampler(families=1, numOffspring=(2,4), affectedParents=(1,2), affectedOffspring=(1,3))
])
-plotPedigree(sample, filename='log/combinedSampling.png')
+#plotPedigree(sample, filename='log/combinedSampling.png')
#end_file
Modified: trunk/src/plotter.py
===================================================================
--- trunk/src/plotter.py 2015-10-30 17:42:45 UTC (rev 4994)
+++ trunk/src/plotter.py 2015-10-31 00:32:59 UTC (rev 4995)
@@ -151,15 +151,9 @@
if plotter == 'rpy':
# open a new window
try:
- # 46754 is the revision number for R 2.8.0
- if int(r.R_Version()['svn rev']) < 46754:
- # For R < 2.8.0, getOption('device') returns a string (such as 'X11')
- r(r.getOption('device') + '()')
- else:
- # For R >= 2.8.0, getOption('device') returns a function
- r('getOption("device")()')
- except:
- raise RuntimeError("Failed to get R version to start a graphical device");
+ r('getOption("device")()')
+ except Exception as e:
+ raise RuntimeError("Failed to get R version to start a graphical device: {}".format(e));
# get device number
device = r.dev_cur()
if device == 0:
@@ -169,14 +163,10 @@
# open a new window
try:
# 46754 is the revision number for R 2.8.0
- if int(r('R.Version()$"svn rev"')[0]) < 46754:
- # For R < 2.8.0, getOption('device') returns a string (such as 'X11')
- r(r.getOption('device') + '()')
- else:
- # For R >= 2.8.0, getOption('device') returns a function
- r('getOption("device")()')
- except:
- raise RuntimeError("Failed to get R version to start a graphical device");
+ # For R >= 2.8.0, getOption('device') returns a function
+ r('getOption("device")()')
+ except Exception as e:
+ raise RuntimeError("Failed to get R version to start a graphical device: {}".format(e));
# get device number
device = r['dev.cur']()
if device == 0:
Modified: trunk/src/sampling.py
===================================================================
--- trunk/src/sampling.py 2015-10-30 17:42:45 UTC (rev 4994)
+++ trunk/src/sampling.py 2015-10-31 00:32:59 UTC (rev 4995)
@@ -157,14 +157,15 @@
plotter.newDevice()
#
try:
+ from rpy2.robjects.packages import importr
# try to use the kinship2 library
- plotter.r.library('kinship2')
- except:
+ importr('kinship2')
+ except Exception as e:
# if not found, try the older version
try:
plotter.r.library('kinship')
except:
- raise ImportError('Failed to load R library kinship2.')
+ raise ImportError('Failed to load R library kinship2. {}'.format(e))
id = []
dadid = []
momid = []
@@ -207,8 +208,9 @@
momid[idx] = 0
dadid[idx] = 0
# create an object of ped structure recognizable by R library
- ptemp = plotter.with_mode(plotter.NO_CONVERSION, plotter.r.pedigree)(
- id=id, dadid=dadid, momid=momid, sex=sex, affected=aff)
+ #ptemp = plotter.with_mode(plotter.NO_CONVERSION, plotter.r.pedigree)(
+ # id=id, dadid=dadid, momid=momid, sex=sex, affected=aff)
+ ptemp = plotter.r.pedigree(id=id, dadid=dadid, momid=momid, sex=sex, affected=aff)
# plot the ped structure
plotter.r.par(**args.getArgs('par', None))
plotter.r.plot(ptemp, **args.getArgs('plot', None))
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2015-10-30 17:42:47
|
Revision: 4994
http://sourceforge.net/p/simupop/code/4994
Author: bpeng2000
Date: 2015-10-30 17:42:45 +0000 (Fri, 30 Oct 2015)
Log Message:
-----------
Remove __builtins__ before a population is saved because it contains non-pickable items
Modified Paths:
--------------
trunk/src/utility.cpp
Modified: trunk/src/utility.cpp
===================================================================
--- trunk/src/utility.cpp 2015-10-30 04:42:01 UTC (rev 4993)
+++ trunk/src/utility.cpp 2015-10-30 17:42:45 UTC (rev 4994)
@@ -2406,6 +2406,13 @@
if (! pickle)
throw RuntimeError("Failed to import module pickle to serialize population variables.");
+ // Some items in __builtins__ are not pickable so we will have to remove them
+ // before picking.
+ PyObject * builtins = PyString_FromString("__builtins__");
+ if (PyDict_Contains(m_dict, builtins))
+ PyDict_DelItem(m_dict, builtins);
+ Py_DECREF(builtins);
+
// here we use version 2 because this is the latest version that supported by
// both python 2 and python 3, also because it is the one that handles simuPOP's
// defdict type using its __reduce__ interface.
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2015-10-30 04:42:04
|
Revision: 4993
http://sourceforge.net/p/simupop/code/4993
Author: bpeng2000
Date: 2015-10-30 04:42:01 +0000 (Fri, 30 Oct 2015)
Log Message:
-----------
Add an example of dynamically determined offspring population size
Modified Paths:
--------------
trunk/doc/userGuide.lyx
trunk/doc/userGuide.py
Modified: trunk/doc/userGuide.lyx
===================================================================
--- trunk/doc/userGuide.lyx 2015-09-21 15:46:32 UTC (rev 4992)
+++ trunk/doc/userGuide.lyx 2015-10-30 04:42:01 UTC (rev 4993)
@@ -1,5 +1,5 @@
-#LyX 2.0 created this file. For more info see http://www.lyx.org/
-\lyxformat 413
+#LyX 2.1 created this file. For more info see http://www.lyx.org/
+\lyxformat 474
\begin_document
\begin_header
\textclass manual
@@ -52,13 +52,13 @@
\font_roman default
\font_sans default
\font_typewriter beramono
+\font_math auto
\font_default_family default
\use_non_tex_fonts false
\font_sc false
\font_osf false
\font_sf_scale 100
\font_tt_scale 70
-
\graphics default
\default_output_format default
\output_sync 0
@@ -83,15 +83,24 @@
\pdf_quoted_options "linkcolor=TitleColor,urlcolor=LinkColor"
\papersize default
\use_geometry false
-\use_amsmath 1
-\use_esint 0
-\use_mhchem 1
-\use_mathdots 1
-\cite_engine natbib_authoryear
+\use_package amsmath 1
+\use_package amssymb 1
+\use_package cancel 1
+\use_package esint 0
+\use_package mathdots 1
+\use_package mathtools 1
+\use_package mhchem 1
+\use_package stackrel 1
+\use_package stmaryrd 1
+\use_package undertilde 1
+\cite_engine natbib
+\cite_engine_type authoryear
+\biblio_style plainnat
\use_bibtopic false
\use_indices false
\paperorientation portrait
\suppress_date false
+\justification true
\use_refstyle 0
\index Index
\shortcut idx
@@ -473,7 +482,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -1669,7 +1678,7 @@
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -2665,7 +2674,7 @@
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -3101,7 +3110,7 @@
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
Turn on/off debug information
@@ -3964,7 +3973,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -4235,7 +4244,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -18760,7 +18769,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -19401,6 +19410,87 @@
\end_layout
\begin_layout Subsection
+Dynamic population size determined by number of offspring *
+\end_layout
+
+\begin_layout Standard
+What we have described so far requires you to determine the size of offspring
+ population in advance.
+ Each mating event produces a number of offspring that is determined by
+ parameter
+\family typewriter
+NumOffspring
+\family default
+.
+ The mating process stops when the offspring population is filled.
+ This works for most scenarios but there are cases where the offspring populatio
+n size is determined dynamically from a fixed number of mating events with
+ random number of offspring.
+ For example, you might design a mating scheme where all males in a population
+ mate only once and produce random number of offspring.
+\end_layout
+
+\begin_layout Standard
+These kind of mating schemes can be simulated using a demographic model
+ that calculates offspring population size from pre-simulated number of
+ offspring for each family.
+ More specifically, we
+\end_layout
+
+\begin_layout Itemize
+Define a demogrphic function (model) that will be called before mating happens.
+\end_layout
+
+\begin_layout Itemize
+This function determines and save the number of offspring for each mating
+ event, and return the total number of offspring as offspring population
+ size.
+\end_layout
+
+\begin_layout Itemize
+Pass a function or generator to parameter numOffspring to pass pre-determined
+ number of offspring.
+ This function will be called each time when number of offspring is needed.
+\end_layout
+
+\begin_layout Standard
+The number of offspring could be saved and retrieved as global variable
+ but a more clever method is to store the numbers of offspring in a demographic
+ model (class).
+ Example
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "dynamicNumOff"
+
+\end_inset
+
+ demonstrates this method by implementing a demographic model that simulate,
+ save, and return the number of offspring.
+ Note that although we determine the number of mating events from number
+ of males in the parental population, a random mating scheme will choose
+ parents with replacement so it is likely that some parents will be chosen
+ multiple times while some others are not chosen at all.
+ Please refer to section
+\begin_inset Quotes eld
+\end_inset
+
+Non-random and customized mating schemes
+\begin_inset Quotes erd
+\end_inset
+
+ to learn how to define a mating scheme that picks parents without replacement.
+
+\begin_inset CommandInset include
+LatexCommand lstinputlisting
+filename "log/dynamicNumOff.log"
+lstparams "caption={Dynamic population size determined by number of offspring},label={dynamicNumOff}"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsection
Determine sex of offspring
\begin_inset CommandInset label
LatexCommand label
@@ -19956,7 +20046,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -20846,7 +20936,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -23870,7 +23960,7 @@
\end_layout
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -24325,7 +24415,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -24424,7 +24514,7 @@
status collapsed
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -24589,7 +24679,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -25873,7 +25963,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -25944,7 +26034,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -26244,7 +26334,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -26324,7 +26414,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -26395,7 +26485,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -26457,7 +26547,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -27006,7 +27096,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -27164,7 +27254,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -27288,7 +27378,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -27549,7 +27639,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -27734,7 +27824,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -27775,7 +27865,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -27813,7 +27903,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -28211,7 +28301,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -28316,7 +28406,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -28423,7 +28513,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -28542,7 +28632,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -28653,7 +28743,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
@@ -29956,7 +30046,7 @@
status open
\begin_layout Plain Layout
-\begin_inset Caption
+\begin_inset Caption Standard
\begin_layout Plain Layout
\begin_inset CommandInset label
Modified: trunk/doc/userGuide.py
===================================================================
--- trunk/doc/userGuide.py 2015-09-21 15:46:32 UTC (rev 4992)
+++ trunk/doc/userGuide.py 2015-10-30 04:42:01 UTC (rev 4993)
@@ -1389,6 +1389,60 @@
#end_file
+
+#begin_file log/dynamicNumOff.py
+#begin_ignore
+import simuOpt
+simuOpt.setOptions(quiet=True)
+#end_ignore
+import simuPOP as sim
+#begin_ignore
+sim.setRNG(seed=12345)
+#end_ignore
+
+import random
+
+class RandomNumOff:
+ # a demographic model
+ def __init__(self):
+ self.numOff = []
+
+ def getNumOff(self):
+ # return the pre-simulated number of offspring as a generator function
+ for item in self.numOff:
+ yield item
+
+ def __call__(self, pop):
+ # define __call__ so that a RandomNumOff object is callable.
+ #
+ # Each male produce from 1 to 3 offspring. For large population, get the
+ # number of males instead of checking the sex of each individual
+ self.numOff = [random.randint(1, 3) for ind in pop.individuals() if ind.sex() == sim.MALE]
+ # return the total population size
+ print('{} mating events with number of offspring {}'.format(len(self.numOff), self.numOff))
+ return sum(self.numOff)
+
+
+pop = sim.Population(10)
+
+# create a demogranic model
+numOffModel = RandomNumOff()
+
+pop.evolve(
+ preOps=sim.InitSex(),
+ matingScheme=sim.RandomMating(
+ # the model will be called before mating to deteremine
+ # family and population size
+ subPopSize=numOffModel,
+ # the getNumOff function (generator) returns number of offspring
+ # for each mating event
+ numOffspring=numOffModel.getNumOff
+ ),
+ gen=3
+)
+
+#end_file
+
#begin_file log/sexMode.py
#begin_ignore
import simuOpt
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2015-09-21 15:46:35
|
Revision: 4992
http://sourceforge.net/p/simupop/code/4992
Author: bpeng2000
Date: 2015-09-21 15:46:32 +0000 (Mon, 21 Sep 2015)
Log Message:
-----------
Fix test_03_operator for python2/3 compatibility
Modified Paths:
--------------
trunk/test/test_03_operator.py
Modified: trunk/test/test_03_operator.py
===================================================================
--- trunk/test/test_03_operator.py 2015-08-31 15:19:38 UTC (rev 4991)
+++ trunk/test/test_03_operator.py 2015-09-21 15:46:32 UTC (rev 4992)
@@ -266,14 +266,15 @@
with open('test1.txt') as test1txt:
self.assertEqual(test1txt.read(), 'func2'*50)
# runtime error only raise in python 3
- simu = Simulator(Population(), rep=5)
- with open('test1.txt', 'wb') as out:
- # each replicate
- self.assertRaises(RuntimeError, simu.evolve,
- matingScheme=CloneMating(),
- postOps = [
- PyOutput("func2", output=out),
- ], gen=10)
+ if sys.version_info.major > 2:
+ simu = Simulator(Population(), rep=5)
+ with open('test1.txt', 'wb') as out:
+ # each replicate
+ self.assertRaises(RuntimeError, simu.evolve,
+ matingScheme=CloneMating(),
+ postOps = [
+ PyOutput("func2", output=out),
+ ], gen=10)
#
simu = Simulator(Population(), rep=5)
with open('test1.txt', 'wb') as out:
@@ -903,11 +904,14 @@
postOps=PyOperator(func=Functor(), param=self),
gen=10
)
- self.assertEqual(Functor.instance_count, 0)
+ self.assertLessEqual(Functor.instance_count, 1)
def testWrapperOp(self):
- # func
+ # test memory leak of Function or class that involves self
+ # simuPOP automatically clears memory but the first instance
+ # is not cleared (due to the fact that we do not know the
+ # life of the object at the C++ level.
for i in range(100):
pop = Population(10)
pop.evolve(
@@ -915,7 +919,7 @@
postOps=WrapperOpFunc(param=self),
gen=10
)
- self.assertEqual(WrapperOpFunc.instance_count, 0)
+ self.assertLessEqual(WrapperOpFunc.instance_count, 1)
#
for i in range(100):
pop = Population(10)
@@ -924,7 +928,7 @@
postOps=WrapperOpCall(param=self),
gen=10
)
- self.assertEqual(WrapperOpCall.instance_count, 0)
+ self.assertLessEqual(WrapperOpCall.instance_count, 1)
#
for i in range(100):
pop = Population(10)
@@ -933,7 +937,7 @@
postOps=WrapperOpSelf(param=self),
gen=10
)
- self.assertEqual(WrapperOpSelf.instance_count, 0)
+ self.assertLessEqual(WrapperOpSelf.instance_count, 1)
#
class lociFunc:
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <bpe...@us...> - 2015-08-31 15:19:41
|
Revision: 4991
http://sourceforge.net/p/simupop/code/4991
Author: bpeng2000
Date: 2015-08-31 15:19:38 +0000 (Mon, 31 Aug 2015)
Log Message:
-----------
Remove SContruct because it is no longer used
Added Paths:
-----------
trunk/test/sample_64_ba_v3.pop
Removed Paths:
-------------
trunk/SConstruct
Deleted: trunk/SConstruct
===================================================================
--- trunk/SConstruct 2015-08-31 04:46:59 UTC (rev 4990)
+++ trunk/SConstruct 2015-08-31 15:19:38 UTC (rev 4991)
@@ -1,320 +0,0 @@
-#
-# $File: SConstruct $
-# $LastChangedDate: 2009-02-24 15:50:38 -0600 (Tue, 24 Feb 2009) $
-# $Rev: 2491 $
-#
-# This file is part of simuPOP, a forward-time population genetics
-# simulation environment. Please visit http://simupop.sourceforge.net
-# for details.
-#
-# Copyright (C) 2004 - 2009 Bo Peng (bp...@md...)
-#
-# This program is free software: you can redistribute it and/or modify
-# it under the terms of the GNU General Public License as published by
-# the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# This program is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU General Public License for more details.
-#
-# You should have received a copy of the GNU General Public License
-# along with this program. If not, see <http://www.gnu.org/licenses/>.
-
-# A scons based build system for simuPOP. It is better for development
-# usage because it is a real multi-thread, dependency-based build system
-# that can build part of the simuPOP modules..
-#
-# Usage:
-# scons [scons-options] [options] [targets]
-# where:
-# scons-options:
-# standard scons options like -j (number of threads)
-#
-# options:
-# prefix=/path/to/prefix:
-# prefix of installation, equivalent to the --prefix option of python setup.py
-# include-dirs=/path/to/include:/path/to/includes2
-# extra include directories, usually to boost library. The path separator is ; under windows.
-# library-dirs=/path/to/lib:/path/to/lib2
-# extra library directories, usually to boost library. The path separator is ; under windows.
-#
-# targets: one of more of
-# std op la laop ba baop mu muop: individual module
-# all = all modules
-# default to all modules
-# install: install specified targets
-# std_exe op_exe ....: individual executable for debugging purposes.
-# executable: std_exe, op_exe etc
-#
-#
-import os, sys
-import distutils.sysconfig
-
-# load all the module information from setup.py
-from setup import *
-# do not update version for this development version to avoid rebuild
-SIMUPOP_VER, SIMUPOP_REV = simuPOP_version()
-print 'Using swig version %d.%d.%d' % tuple(swig_version())
-# get information from python distutils
-vars = distutils.sysconfig.get_config_vars('CC', 'CXX', 'OPT', 'BASECFLAGS',
- 'CCSHARED', 'LDSHARED', 'SO', 'LIBDEST', 'INCLUDEPY')
-for i in range(len(vars)):
- if vars[i] is None:
- vars[i] = ""
-(cc, cxx, opt, basicflags, ccshared, ldshared, so_ext, lib_dest, python_inc_dir) = vars
-# C++ does not need this option, remove it to avoid annoying warning messages.
-opt = opt.replace('-Wstrict-prototypes', '').replace('-Wall', '')
-
-# get compiler and hack options
-from distutils.ccompiler import new_compiler
-comp = new_compiler()
-if comp.__dict__.has_key('initialized'):
- comp.initialize()
-if not comp.__dict__.has_key('ldflags_shared'):
- comp.ldflags_shared = ''
-if not comp.__dict__.has_key('compile_options'):
- comp.compile_options = []
-
-opts = Variables()
-opts.AddVariables(
- PathVariable('prefix', 'Where to install. see "python setup.py install --prefix"', None),
- PathVariable('include-dirs', 'Extra include directories, see "python setup.py build_ext --help"', None),
- PathVariable('library-dirs', 'Extra library directories, see "python setup.py build_ext --help"', None),
-)
-
-# common environment
-env = Environment(
- options=opts,
- # pass all environment variables because MSVC needs INCLUDE and LIB
- # but they may not exist on other platforms
- ENV=os.environ,
- tools=['default', 'swig'])
-
-# try to use the compiler used to build python
-if cc != '':
- env['CC'] = cc
-if cxx != '':
- env['CXX'] = cxx
-if ldshared != '':
- env['SHLINK'] = ldshared
-
-if env.has_key('prefix') and env['prefix'] is not None:
- pylib_dir = distutils.sysconfig.get_python_lib(plat_specific=1, prefix=env['prefix'])
- dest_dir = os.path.join(pylib_dir, 'simuPOP')
- prefix = env['prefix']
- print "Installing to", dest_dir
-else:
- pylib_dir = distutils.sysconfig.get_python_lib(plat_specific=1, prefix=None)
- dest_dir = os.path.join(pylib_dir, 'simuPOP')
- prefix = sys.prefix
-
-if os.name == 'nt':
- wrapper_flags = []
- extra_path = [os.path.join(sys.exec_prefix, 'libs'), 'win32']
- # treat warning as error. This is only set for development builds (scons)
- # because setup.py builds might be built by imcompatible version of compilers
- extra_flags = ['/Zi', '/WX']
-else:
- wrapper_flags = ['-Wno-unused-parameter', '-Wno-missing-field-initializers']
- extra_path = [os.path.join(os.path.split(lib_dest)[:-1])]
- # during development, output symbols (for profiling and debugging)
- # treat warning as error... (although an unused parameter warning from boost has to be allowed)
- extra_flags = ['-g', '-Wno-unused-parameter', '-Wno-char-subscripts', '-Wno-unused-variable',
- '-Wno-unused-function']
- #
- # GNU gcc:
- #
- # -Wconversion will produce a lot of warnings about type conversion. It is
- # useful but I cannot always turn it on because many warnings from boost will
- # overwhelm the output, and stop compiling if '-Werror' is defined. Because
- # MSVC does a reasonable job in detecting conversion problems and we have /WX
- # turned on for MSVC, we can turn on this option (and disable Werror) from
- # time to time to pick out warnings not detected by MSVC.
- #
- # Intel C++:
- #
- # -Wcheck and -Weffc++ give more warnings about coding style etc. They
- # also should not be turned on by default.
- #
- # -vec-report2 reports loops that are automatically vectorized, and remarks
- # for loops that cannot be vectorized.
- #
- #if USE_ICC:
- # extra_flags.extend(['-Wcheck', '-Weffc++', '-vec-report2'])
- #else:
- # extra_flags.append('-Wconversion')
- #
-
-def convert_def(defines):
- new_list = []
- for d in defines:
- if d[1] is not None:
- new_list.append(d)
- else:
- new_list.append(d[0])
- return new_list
-#
-# Building library gsl
-#
-gsl_env = env.Clone()
-gsl_env.VariantDir('build/gsl', '.')
-gsl_env['SWIGOUTDIR'] = 'build/gsl/src'
-gsl_env['SWIGFLAGS'] = SWIG_CC_FLAGS
-gsl = gsl_env.SharedLibrary(
- target = 'build/gsl/_gsl%s' % so_ext,
- source = ['build/gsl/' + x for x in GSL_FILES] + ['build/gsl/src/gsl.i'],
- SHLIBPREFIX = "",
- SHLIBSUFFIX = so_ext,
- SHLINKFLAGS = comp.ldflags_shared,
- LIBPATH = extra_path,
- # build for config.h
- CPPPATH = ['.', 'gsl', 'gsl/cdf', 'gsl/specfunc', 'build', python_inc_dir],
- CCFLAGS = comp.compile_options,
- CPPFLAGS = ' '.join([basicflags, ccshared, opt])
-)
-Alias('gsl', gsl)
-Alias('all', gsl)
-Alias('install', gsl_env.InstallAs(os.path.join(dest_dir, '_gsl%s' % so_ext), gsl[0]))
-Alias('install', gsl_env.InstallAs(os.path.join(dest_dir, 'gsl.py'), 'build/gsl/src/gsl.py'))
-#
-#
-# Building a library for common files
-#
-# Note that we suppress warnings when building these external source files
-# using -wd options
-std_common_env = env.Clone()
-std_common_env.VariantDir('build/std_common', '.')
-std_common_lib = std_common_env.StaticLibrary(
- target = 'build/std_common/common',
- source = ['build/std_common/' + x for x in LIB_FILES],
- CCFLAGS = ModuInfo('std', SIMUPOP_VER, SIMUPOP_REV)['extra_compile_args'] +\
- comp.compile_options + (['-wd1125', '-wd2259'] if USE_ICC else []),
- CPPPATH = ['.', 'gsl', 'build', ModuInfo('std', SIMUPOP_VER, SIMUPOP_REV)['include_dirs']],
- CPPDEFINES = convert_def(ModuInfo('std', SIMUPOP_VER, SIMUPOP_REV)['define_macros']),
- CPPFLAGS = ccshared + ' ' + opt,
-)
-op_common_env = env.Clone()
-op_common_env.VariantDir('build/op_common', '.')
-op_common_lib = op_common_env.StaticLibrary(
- target = 'build/op_common/common',
- source = ['build/op_common/' + x for x in LIB_FILES],
- CCFLAGS = ModuInfo('op', SIMUPOP_VER, SIMUPOP_REV)['extra_compile_args'] + \
- comp.compile_options + (['-wd1125'] if USE_ICC else []),
- CPPPATH = ['.', 'gsl', 'build', ModuInfo('op', SIMUPOP_VER, SIMUPOP_REV)['include_dirs']],
- CPPDEFINES = convert_def(ModuInfo('op', SIMUPOP_VER, SIMUPOP_REV)['define_macros']),
- CPPFLAGS = ccshared + ' ' + opt,
-)
-Alias('common', (std_common_lib, op_common_lib))
-
-#
-# Building modules
-#
-targets = []
-exe_targets = []
-all_modu = ['std', 'op', 'la', 'laop', 'ba', 'baop', 'mu', 'muop', 'lin', 'linop']
-all_exe = [x + '_exe' for x in all_modu]
-for key in all_modu:
- if key in BUILD_TARGETS:
- targets.append(key)
-for key in all_exe:
- if key in BUILD_TARGETS:
- exe_targets.append(key)
-
-if targets == [] and 'gsl' not in BUILD_TARGETS and len(exe_targets) == 0:
- targets = all_modu
-if 'all' in BUILD_TARGETS:
- targets = all_modu
-
-for mod in targets:
- mod_env = env.Clone()
- mod_env.VariantDir('build/' + mod, '.')
- mod_env['SWIGFLAGS'] = SWIG_CPP_FLAGS + ' -Isrc' # -Isrc for %include interface files under src
- mod_env['SWIGOUTDIR'] = 'build/%s/src' % mod
- info = ModuInfo(mod, SIMUPOP_VER, SIMUPOP_REV)
- mod_env.Command('build/%s/src/swigpyrun.h' % mod, None, ['swig %s $TARGET' % SWIG_RUNTIME_FLAGS])
- if 'op' in mod:
- common_lib = op_common_lib
- else:
- common_lib = std_common_lib
- mod_wrapper = mod_env.StaticLibrary(
- target = 'build/%s/src/_simuPOP_wrap_%s' % (mod, mod),
- source = ['build/%s/src/simuPOP_%s.i' % (mod, mod)],
- SHLIBPREFIX = "",
- SHLIBSUFFIX = so_ext,
- SHLINKFLAGS = comp.ldflags_shared,
- LIBPATH = info['library_dirs'] + extra_path,
- CPPPATH = [python_inc_dir, '.', 'src', 'build'] + info['include_dirs'],
- CPPDEFINES = convert_def(info['define_macros']),
- # No extra-flags because it can have -Werror and we cannot guarantee that there is
- # no warning coming out of the wrapper code.
- CCFLAGS = opt.split() + info['extra_compile_args'] + comp.compile_options + wrapper_flags,
- CPPFLAGS = ' '.join([basicflags, ccshared])
- )
- mod_lib = mod_env.SharedLibrary(
- target = 'build/%s/_simuPOP_%s' % (mod, mod),
- source = ['build/%s/src/%s' % (mod, x) for x in SOURCE_FILES],
- LIBS = info['libraries'] + [mod_wrapper, common_lib],
- SHLIBPREFIX = "",
- SHLIBSUFFIX = so_ext,
- SHLINKFLAGS = comp.ldflags_shared,
- LIBPATH = info['library_dirs'] + extra_path,
- CPPPATH = [python_inc_dir, '.', 'src', 'build'] + info['include_dirs'],
- CPPDEFINES = convert_def(info['define_macros']),
- CCFLAGS = info['extra_compile_args'] + comp.compile_options + extra_flags,
- CPPFLAGS = ' '.join([basicflags, ccshared, opt])
- )
- env.Depends('build/%s/src/swigpyrun.h' % mod, 'build/%s/src/utility.cpp' % mod)
- Alias(mod, mod_lib)
- Alias('all', mod_lib)
- pyInstall = env.InstallAs(os.path.join(dest_dir, 'simuPOP_%s.py' % mod),
- 'build/%s/src/simuPOP_%s.py' % (mod, mod))
- libInstall = env.InstallAs(os.path.join(dest_dir, '_simuPOP_%s%s' % (mod, so_ext)),
- mod_lib[0])
- mod_env.Depends(pyInstall, libInstall)
- Alias('install', libInstall)
- Alias('install', pyInstall)
-
-# module executable, for testing purposes.
-if os.name == 'nt':
- pylibs = ['python%d%d' % (sys.version_info[0], sys.version_info[1])]
-else:
- pylibs = ['python%d.%d' % (sys.version_info[0], sys.version_info[1]), 'util']
-
-for mod in exe_targets:
- mod_name = mod.split('_')[0]
- exe_env = env.Clone()
- exe_env.VariantDir('build/' + mod, '.')
- exe_env['SWIGFLAGS'] = SWIG_CPP_FLAGS + ' -Isrc' # -Isrc for %include interface files under src
- exe_env['SWIGOUTDIR'] = 'build/%s/src' % mod
- info = ModuInfo(mod_name, SIMUPOP_VER, SIMUPOP_REV)
- if 'op' in mod:
- common_lib = op_common_lib
- else:
- common_lib = std_common_lib
- mod_executable = exe_env.Program(
- target = 'build/simuPOP_%s' % mod,
- source = ['build/%s/src/%s' % (mod, x) for x in SOURCE_FILES],
- LIBS = info['libraries'] + [common_lib] + pylibs,
- LIBPATH = info['library_dirs'] + extra_path,
- CPPPATH = [python_inc_dir, '.', 'src', 'build'] + info['include_dirs'],
- CPPDEFINES = convert_def(info['define_macros'] + [('STANDALONE_EXECUTABLE', None)]),
- CCFLAGS = info['extra_compile_args'] + comp.compile_options + extra_flags,
- CPPFLAGS = ' '.join([basicflags, ccshared, opt])
- )
- Alias(mod, mod_executable)
- Alias('debug', mod_executable)
-
-
-env.Install(pylib_dir, 'simuOpt.py')
-Alias('install', pylib_dir)
-for pyfile in ['__init__.py', 'utils.py', 'plotter.py', 'sampling.py', 'sandbox.py', 'demography.py']:
- env.Install(dest_dir, 'src/%s' % pyfile)
- Alias('install', dest_dir)
-
-for datafile in PACKAGE_DATA:
- env.Install(dest_dir, 'src/%s' % datafile)
- Alias('install', dest_dir)
-
-Default('install')
Added: trunk/test/sample_64_ba_v3.pop
===================================================================
(Binary files differ)
Index: trunk/test/sample_64_ba_v3.pop
===================================================================
--- trunk/test/sample_64_ba_v3.pop 2015-08-31 04:46:59 UTC (rev 4990)
+++ trunk/test/sample_64_ba_v3.pop 2015-08-31 15:19:38 UTC (rev 4991)
Property changes on: trunk/test/sample_64_ba_v3.pop
___________________________________________________________________
Added: svn:mime-type
## -0,0 +1 ##
+application/octet-stream
\ No newline at end of property
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
2 messages has been excluded from this view by a project administrator.