content/08_phylogenetic_trees/birth_death_models.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Implementing a Birth-Death model of Speciation and Extinction\n",
    "\n",
    "Simple models are extremely useful in bioinformatics. Although a simple model will not generally capture all aspects of a biological process, it can provide a very useful baseline: if available biological data is well explained by the simple model, then no further explanation is demanded. On the other hand, if the simple model does a poor job of explaining the observed data, then further investigations, and perhaps a more complex model may be justified. It is by digging in deeply to such cases - where something 'weird' or surprising is going on - that many interesting discoveries in biology have been uncovered\n",
    "\n",
    "A simple way to model the evolution of species diversity is through a birth-death process. In such a process, species diverge (forming two new species) or go extinct at some rate. In this section we will implement a birth-death model of species diversity using our PhyloNode object.\n",
    "\n",
    "In our model we start with one initial species representing the common ancestor of the tree. New species form by speciation events which cause two child species to appear. Species can also be removed by extinction events. We can then add a chance of speciating or going extinct at each point in time.\n",
    "\n",
    "To do so we'll have to add a speciation function that automatically splits a lineage in two, a property to keep track of which nodes are extinct, and an update function that updates a node on the tree, causing it to possibly speciate or go extinct.\n",
    "\n",
    "Finally, outside of our class we'll want to write a for loop to run the simulation for a given number of iterations and report the results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from random import random, choice\n",
    "\n",
    "class PhyloNode():\n",
    "    \"\"\"A node on a phylogenetic tree\"\"\"\n",
    "\n",
    "    def __init__(self,children=None,parent=None,\\\n",
    "      name = None):\n",
    "      \"\"\"Initiate a node on a phylogenetic tree\n",
    "      children -- a list of PhyloNode objects\n",
    "        descending immediately from this node\n",
    "      name -- a string for the name of the node\n",
    "      parent -- a single PhyloNode object for the parent\n",
    "      \"\"\"\n",
    "      self.Name = name\n",
    "      self.Children = []\n",
    "      if children:\n",
    "        self.Children.extend(children)\n",
    "      self.Parent = parent\n",
    "      self.Extinct = False\n",
    "\n",
    "    def isTip(self):\n",
    "      \"\"\"Return True if the node is a tip\"\"\"\n",
    "      if not self.Children: #capture None or []\n",
    "        return True\n",
    "      else:\n",
    "        return False\n",
    "\n",
    "    def isRoot(self):\n",
    "      \"\"\"Return True if the node is the root of the whole tree\"\"\"\n",
    "      if not self.Parent:\n",
    "        return True\n",
    "      else:\n",
    "        return False\n",
    "\n",
    "    def getDescendants(self):\n",
    "        \"\"\"Return a list of PhyloNodes descending from the current node\"\"\"\n",
    "\n",
    "        if self.isTip():\n",
    "            print(self.Name,\" is a tip ... returning []\")\n",
    "            return []\n",
    "\n",
    "        descendants = self.Children or []\n",
    "        for c in self.Children:\n",
    "            #The set of descendants is described\n",
    "            #by the descendants of all the nodes\n",
    "            #immediate children.\n",
    "            if not c.isTip():\n",
    "                child_descendants = c.getDescendants()\n",
    "                descendants.extend([c for c in child_descendants if c not in descendants])\n",
    "\n",
    "            #Side note: this will fail on enormous trees\n",
    "            #due to the recursion limit. Not normally a problem though.\n",
    "        return descendants\n",
    "\n",
    "    def getAncestors(self):\n",
    "        \"\"\"Return the ancestors of the given node\"\"\"\n",
    "        if self.isRoot():\n",
    "            return None\n",
    "\n",
    "        ancestors = [self.Parent]\n",
    "        parents_ancestors = self.Parent.getAncestors()\n",
    "        if parents_ancestors:\n",
    "            ancestors.extend(parents_ancestors)\n",
    "        return ancestors\n",
    "\n",
    "    def getRoot(self):\n",
    "        \"\"\"Return the root node\"\"\"\n",
    "        curr_node = self\n",
    "        while not curr_node.isRoot():\n",
    "            curr_node = curr_node.Parent\n",
    "        return curr_node\n",
    "\n",
    "    def addChild(self,child):\n",
    "        \"\"\"Attach a child node\"\"\"\n",
    "        if child not in self.Children:\n",
    "            self.Children.append(child)\n",
    "        child.Parent = self\n",
    " \n",
    "    def addParent(self,parent):\n",
    "        \"\"\"Attach a parent node\"\"\"\n",
    "        self.Parent = parent\n",
    "        parent.Children.append(self)\n",
    "    \n",
    "    def getMRCA(self,other):\n",
    "        \"\"\"Return PhyloNode for the most recent common ancestor with other\n",
    "        other -- another PhyloNode object in the same tree\n",
    "        \"\"\"\n",
    "        most_recent_common_ancestor = None\n",
    "        #Cache the ancestors of other since\n",
    "        #we'll refer to it often\n",
    "        other_ancestors = other.getAncestors()\n",
    "        self_ancestors = self.getAncestors()\n",
    "        #First check for the trivial case\n",
    "        #where one node is root (but be sure they're on the same tree)\n",
    "        if self.isRoot() and self in other_ancestors:\n",
    "            return self\n",
    "\n",
    "        if other.isRoot() and other in self_ancestors:\n",
    "            return other\n",
    "\n",
    "        for a in self_ancestors:\n",
    "            #Note these will be in order of relatedness\n",
    "            if a in other_ancestors:\n",
    "                #End the loop the first time this happens\n",
    "                #since we want the most\n",
    "                #recent common ancestor\n",
    "                most_recent_common_ancestor = a\n",
    "                break\n",
    "\n",
    "        if not most_recent_common_ancestor:\n",
    "            raise ValueError(\"No common ancestor found for \",self.Name,\" and \",other.Name,\\\n",
    "            \". Are they on the same tree?\")\n",
    "\n",
    "        return most_recent_common_ancestor\n",
    "    \n",
    "    def speciate(self):\n",
    "        \"\"\"Add two children descending from this node\"\"\"\n",
    "        if self.Children:\n",
    "            raise ValueError(\"Internal nodes can't speciate\")\n",
    "        if self.Extinct:\n",
    "            raise ValueError(\"Extinct nodes can't speciate\")\n",
    "\n",
    "        child1_name = self.Name + \"A\"\n",
    "        child1=PhyloNode(name=child1_name)\n",
    "        self.addChild(child1)\n",
    "        child2_name = self.Name + \"B\"\n",
    "        child2=PhyloNode(name=child2_name)\n",
    "        self.addChild(child2)\n",
    "\n",
    "\n",
    "    def update(self,speciation_chance=0.25,\\\n",
    "        extinction_chance=0.25):\n",
    "        \"\"\"Update the node by speciation or extinction\"\"\"\n",
    "        if not self.isTip():\n",
    "            print(\"Not updating node {name} - it's not a tip\".format(name=self.Name))\n",
    "            return None\n",
    "\n",
    "        if random() < extinction_chance:\n",
    "            print(self.Name,\" goes extinct!\")\n",
    "            self.Extinct = True\n",
    "\n",
    "        if self.Extinct:\n",
    "            print(\"Not updating node {name} - it's extinct\".format(name=self.Name))\n",
    "            return None\n",
    "        print(\"Updating node:{name}\".format(name=self.Name))\n",
    "\n",
    "        if random() < speciation_chance:\n",
    "            #Speciate\n",
    "            print(\"Node {name} speciates!\".format(name=self.Name))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Running the birth-death model\n",
    "\n",
    "Now that we have a way to update the tree at each timestep, we can build a simple tree and use a for-loop to run the birth-death model in a loop."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Root node name: root\n",
      "================================================================================\n",
      "Beginning simulation of time: 1\n",
      "================================================================================\n",
      "About to update node: A\n",
      "A  goes extinct!\n",
      "Not updating node A - it's extinct\n",
      "About to update node: B\n",
      "Updating node:B\n",
      "Nodes in tree at time  1 : 2\n",
      "Tree tips: ['A', 'B']\n",
      "================================================================================\n",
      "Beginning simulation of time: 2\n",
      "================================================================================\n",
      "About to update node: A\n",
      "A  goes extinct!\n",
      "Not updating node A - it's extinct\n",
      "About to update node: B\n",
      "Updating node:B\n",
      "Nodes in tree at time  2 : 2\n",
      "Tree tips: ['A', 'B']\n",
      "================================================================================\n",
      "Beginning simulation of time: 3\n",
      "================================================================================\n",
      "About to update node: A\n",
      "Not updating node A - it's extinct\n",
      "About to update node: B\n",
      "Updating node:B\n",
      "Nodes in tree at time  3 : 2\n",
      "Tree tips: ['A', 'B']\n",
      "================================================================================\n",
      "Beginning simulation of time: 4\n",
      "================================================================================\n",
      "About to update node: A\n",
      "A  goes extinct!\n",
      "Not updating node A - it's extinct\n",
      "About to update node: B\n",
      "Updating node:B\n",
      "Nodes in tree at time  4 : 2\n",
      "Tree tips: ['A', 'B']\n"
     ]
    }
   ],
   "source": [
    "\n",
    "#Build the tree!\n",
    "root = PhyloNode(name=\"root\")\n",
    "A = PhyloNode(name=\"A\")\n",
    "B = PhyloNode(name=\"B\")\n",
    "root.addChild(A)\n",
    "root.addChild(B)\n",
    "\n",
    "time = 5\n",
    "#Demo our tree functions\n",
    "print(\"Root node name:\",root.Name)\n",
    "\n",
    "for t in range(1,time):\n",
    "    print(\"=\"*80)\n",
    "    print(\"Beginning simulation of time:\",t)\n",
    "    print(\"=\"*80)\n",
    "    #I want to keep the original set of nodes to update\n",
    "    #so that new species don't immediately get simulated\n",
    "    nodes_to_update = [n for n in root.getDescendants()]\n",
    "    for node in nodes_to_update:\n",
    "\n",
    "        if not node.isTip():\n",
    "            continue\n",
    "        print(\"About to update node:\",node.Name)\n",
    "        node.update()\n",
    "    print(\"Nodes in tree at time \",t,\":\",len([n.Name for n in root.getDescendants()]))\n",
    "    print(\"Tree tips:\",[n.Name for n in root.getDescendants() if n.isTip()])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Simulating Trait Evolution\n",
    "\n",
    "We can add traits to our birth-death model of tree evolution by adding a dict of BinaryTraits to each tree node, and then modifying it as necessary over time. We need to set up the __init__ method of the PhyloNode object to accept trait information. We can then pass on the traits of ancestral organisms to their descendants when nodes speciate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "from copy import deepcopy\n",
    "from random import random,choice\n",
    "\n",
    "class PhyloNode():\n",
    "    \"\"\"A node on a phylogenetic tree\"\"\"\n",
    "\n",
    "    def __init__(self,children=None,parent=None,\\\n",
    "      name = None,binary_traits=None):\n",
    "      \"\"\"Initiate a node on a phylogenetic tree\n",
    "      children -- a list of PhyloNode objects\n",
    "        descending immediately from this node\n",
    "      name -- a string for the name of the node\n",
    "      parent -- a single PhyloNode object for the parent\n",
    "      binary_traits -- a dict of True/False traits\n",
    "      \"\"\"\n",
    "      self.Name = name\n",
    "      self.Children = []\n",
    "      if children:\n",
    "        self.Children.extend(children)\n",
    "      self.Parent = parent\n",
    "      self.BinaryTraits = deepcopy(binary_traits) or {}\n",
    "      self.Extinct = False\n",
    "\n",
    "    def isTip(self):\n",
    "      \"\"\"Return True if the node is a tip\"\"\"\n",
    "      if not self.Children: #capture None or []\n",
    "        return True\n",
    "      else:\n",
    "        return False\n",
    "\n",
    "    def isRoot(self):\n",
    "      \"\"\"Return True if the node is the root of the whole tree\"\"\"\n",
    "      if not self.Parent:\n",
    "        return True\n",
    "      else:\n",
    "        return False\n",
    "\n",
    "    def getDescendants(self):\n",
    "        \"\"\"Return a list of PhyloNodes descending from the current node\"\"\"\n",
    "\n",
    "        if self.isTip():\n",
    "            print(self.Name,\" is a tip ... returning []\")\n",
    "            return []\n",
    "\n",
    "        descendants = self.Children or []\n",
    "        for c in self.Children:\n",
    "            #The set of descendants is described\n",
    "            #by the descendants of all the nodes\n",
    "            #immediate children.\n",
    "            if not c.isTip():\n",
    "                child_descendants = c.getDescendants()\n",
    "                descendants.extend([c for c in child_descendants if c not in descendants])\n",
    "\n",
    "            #Side note: this will fail on enormous trees\n",
    "            #due to the recursion limit. Not normally a problem though.\n",
    "        return descendants\n",
    "\n",
    "    def getAncestors(self):\n",
    "        \"\"\"Return the ancestors of the given node\"\"\"\n",
    "        if self.isRoot():\n",
    "            return None\n",
    "\n",
    "        ancestors = [self.Parent]\n",
    "        parents_ancestors = self.Parent.getAncestors()\n",
    "        if parents_ancestors:\n",
    "            ancestors.extend(parents_ancestors)\n",
    "        return ancestors\n",
    "\n",
    "    def getRoot(self):\n",
    "        \"\"\"Return the root node\"\"\"\n",
    "        curr_node = self\n",
    "        while not curr_node.isRoot():\n",
    "            curr_node = curr_node.Parent\n",
    "        return curr_node\n",
    "\n",
    "    def addChild(self,child):\n",
    "        \"\"\"Attach a child node\"\"\"\n",
    "        if child not in self.Children:\n",
    "            self.Children.append(child)\n",
    "        child.Parent = self\n",
    " \n",
    "    def addParent(self,parent):\n",
    "        \"\"\"Attach a parent node\"\"\"\n",
    "        self.Parent = parent\n",
    "        parent.Children.append(self)\n",
    "    \n",
    "    def getMRCA(self,other):\n",
    "        \"\"\"Return PhyloNode for the most recent common ancestor with other\n",
    "        other -- another PhyloNode object in the same tree\n",
    "        \"\"\"\n",
    "        most_recent_common_ancestor = None\n",
    "        #Cache the ancestors of other since\n",
    "        #we'll refer to it often\n",
    "        other_ancestors = other.getAncestors()\n",
    "        self_ancestors = self.getAncestors()\n",
    "        #First check for the trivial case\n",
    "        #where one node is root (but be sure they're on the same tree)\n",
    "        if self.isRoot() and self in other_ancestors:\n",
    "            return self\n",
    "\n",
    "        if other.isRoot() and other in self_ancestors:\n",
    "            return other\n",
    "\n",
    "        for a in self_ancestors:\n",
    "            #Note these will be in order of relatedness\n",
    "            if a in other_ancestors:\n",
    "                #End the loop the first time this happens\n",
    "                #since we want the most\n",
    "                #recent common ancestor\n",
    "                most_recent_common_ancestor = a\n",
    "                break\n",
    "\n",
    "        if not most_recent_common_ancestor:\n",
    "            raise ValueError(\"No common ancestor found for \",self.Name,\" and \",other.Name,\\\n",
    "            \". Are they on the same tree?\")\n",
    "\n",
    "        return most_recent_common_ancestor\n",
    "    \n",
    "    def speciate(self):\n",
    "        \"\"\"Add two children descending from this node\"\"\"\n",
    "        if self.Children:\n",
    "            raise ValueError(\"Internal nodes can't speciate\")\n",
    "        if self.Extinct:\n",
    "            raise ValueError(\"Extinct nodes can't speciate\")\n",
    "\n",
    "        child1_name = self.Name + \"A\"\n",
    "        child1=PhyloNode(name=child1_name,binary_traits=self.BinaryTraits)\n",
    "        self.addChild(child1)\n",
    "        child2_name = self.Name + \"B\"\n",
    "        child2=PhyloNode(name=child2_name,binary_traits=self.BinaryTraits)\n",
    "        self.addChild(child2)\n",
    "\n",
    "\n",
    "    def update(self,speciation_chance=0.25,trait_change_chance=0.25,\\\n",
    "        extinction_chance=0.25):\n",
    "        \"\"\"Update the node by speciation, extinction or trait gain or loss\"\"\"\n",
    "        if not self.isTip():\n",
    "            print(\"Not updating node {name} - it's not a tip\".format(name=self.Name))\n",
    "            return None\n",
    "\n",
    "        if random() < extinction_chance:\n",
    "            print(self.Name,\" goes extinct!\")\n",
    "            self.Extinct = True\n",
    "\n",
    "        if self.Extinct:\n",
    "            print(\"Not updating node {name} - it's extinct\".format(name=self.Name))\n",
    "            return None\n",
    "        print(\"Updating node:{name}\".format(name=self.Name))\n",
    "\n",
    "\n",
    "        for trait,value in self.BinaryTraits.items():\n",
    "            if random() < trait_change_chance:\n",
    "                self.BinaryTraits[trait] = not self.BinaryTraits[trait]\n",
    "                print(self.Name,\" has a new trait value for \",\\\n",
    "                  trait,\": \",self.BinaryTraits[trait],\"!\")\n",
    "\n",
    "        if random() < speciation_chance:\n",
    "            #Speciate\n",
    "            print(\"Node {name} speciates!\".format(name=self.Name))\n",
    "            self.speciate()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's use our updated, trait-aware PhyloNode object to simulate tree evolution:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Root node name: root\n",
      "================================================================================\n",
      "Beginning simulation of time: 1\n",
      "================================================================================\n",
      "About to update node: A\n",
      "Updating node:A\n",
      "About to update node: B\n",
      "B  goes extinct!\n",
      "Not updating node B - it's extinct\n",
      "Nodes in tree at time  1 : 2\n",
      "Tree tips: ['A', 'B']\n",
      "================================================================================\n",
      "Beginning simulation of time: 2\n",
      "================================================================================\n",
      "About to update node: A\n",
      "Updating node:A\n",
      "A  has a new trait value for  Horn :  True !\n",
      "About to update node: B\n",
      "B  goes extinct!\n",
      "Not updating node B - it's extinct\n",
      "Nodes in tree at time  2 : 2\n",
      "Tree tips: ['A', 'B']\n",
      "================================================================================\n",
      "Beginning simulation of time: 3\n",
      "================================================================================\n",
      "About to update node: A\n",
      "Updating node:A\n",
      "A  has a new trait value for  Red coat :  False !\n",
      "A  has a new trait value for  Big Ears :  True !\n",
      "A  has a new trait value for  Horn :  False !\n",
      "About to update node: B\n",
      "Not updating node B - it's extinct\n",
      "Nodes in tree at time  3 : 2\n",
      "Tree tips: ['A', 'B']\n",
      "================================================================================\n",
      "Beginning simulation of time: 4\n",
      "================================================================================\n",
      "About to update node: A\n",
      "A  goes extinct!\n",
      "Not updating node A - it's extinct\n",
      "About to update node: B\n",
      "B  goes extinct!\n",
      "Not updating node B - it's extinct\n",
      "Nodes in tree at time  4 : 2\n",
      "Tree tips: ['A', 'B']\n"
     ]
    }
   ],
   "source": [
    "#Define traits\n",
    "binary_traits = {\"Blubber\":False,\"Red coat\":True,\"Big Ears\":False,\\\n",
    "  \"Horn\":False,\"Wings\":False}\n",
    "#Build the tree!\n",
    "root = PhyloNode(name=\"root\",binary_traits=binary_traits)\n",
    "A = PhyloNode(name=\"A\",binary_traits=binary_traits)\n",
    "B = PhyloNode(name=\"B\",binary_traits=binary_traits)\n",
    "root.addChild(A)\n",
    "root.addChild(B)\n",
    "\n",
    "time = 5\n",
    "#Demo our tree functions\n",
    "print(\"Root node name:\",root.Name)\n",
    "\n",
    "for t in range(1,time):\n",
    "    print(\"=\"*80)\n",
    "    print(\"Beginning simulation of time:\",t)\n",
    "    print(\"=\"*80)\n",
    "    #I want to keep the original set of nodes to update\n",
    "    #so that new species don't immediately get simulated\n",
    "    nodes_to_update = [n for n in root.getDescendants()]\n",
    "    for node in nodes_to_update:\n",
    "\n",
    "        if not node.isTip():\n",
    "            continue\n",
    "        print(\"About to update node:\",node.Name)\n",
    "        node.update()\n",
    "    print(\"Nodes in tree at time \",t,\":\",len([n.Name for n in root.getDescendants()]))\n",
    "    print(\"Tree tips:\",[n.Name for n in root.getDescendants() if n.isTip()])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercises\n",
    "\n",
    "**Exercise 1**. Quantify the range of likely species diversity for given speciation/extinction rates\n",
    "\n",
    "Rerun the simulation many times. Each time, count up the number of surviving tree tips by time 5. You'll notice it's not always the same number! The number of surviving species (tips) is a measure of the range of species diversity we might observe in a particular taxon with the speciation and extinction rates we specified up above. What are the minimum,maximum,mean, and standard deviation for species diversity with the parameters given up above?  Now try changing the speciation or extinction rate. How does that change species diversity?\n",
    "\n",
    "*Hint*: it will be less laborious to do this is you make a for loop to run the simulation many times and count up results.\n",
    "\n",
    "*Hint 2*: if you like, you can also plot the results in a histogram. See for example how to do so in matplotlib (a graphing library that comes with Anaconda python): [histogram example](https://matplotlib.org/3.1.1/gallery/statistics/hist.html)\n",
    "\n",
    "**Exercise 2** Add continuous traits. Modify our birth-death simulation to allow for continuous traits. You will likely want to define a new class property called self.ContinuousTraits that holds a dict with keys equal to the names of continuous traits, and their values represented by floats. Trait changes will be represented by modifying the trait value. You will have to make some subjective choices, but that is OK. However, do be careful to record which choices you make and consider what assumptions, if any, they encode about the process of evolution for your traits. Examples include whether you modify the value of your trait by addition vs. multiplication,  and whether you want to bound your traits to values >=0. \n",
    "\n",
    "**Exercise 3**. Allow for unequal rates of evolution. Our birth-death simulation as written currently assumes that all traits evolve at equal rates. Modify it to allow for each trait to evolve at a different rate (e.g. to have a different chance of changing on each branch). *Hint*: you may wish to store the rate of evolution of each trait in a dict."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Project Ideas"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### How do the rates of speciation, extinction, and trait change affect the occurence of homology vs. convergent evolution?\n",
    "\n",
    "*Rationale*: When species share a trait, a fundamental evolutionary question is whether that trait arose from a common evolutionary origin or if the trait evolved separately in the lineages leading to the two species. Presumably, rates of speciation vs. rates of trait change should affect how often shared traits are homologous vs. convergently evolved. Here are some example predictions one might make:\n",
    "\n",
    "* if speciation is fast and trait change is slow, many species will tend to share common trait values that they inherit from their common ancestor (homologous traits)\n",
    "* if trait change is extremely fast and speciation is slow, even sister taxa may have essentially random trait values and nearly all shared traits will be convergently evolved.\n",
    "\n",
    "*Formal problem* For a tree with N nodes and a single binary trait, use simulation to find the probability that two random nodes on the tree share the trait by homology. Also find the probability that two random nodes on the tree share the trait by convergent evolution. How does the ratio of these probabilities change with varying sizes of tree, and rates of speciation, extinction, and trait change?\n",
    "\n",
    "*Relevant methods*. Testing the probability above is a form of Monte Carlo method. Reviewing or learning for the first time about Monte Carlo methods (e.g. in the statistics section of this text) may provide useful background."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}