tests/cookbook-baba-simple.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## (&#x1F4D7;) Cookbook: `ipyrad.analysis.baba` \n",
    "\n",
    "This notebook demonstrates how to use the `baba` module in `ipyrad` to test for admixture and introgression. The code is written in Python and is best implemented in a Jupyter-notebook like this one. Analyses can be split up among many computing cores. Finally, there are some simple plotting methods to accompany tests which I show below. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### import packages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  ipyrad v.0.5.15\n",
      "  ipyparallel v.5.2.0\n",
      "  numpy v.1.11.2\n"
     ]
    }
   ],
   "source": [
    "## imports\n",
    "import numpy as np\n",
    "import ipyrad as ip\n",
    "import ipyrad.analysis.baba as baba\n",
    "import ipyparallel as ipp\n",
    "\n",
    "## print versions\n",
    "print \"  ipyrad v.{}\".format(ip.__version__)\n",
    "print \"  ipyparallel v.{}\".format(ipp.__version__)\n",
    "print \"  numpy v.{}\".format(np.__version__)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### ipcluster setup to run parallel code\n",
    "You can find more details about this in our [ipyparallel tutorial]."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "## start 4 engines by running the commented line below in a separate bash terminal. \n",
    "##\n",
    "## ipcluster start --n=4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  host compute node: [4 cores] on oud\n"
     ]
    }
   ],
   "source": [
    "## connect to client and print connection summary\n",
    "ipyclient = ipp.Client()\n",
    "print ip.cluster_info(client=ipyclient)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Get default (example) 12 taxon tree \n",
    "The function  `baba.Tree()`  takes a  `newick`  string as an argument, or, if no value is passed, it constructs the default 12 taxon tree shown below. The  `Tree`  Class object holds a representation of the tree, which can be used for plotting, as well as for describing a demographic model for simulating data, as we'll show below. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div align=\"center\" class=\"toyplot\" id=\"t3296c26a307e4f83891a63a0c24429c6\"><svg class=\"toyplot-canvas-Canvas\" height=\"300.0px\" id=\"t5df3d57c59ba4afcb02c4b65b7439a67\" preserveAspectRatio=\"xMidYMid meet\" style=\"background-color:transparent;fill:rgb(16.1%,15.3%,14.1%);fill-opacity:1.0;font-family:Helvetica;font-size:12px;opacity:1.0;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:1.0\" viewBox=\"0 0 400.0 300.0\" width=\"400.0px\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:toyplot=\"http://www.sandia.gov/toyplot\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><g class=\"toyplot-coordinates-Cartesian\" id=\"t38979d1fcd21497a9ce28c74d1f3cc48\"><clipPath id=\"t52f6f6688d95409ea5d8c335c31724da\"><rect height=\"220.0\" width=\"320.0\" x=\"40.0\" y=\"40.0\"></rect></clipPath><g clip-path=\"url(#t52f6f6688d95409ea5d8c335c31724da)\"><g class=\"toyplot-mark-Graph\" id=\"t0e71098d6dd54d90b7605f650482b836\"><g class=\"toyplot-Edges\"><path d=\"M 302.272727273 170.0 L 322.727272727 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 281.818181818 210.0 L 295.454545455 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 281.818181818 210.0 L 268.181818182 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 302.272727273 170.0 L 281.818181818 210.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 326.136363636 130.0 L 302.272727273 170.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 244.318181818 50.0 L 326.136363636 130.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 217.045454545 130.0 L 240.909090909 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 193.181818182 170.0 L 213.636363636 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 172.727272727 210.0 L 186.363636364 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 172.727272727 210.0 L 159.090909091 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 193.181818182 170.0 L 172.727272727 210.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 217.045454545 130.0 L 193.181818182 170.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 162.5 90.0 L 217.045454545 130.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 107.954545455 130.0 L 131.818181818 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 84.0909090909 170.0 L 104.545454545 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 63.6363636364 210.0 L 77.2727272727 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 63.6363636364 210.0 L 50.0 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 84.0909090909 170.0 L 63.6363636364 210.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 107.954545455 130.0 L 84.0909090909 170.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 162.5 90.0 L 107.954545455 130.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 244.318181818 50.0 L 162.5 90.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 326.136363636 130.0 L 350.0 250.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path></g><g class=\"toyplot-Vertices\"><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"244.31818181818181\" cy=\"50.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"326.13636363636363\" cy=\"130.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"302.27272727272725\" cy=\"170.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"281.81818181818181\" cy=\"210.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"162.5\" cy=\"89.999999999999986\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"217.04545454545453\" cy=\"130.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"193.18181818181819\" cy=\"170.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"172.72727272727272\" cy=\"210.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"107.95454545454545\" cy=\"130.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"84.090909090909093\" cy=\"170.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"63.63636363636364\" cy=\"210.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"295.45454545454544\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"268.18181818181819\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"322.72727272727275\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"350.0\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"186.36363636363637\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"159.09090909090909\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"213.63636363636363\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"240.90909090909091\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"77.27272727272728\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"50.0\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"104.54545454545455\" cy=\"250.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"131.81818181818181\" cy=\"250.0\" r=\"10.0\"></circle></g></g><g class=\"toyplot-Labels\"><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(244.31818181818181,50.0)\"><tspan style=\"dominant-baseline:inherit\">0</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(326.13636363636363,130.0)\"><tspan style=\"dominant-baseline:inherit\">1</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(302.27272727272725,170.0)\"><tspan style=\"dominant-baseline:inherit\">2</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(281.81818181818181,210.0)\"><tspan style=\"dominant-baseline:inherit\">3</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(162.5,89.999999999999986)\"><tspan style=\"dominant-baseline:inherit\">4</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(217.04545454545453,130.0)\"><tspan style=\"dominant-baseline:inherit\">5</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(193.18181818181819,170.0)\"><tspan style=\"dominant-baseline:inherit\">6</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(172.72727272727272,210.0)\"><tspan style=\"dominant-baseline:inherit\">7</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(107.95454545454545,130.0)\"><tspan style=\"dominant-baseline:inherit\">8</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(84.090909090909093,170.0)\"><tspan style=\"dominant-baseline:inherit\">9</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(63.63636363636364,210.0)\"><tspan style=\"dominant-baseline:inherit\">10</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(295.45454545454544,250.0)\"><tspan style=\"dominant-baseline:inherit\">a</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(268.18181818181819,250.0)\"><tspan style=\"dominant-baseline:inherit\">b</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(322.72727272727275,250.0)\"><tspan style=\"dominant-baseline:inherit\">c</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(350.0,250.0)\"><tspan style=\"dominant-baseline:inherit\">d</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(186.36363636363637,250.0)\"><tspan style=\"dominant-baseline:inherit\">e</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(159.09090909090909,250.0)\"><tspan style=\"dominant-baseline:inherit\">f</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(213.63636363636363,250.0)\"><tspan style=\"dominant-baseline:inherit\">g</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(240.90909090909091,250.0)\"><tspan style=\"dominant-baseline:inherit\">h</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(77.27272727272728,250.0)\"><tspan style=\"dominant-baseline:inherit\">i</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(50.0,250.0)\"><tspan style=\"dominant-baseline:inherit\">j</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(104.54545454545455,250.0)\"><tspan style=\"dominant-baseline:inherit\">k</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(131.81818181818181,250.0)\"><tspan style=\"dominant-baseline:inherit\">l</tspan></text></g></g></g></g></svg><div class=\"toyplot-interactive\"></div></div>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "## create a Tree class object with the default tree\n",
    "tre = baba.Tree()\n",
    "tre.show(width=400, height=300);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Or get any arbitrary tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div align=\"center\" class=\"toyplot\" id=\"t42a429e521e143508938264e494d3638\"><svg class=\"toyplot-canvas-Canvas\" height=\"200.0px\" id=\"t9377d192f94a427ca8c95a7b73c164d8\" preserveAspectRatio=\"xMidYMid meet\" style=\"background-color:transparent;fill:rgb(16.1%,15.3%,14.1%);fill-opacity:1.0;font-family:Helvetica;font-size:12px;opacity:1.0;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:1.0\" viewBox=\"0 0 200.0 200.0\" width=\"200.0px\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:toyplot=\"http://www.sandia.gov/toyplot\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><g class=\"toyplot-coordinates-Cartesian\" id=\"t6cfea81f17bf4853bed923c9fd4d0af2\"><clipPath id=\"t0d568f03427340f19baf2ecc21805b36\"><rect height=\"120.0\" width=\"120.0\" x=\"40.0\" y=\"40.0\"></rect></clipPath><g clip-path=\"url(#t0d568f03427340f19baf2ecc21805b36)\"><g class=\"toyplot-mark-Graph\" id=\"t007f11bcd884459d91b3866a7be80b07\"><g class=\"toyplot-Edges\"><path d=\"M 75.0 100.0 L 100.0 150.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 75.0 100.0 L 50.0 150.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 112.5 50.0 L 75.0 100.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 112.5 50.0 L 150.0 150.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path></g><g class=\"toyplot-Vertices\"><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"112.5\" cy=\"50.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"75.0\" cy=\"100.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"100.0\" cy=\"150.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"50.0\" cy=\"150.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"150.0\" cy=\"150.0\" r=\"10.0\"></circle></g></g><g class=\"toyplot-Labels\"><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(112.5,50.0)\"><tspan style=\"dominant-baseline:inherit\">0</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(75.0,100.0)\"><tspan style=\"dominant-baseline:inherit\">1</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(100.0,150.0)\"><tspan style=\"dominant-baseline:inherit\">a</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(50.0,150.0)\"><tspan style=\"dominant-baseline:inherit\">b</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(150.0,150.0)\"><tspan style=\"dominant-baseline:inherit\">c</tspan></text></g></g></g></g></svg><div class=\"toyplot-interactive\"></div></div>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "## pass a newick string to the Tree class object\n",
    "tre = baba.Tree(newick=\"((a,b),c);\")\n",
    "tre.show(width=200, height=200);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### get a tree with an admixture edge\n",
    "Shown by an arrow pointing backwards from source to sink (indicating that geneflow should be entered as occurring backwards in time). This is how the program `msprime` will interpret it. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div align=\"center\" class=\"toyplot\" id=\"tae7587270788481b892f8f01e655bbc3\"><svg class=\"toyplot-canvas-Canvas\" height=\"250.0px\" id=\"t0ccec0fc3d2b46c1ac05081d731ce6fd\" preserveAspectRatio=\"xMidYMid meet\" style=\"background-color:transparent;fill:rgb(16.1%,15.3%,14.1%);fill-opacity:1.0;font-family:Helvetica;font-size:12px;opacity:1.0;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:1.0\" viewBox=\"0 0 250.0 250.0\" width=\"250.0px\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:toyplot=\"http://www.sandia.gov/toyplot\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><g class=\"toyplot-coordinates-Cartesian\" id=\"t45f8f03c09c44b5fa6c001f5e1f4ebc2\"><clipPath id=\"t64c9948f2d0c49f2bd0b582f81a96aef\"><rect height=\"170.0\" width=\"170.0\" x=\"40.0\" y=\"40.0\"></rect></clipPath><g clip-path=\"url(#t64c9948f2d0c49f2bd0b582f81a96aef)\"><g class=\"toyplot-mark-Graph\" id=\"t2c386d1e9a4f4720bcf76a8ada64f1ce\"><g class=\"toyplot-Edges\"><path d=\"M 181.25 150.0 L 162.5 200.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 139.0625 50.0 L 181.25 150.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 96.875 100.0 L 125.0 200.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 68.75 150.0 L 87.5 200.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 68.75 150.0 L 50.0 200.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 96.875 100.0 L 68.75 150.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 139.0625 50.0 L 96.875 100.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path><path d=\"M 181.25 150.0 L 200.0 200.0\" style=\"fill:none;stroke:rgb(16.1%,15.3%,14.1%);stroke-opacity:1.0;stroke-width:3.0\"></path></g><g class=\"toyplot-Vertices\"><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"139.0625\" cy=\"50.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"181.25\" cy=\"150.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"96.875\" cy=\"100.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"68.75\" cy=\"150.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"87.5\" cy=\"200.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"50.0\" cy=\"200.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"125.0\" cy=\"200.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"200.0\" cy=\"200.0\" r=\"10.0\"></circle></g><g class=\"toyplot-Datum\" style=\"fill:rgb(40%,76.1%,64.7%);fill-opacity:1.0;opacity:1.0;stroke:rgb(40%,76.1%,64.7%);stroke-opacity:1.0\"><circle cx=\"162.5\" cy=\"200.0\" r=\"10.0\"></circle></g></g><g class=\"toyplot-Labels\"><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(139.0625,50.0)\"><tspan style=\"dominant-baseline:inherit\">0</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(181.25,150.0)\"><tspan style=\"dominant-baseline:inherit\">1</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(96.875,100.0)\"><tspan style=\"dominant-baseline:inherit\">2</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(68.75,150.0)\"><tspan style=\"dominant-baseline:inherit\">3</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(87.5,200.0)\"><tspan style=\"dominant-baseline:inherit\">a</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(50.0,200.0)\"><tspan style=\"dominant-baseline:inherit\">b</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(125.0,200.0)\"><tspan style=\"dominant-baseline:inherit\">c</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(200.0,200.0)\"><tspan style=\"dominant-baseline:inherit\">d</tspan></text><text class=\"toyplot-Datum\" style=\"dominant-baseline:middle;font-size:18px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(162.5,200.0)\"><tspan style=\"dominant-baseline:inherit\">e</tspan></text></g></g><g class=\"toyplot-mark-Plot\" id=\"te318f22b8b5c43d099e623a417cd214c\" style=\"fill:none;stroke-dasharray:2, 2;stroke-width:3\"><g class=\"toyplot-Series\"><path d=\"M 171.875 175.00000000000003 L 110.9375 150.0\" style=\"stroke:rgb(98.8%,55.3%,38.4%);stroke-dasharray:2, 2;stroke-opacity:1.0;stroke-width:3\"></path></g></g><g class=\"toyplot-mark-Plot\" id=\"t6b7af0b3d6934f76b786105ceb44a536\" style=\"fill:none;stroke-dasharray:2, 2;stroke-width:3\"><g class=\"toyplot-Series\"><path d=\"M 82.8125 125.0 L 160.15625 100.0\" style=\"stroke:rgb(98.8%,55.3%,38.4%);stroke-dasharray:2, 2;stroke-opacity:1.0;stroke-width:3\"></path></g></g></g><g class=\"toyplot-coordinates-Axis\" id=\"tfbf8c40a478b474686d81fca29ab72e5\" transform=\"translate(50.0,200.0)rotate(-90.0)translate(0,-10.0)\"><line style=\"\" x1=\"0\" x2=\"150.0\" y1=\"0\" y2=\"0\"></line><g><text style=\"dominant-baseline:alphabetic;font-size:10px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(0.0,-6)\"><tspan style=\"dominant-baseline:inherit\">0</tspan></text><text style=\"dominant-baseline:alphabetic;font-size:10px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(50.0,-6)\"><tspan style=\"dominant-baseline:inherit\">1</tspan></text><text style=\"dominant-baseline:alphabetic;font-size:10px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(100.0,-6)\"><tspan style=\"dominant-baseline:inherit\">2</tspan></text><text style=\"dominant-baseline:alphabetic;font-size:10px;font-weight:normal;stroke:none;text-anchor:middle\" transform=\"translate(150.0,-6)\"><tspan style=\"dominant-baseline:inherit\">3</tspan></text></g><g class=\"toyplot-coordinates-Axis-coordinates\" style=\"visibility:hidden\" transform=\"\"><line style=\"stroke:rgb(43.9%,50.2%,56.5%);stroke-opacity:1.0;stroke-width:1.0\" x1=\"0\" x2=\"0\" y1=\"3.0\" y2=\"-4.5\"></line><text style=\"alignment-baseline:hanging;fill:rgb(43.9%,50.2%,56.5%);fill-opacity:1.0;font-size:10px;font-weight:normal;stroke:none;text-anchor:middle\" x=\"0\" y=\"6\"></text></g></g></g></svg><div class=\"toyplot-interactive\"><ul class=\"toyplot-mark-popup\" onmouseleave=\"this.style.visibility='hidden'\" style=\"background:rgba(0%,0%,0%,0.75);border:0;border-radius:6px;color:white;cursor:default;list-style:none;margin:0;padding:5px;position:fixed;visibility:hidden\">\n",
       "            <li class=\"toyplot-mark-popup-title\" style=\"color:lightgray;cursor:default;padding:5px;list-style:none;margin:0\"></li>\n",
       "            <li class=\"toyplot-mark-popup-save-csv\" onmouseout=\"this.style.color='white';this.style.background='steelblue'\" onmouseover=\"this.style.color='steelblue';this.style.background='white'\" style=\"border-radius:3px;padding:5px;list-style:none;margin:0\">\n",
       "                Save as .csv\n",
       "            </li>\n",
       "        </ul><script>\n",
       "        (function()\n",
       "        {\n",
       "          var data_tables = [{\"title\": \"Plot Data\", \"names\": [\"x\", \"y0\"], \"id\": \"te318f22b8b5c43d099e623a417cd214c\", \"columns\": [[3.25, 1.625], [0.5, 1.0]], \"filename\": \"toyplot\"}, {\"title\": \"Plot Data\", \"names\": [\"x\", \"y0\"], \"id\": \"t6b7af0b3d6934f76b786105ceb44a536\", \"columns\": [[0.875, 2.9375], [1.5, 2.0]], \"filename\": \"toyplot\"}];\n",
       "\n",
       "          function save_csv(data_table)\n",
       "          {\n",
       "            var uri = \"data:text/csv;charset=utf-8,\";\n",
       "            uri += data_table.names.join(\",\") + \"\\n\";\n",
       "            for(var i = 0; i != data_table.columns[0].length; ++i)\n",
       "            {\n",
       "              for(var j = 0; j != data_table.columns.length; ++j)\n",
       "              {\n",
       "                if(j)\n",
       "                  uri += \",\";\n",
       "                uri += data_table.columns[j][i];\n",
       "              }\n",
       "              uri += \"\\n\";\n",
       "            }\n",
       "            uri = encodeURI(uri);\n",
       "\n",
       "            var link = document.createElement(\"a\");\n",
       "            if(typeof link.download != \"undefined\")\n",
       "            {\n",
       "              link.href = uri;\n",
       "              link.style = \"visibility:hidden\";\n",
       "              link.download = data_table.filename + \".csv\";\n",
       "\n",
       "              document.body.appendChild(link);\n",
       "              link.click();\n",
       "              document.body.removeChild(link);\n",
       "            }\n",
       "            else\n",
       "            {\n",
       "              window.open(uri);\n",
       "            }\n",
       "          }\n",
       "\n",
       "          function open_popup(data_table)\n",
       "          {\n",
       "            return function(e)\n",
       "            {\n",
       "              var popup = document.querySelector(\"#tae7587270788481b892f8f01e655bbc3 .toyplot-mark-popup\");\n",
       "              popup.querySelector(\".toyplot-mark-popup-title\").innerHTML = data_table.title;\n",
       "              popup.querySelector(\".toyplot-mark-popup-save-csv\").onclick = function() { popup.style.visibility = \"hidden\"; save_csv(data_table); }\n",
       "              popup.style.left = (e.clientX - 50) + \"px\";\n",
       "              popup.style.top = (e.clientY - 20) + \"px\";\n",
       "              popup.style.visibility = \"visible\";\n",
       "              e.stopPropagation();\n",
       "              e.preventDefault();\n",
       "            }\n",
       "\n",
       "          }\n",
       "\n",
       "          for(var i = 0; i != data_tables.length; ++i)\n",
       "          {\n",
       "            var data_table = data_tables[i];\n",
       "            var event_target = document.querySelector(\"#\" + data_table.id);\n",
       "            event_target.oncontextmenu = open_popup(data_table);\n",
       "          }\n",
       "        })();\n",
       "        </script><script>\n",
       "        (function()\n",
       "        {\n",
       "            function _sign(x)\n",
       "            {\n",
       "                return x < 0 ? -1 : x > 0 ? 1 : 0;\n",
       "            }\n",
       "\n",
       "            function _mix(a, b, amount)\n",
       "            {\n",
       "                return ((1.0 - amount) * a) + (amount * b);\n",
       "            }\n",
       "\n",
       "            function _log(x, base)\n",
       "            {\n",
       "                return Math.log(Math.abs(x)) / Math.log(base);\n",
       "            }\n",
       "\n",
       "            function _in_range(a, x, b)\n",
       "            {\n",
       "                var left = Math.min(a, b);\n",
       "                var right = Math.max(a, b);\n",
       "                return left <= x && x <= right;\n",
       "            }\n",
       "\n",
       "            function inside(range, projection)\n",
       "            {\n",
       "                for(var i = 0; i != projection.length; ++i)\n",
       "                {\n",
       "                    var segment = projection[i];\n",
       "                    if(_in_range(segment.range.min, range, segment.range.max))\n",
       "                        return true;\n",
       "                }\n",
       "                return false;\n",
       "            }\n",
       "\n",
       "            function to_domain(range, projection)\n",
       "            {\n",
       "                for(var i = 0; i != projection.length; ++i)\n",
       "                {\n",
       "                    var segment = projection[i];\n",
       "                    if(_in_range(segment.range.bounds.min, range, segment.range.bounds.max))\n",
       "                    {\n",
       "                        if(segment.scale == \"linear\")\n",
       "                        {\n",
       "                            var amount = (range - segment.range.min) / (segment.range.max - segment.range.min);\n",
       "                            return _mix(segment.domain.min, segment.domain.max, amount)\n",
       "                        }\n",
       "                        else if(segment.scale[0] == \"log\")\n",
       "                        {\n",
       "                            var amount = (range - segment.range.min) / (segment.range.max - segment.range.min);\n",
       "                            var base = segment.scale[1];\n",
       "                            return _sign(segment.domain.min) * Math.pow(base, _mix(_log(segment.domain.min, base), _log(segment.domain.max, base), amount));\n",
       "                        }\n",
       "                    }\n",
       "                }\n",
       "            }\n",
       "\n",
       "            function display_coordinates(e)\n",
       "            {\n",
       "                var current = svg.createSVGPoint();\n",
       "                current.x = e.clientX;\n",
       "                current.y = e.clientY;\n",
       "\n",
       "                for(var axis_id in axes)\n",
       "                {\n",
       "                    var axis = document.querySelector(\"#\" + axis_id);\n",
       "                    var coordinates = axis.querySelector(\".toyplot-coordinates-Axis-coordinates\");\n",
       "                    if(coordinates)\n",
       "                    {\n",
       "                        var projection = axes[axis_id];\n",
       "                        var local = current.matrixTransform(axis.getScreenCTM().inverse());\n",
       "                        if(inside(local.x, projection))\n",
       "                        {\n",
       "                            var domain = to_domain(local.x, projection);\n",
       "                            coordinates.style.visibility = \"visible\";\n",
       "                            coordinates.setAttribute(\"transform\", \"translate(\" + local.x + \")\");\n",
       "                            var text = coordinates.querySelector(\"text\");\n",
       "                            text.textContent = domain.toFixed(2);\n",
       "                        }\n",
       "                        else\n",
       "                        {\n",
       "                            coordinates.style.visibility= \"hidden\";\n",
       "                        }\n",
       "                    }\n",
       "                }\n",
       "            }\n",
       "\n",
       "            var root_id = \"tae7587270788481b892f8f01e655bbc3\";\n",
       "            var axes = {\"tfbf8c40a478b474686d81fca29ab72e5\": [{\"domain\": {\"bounds\": {\"max\": Infinity, \"min\": -Infinity}, \"max\": 3.0, \"min\": 0.0}, \"range\": {\"bounds\": {\"max\": Infinity, \"min\": -Infinity}, \"max\": 150.0, \"min\": 0.0}, \"scale\": \"linear\"}]};\n",
       "\n",
       "            var svg = document.querySelector(\"#\" + root_id + \" svg\");\n",
       "            svg.addEventListener(\"click\", display_coordinates);\n",
       "        })();\n",
       "        </script></div></div>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "## store newick tree and a migration event list [source, sink, start, stop, rate]\n",
    "newick = \"(((a,b),c), (d,e));\"\n",
    "events = [['e', 'c', 0, 1, 1e-6], \n",
    "          ['3', '1', 1, 2, 1e-6]]\n",
    "\n",
    "## add admixture events to tree\n",
    "tre = baba.Tree(newick=newick, admix=events)\n",
    "tre.show(width=250, height=250, yaxis=True);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### simulate data on that tree\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "8 [ 1.  0.]\n",
      "7 [ 2.  0.]\n",
      "6 [ 3.  0.]\n",
      "5 [ 4.  0.]\n",
      "4 [ 0.5  1. ]\n",
      "3 [ 3.5  1. ]\n",
      "2 [ 1.25  2.  ]\n",
      "1 [ 2.375  3.   ]\n",
      "0 [ 0.  0.]\n"
     ]
    }
   ],
   "source": [
    "for idx in xrange(tre.verts.shape[0]-1, -1, -1):\n",
    "    print idx, tre.verts[idx-1]\n",
    "    row = tre.verts[idx-1]\n",
    "    \n",
    "    ## get time of the next event\n",
    "    if row[1] == 0:\n",
    "        pass\n",
    "    else:\n",
    "        tau = row[1]\n",
    "    \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "when [ 3.  2.  1.  1.]\n",
      "who [0 1 2 3]\n"
     ]
    }
   ],
   "source": [
    "## when do coalescent events occur?\n",
    "who = tre.verts[:, 1] > 0\n",
    "print 'when', tre.verts[who, 1]\n",
    "print 'who', np.where(who)[0]\n",
    "\n",
    "## when do migration events occur\n",
    "## ...\n",
    "\n",
    "for time in when:\n",
    "    ## does migration change during the next epoch\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 2.375  3.   ]\n",
      " [ 1.25   2.   ]\n",
      " [ 3.5    1.   ]\n",
      " [ 0.5    1.   ]\n",
      " [ 4.     0.   ]\n",
      " [ 3.     0.   ]\n",
      " [ 2.     0.   ]\n",
      " [ 1.     0.   ]\n",
      " [ 0.     0.   ]]\n",
      "[[1 8]\n",
      " [0 1]\n",
      " [2 6]\n",
      " [3 4]\n",
      " [3 5]\n",
      " [2 3]\n",
      " [0 2]\n",
      " [1 7]]\n"
     ]
    }
   ],
   "source": [
    "for node in tre.tree.traverse(\"postorder\"):\n",
    "    ## check for migration\n",
    "    if node.name in [i[0] for i in events]:\n",
    "        node\n",
    "        #ms.MigrationRateChange(time=0, rate=m_C_B, matrix_index=(1, 2)),\n",
    "        #ms.MigrationRateChange(time=Taus[1], rate=0),"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Simulate data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [conda root]",
   "language": "python",
   "name": "conda-root-py"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}