Introduction

simuPOP is a forward-based population genetics simulation program. Unlike coalescence-based simulation programs, simuPOP evolves population(s) forward in time -- subject to arbitrary number of gentic and environmental forces (mutation, recombination, migration, population size change etc.). Easy simulations like most models in standard population genetics books can be setup easily, whereas very complicated simulations such as spreading of complex diseases, ancient out-of-africa migrations can be built step by step by adding appropriate operators (objects that work on populations).

simuPOP can be used at two levels: a normal user can grab one of the pre-defined simulation templates and run it with appropriate parameters; and an advanced user can work interactively or write simuPOP/Python scripts to construct arbitrarily complex evoluationary processes. For example, you can run

  > simuComplexDisease.py 

from a command line and starting a simulation of the evolution of a complex disease.

On the other hand. simuPOP is a R or matlab programming environment in which users can create, manipulate and evolve populations, calculate and visualize population statistics, and test mapping methods. The full power of simuPOP (and Python) can be used to simulate arbitrarily complex evolutionary scenarios. For example, you can write the following script:

>>> from simuPOP import *

>>> simu = simulator(
...     population(size=1000, ploidy=2, loci=[2]),
...     randomMating(),
...     rep = 3)
>>> simu.evolve(
...     preOp = [initByValue([1,2,2,1]) ],  
...     ops = [
...         recombinator( rate=0.1),
...         alleleCounter( haplotypes=[[0,1]]),
...         pyEval(expr = r"str(haploNum['0-1']['1-2'])+'t'",
...             headers=["hapt"]),
...         output("n", rep=REP_LAST)
...         ],
...     end=5
... )
0_hap	1_hap	2_hap	
888.0	874.0	892.0	
835.0	780.0	819.0	
784.0	724.0	783.0	
762.0	688.0	757.0	
769.0	654.0	723.0	
>>>  

This example demonstrates the change of haplotype number when recombination is in effect.

• The import line import simuPOP module (output suppressed).
• simulator creates a simulator from a population created by population function. The population is diploid (ploidy=2), has 1000 individuals (size=1000) each has two loci on the first chromosome (loci=[2]). The simulator has three copies of this population (rep=3) and will evolve through random mating (randomMating()).
• simu.evolve evolves these populations subject to the following operators.
  • preop=[initByValue]: operators in parameter preop (accept a list of operators) will be applied to the populations at the beginning of evolution. initByValue is an initializer that set the same genotype to all individuals. In this case, everyone will have genotype 12/21 (1 2 on the first chromosome and 2 1 on the second copy of the chromosome) so linkage disequilibrium is 0.25 (maximum possible value).
  • operators in ops parameter will be applied to all populations at each generation. (Not exactly, operators can be inactive at certain generations.)
  • recombinator is a during-mating operator that recombine chromosomes with probability 0.1 (an unrealistically high value) during mating.
  • alleleCounter is a post-mating operator that count haplotype 12 and 21 at locus 0 and locus 1 respectively. The result will be saved in the local namespace of each replicate. (_vars[0], _vars[1], _vars[2] if you are curious.)
  • pyEval accepts any python expression, evaluate it in each replicates' local namespace and return the result. In this example, pyEval get the value of a haplotype '0-1|1-2', print it with a trailing tab.
  • output is a simple operator that output any string it is given. Here, it output a new line at the last replicate. Without this, all output will be in one line.
  • end=5 means evolve 5 generations.

The output is a table of haplotype number for each replicate at each generation. All simuPOP scripts will have similar steps. You can add more operators to the ops list to build more complicated simulations. Obvious choices are mutator, migrator, or some proper visualizer (another pyEval operator to call Python's plotting functions) to plot the dynamics of variables.