Skip to content

Instantly share code, notes, and snippets.

@gtsambos

gtsambos/find-extant-lineages.ipynb

Created August 4, 2019 21:51
Show Gist options
  • Save gtsambos/c18d7e6c8666bcd10854df649da931a5 to your computer and use it in GitHub Desktop.
Save gtsambos/c18d7e6c8666bcd10854df649da931a5 to your computer and use it in GitHub Desktop.
Download ZIP
find extant lineages
Display the source blob
Display the rendered blob
Raw
find-extant-lineages.ipynb
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import msprime, tskit, heapq\n",
"from IPython.display import SVG"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"# Classes.\n",
"class Segment(object):\n",
" \"\"\"An ancestral segment mapping a given node ID to a half-open \n",
" genomic interval [left, right).\"\"\"\n",
" \n",
" def __init__(self, left, right, node):\n",
" assert left < right\n",
" self.left = left\n",
" self.right = right\n",
" self.node = node\n",
"\n",
" def __lt__(self, other):\n",
" return (self.node, self.left, self.right) < (other.node, other.left, other.right)\n",
" \n",
" def __repr__(self):\n",
" return repr((self.left, self.right, self.node))\n",
"\n",
"class Interval(object):\n",
" \"\"\"A genomic interval [left, right].\"\"\"\n",
" \n",
" def __init__(self, left, right):\n",
" assert left < right\n",
" self.left = left\n",
" self.right = right\n",
" \n",
" def __lt__(self, other):\n",
" return (self.left, self.right) < (other.left, other.right)\n",
" \n",
" def __repr__(self):\n",
" return repr((self.left, self.right))"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"def map_ancestors(Ni, Ei, S, ancestral_nodes, L):\n",
"\n",
" # Initialising the algorithm.\n",
" No = tskit.NodeTable()\n",
" Eo = tskit.EdgeTable()\n",
" A = [[] for _ in range(len(Ni))]\n",
"\n",
" # Create No and Eo.\n",
" # A: node is a sample\n",
" for u in S:\n",
" v = No.add_row(time=Ni.time[u], flags=1, population=Ni.population[u])\n",
" A[u] = [Segment(0, L, v)]\n",
"\n",
" # B and C: node not a sample\n",
" Q = []\n",
" v = -1\n",
"\n",
" for ind, e in enumerate(Ei):\n",
" u = e.parent\n",
"\n",
" for x in A[e.child]:\n",
" if x.right > e.left and e.right > x.left:\n",
" y = Segment(max(x.left, e.left), min(x.right, e.right), x.node)\n",
" heapq.heappush(Q, y)\n",
" \n",
" if ind + 1 == len(Ei) or Ei[ind + 1].parent != u:\n",
" # Process queue. \n",
" \n",
" # C. node is ancestral.\n",
" if u in ancestral_nodes:\n",
" \n",
" Q2 = []\n",
" while len(Q) > 0:\n",
" x = heapq.heappop(Q) \n",
" if v == -1:\n",
" v = No.add_row(time = Ni.time[u], population=Ni.population[u])\n",
" child = x.node\n",
" x.node = v\n",
" while len(Q) > 0 and Q[0].left <= x.right and Q[0].node == child:\n",
" r = max(x.right, Q[0].right)\n",
" x.right = r\n",
" heapq.heappop(Q)\n",
"\n",
" Eo.add_row(x.left, x.right, v, child)\n",
" Q2.append(Interval(x.left, x.right)) \n",
"\n",
" Q2.sort()\n",
" while len(Q2) > 0:\n",
" x = heapq.heappop(Q2)\n",
" while len(Q2) > 0 and Q2[0].left <= x.right:\n",
" r = max(x.right, Q2[0].right)\n",
" x.right = r\n",
" heapq.heappop(Q2)\n",
" A[u].append(Segment(x.left, x.right, v))\n",
" \n",
" \n",
" # B. node is internal to a branch\n",
" else:\n",
" while len(Q) > 0:\n",
" x = heapq.heappop(Q)\n",
" while len(Q) > 0 and Q[0].node == x.node and Q[0].left <= x.right:\n",
" r = max(x.right, Q[0].right)\n",
" x.right = r\n",
" heapq.heappop(Q) \n",
" A[u].append(x)\n",
" \n",
" # Reset node index.\n",
" v = -1\n",
"\n",
" new_tables = tskit.TableCollection(sequence_length=L)\n",
" new_tables.nodes.append_columns(time=No.time,flags=No.flags)\n",
" new_tables.edges.append_columns(left=Eo.left, right=Eo.right, \n",
" parent=Eo.parent, child=Eo.child)\n",
" new_tables.sort()\n",
" return(new_tables.tree_sequence())"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"def add_census(Ni, Ei, census_time, L):\n",
" No = tskit.NodeTable()\n",
" No.append_columns(time=Ni.time, flags=Ni.flags)\n",
" Eo = tskit.EdgeTable()\n",
" \n",
" edge_list = []\n",
" for ind, e in enumerate(Ei):\n",
" edge_list.append(e)\n",
" if ind + 1 == len(Ei) or Ei[ind + 1].parent != Ei[ind].parent:\n",
" for e in edge_list:\n",
" if Ni.time[e.parent] > census_time and Ni.time[e.child] < census_time:\n",
" v = No.add_row(time=census_time)\n",
" Eo.add_row(e.left, e.right, v, e.child)\n",
" Eo.add_row(e.left, e.right, e.parent, v)\n",
" else:\n",
" Eo.add_row(e.left, e.right, e.parent, e.child)\n",
" edge_list = []\n",
"\n",
" # Sort output. \n",
" new_tables = tskit.TableCollection(sequence_length=L)\n",
" new_tables.nodes.append_columns(time=No.time,flags=No.flags)\n",
" new_tables.edges.append_columns(left=Eo.left, right=Eo.right, \n",
" parent=Eo.parent, child=Eo.child)\n",
" new_tables.sort()\n",
" return(new_tables.tree_sequence())\n",
" "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Here's our tree sequence."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"<svg baseProfile=\"full\" height=\"200\" version=\"1.1\" width=\"200\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:ev=\"http://www.w3.org/2001/xml-events\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><defs/><g id=\"tree_0\"><g fill=\"none\" id=\"edges\" stroke=\"black\"><path d=\"M 55.0 126.57854736965575 V 30.0 H 88.75\" id=\"edge_0_7\"/><path d=\"M 40.0 170.0 V 126.57854736965575 H 55.0\" id=\"edge_0_1\"/><path d=\"M 70.0 170.0 V 126.57854736965575 H 55.0\" id=\"edge_0_2\"/><path d=\"M 122.5 103.71901009108313 V 30.0 H 88.75\" id=\"edge_0_8\"/><path d=\"M 100.0 170.0 V 103.71901009108313 H 122.5\" id=\"edge_0_0\"/><path d=\"M 145.0 131.2394055670778 V 103.71901009108313 H 122.5\" id=\"edge_0_6\"/><path d=\"M 130.0 170.0 V 131.2394055670778 H 145.0\" id=\"edge_0_3\"/><path d=\"M 160.0 170.0 V 131.2394055670778 H 145.0\" id=\"edge_0_4\"/></g><g id=\"symbols\"><g class=\"nodes\"><circle cx=\"88.75\" cy=\"30.0\" id=\"node_0_9\" r=\"3\"/><circle cx=\"55.0\" cy=\"126.57854736965575\" id=\"node_0_7\" r=\"3\"/><circle cx=\"40.0\" cy=\"170.0\" id=\"node_0_1\" r=\"3\"/><circle cx=\"70.0\" cy=\"170.0\" id=\"node_0_2\" r=\"3\"/><circle cx=\"122.5\" cy=\"103.71901009108313\" id=\"node_0_8\" r=\"3\"/><circle cx=\"100.0\" cy=\"170.0\" id=\"node_0_0\" r=\"3\"/><circle cx=\"145.0\" cy=\"131.2394055670778\" id=\"node_0_6\" r=\"3\"/><circle cx=\"130.0\" cy=\"170.0\" id=\"node_0_3\" r=\"3\"/><circle cx=\"160.0\" cy=\"170.0\" id=\"node_0_4\" r=\"3\"/></g><g class=\"mutations\" fill=\"red\"/></g><g font-size=\"14\" id=\"labels\"><g class=\"nodes\"><g text-anchor=\"start\"><text x=\"127.5\" y=\"98.71901009108313\">8</text><text x=\"150.0\" y=\"126.2394055670778\">6</text></g><g text-anchor=\"middle\"><text x=\"88.75\" y=\"25.0\">9</text><text x=\"40.0\" y=\"190.0\">1</text><text x=\"70.0\" y=\"190.0\">2</text><text x=\"100.0\" y=\"190.0\">0</text><text x=\"130.0\" y=\"190.0\">3</text><text x=\"160.0\" y=\"190.0\">4</text></g><g text-anchor=\"end\"><text x=\"50.0\" y=\"121.57854736965575\">7</text></g></g><g class=\"mutations\" font-style=\"italic\"><g text-anchor=\"start\"/><g text-anchor=\"end\"/></g></g></g></svg>"
],
"text/plain": [
"<IPython.core.display.SVG object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/svg+xml": [
"<svg baseProfile=\"full\" height=\"200\" version=\"1.1\" width=\"200\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:ev=\"http://www.w3.org/2001/xml-events\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><defs/><g id=\"tree_1\"><g fill=\"none\" id=\"edges\" stroke=\"black\"><path d=\"M 66.25 103.71901009108313 V 30.0 H 113.125\" id=\"edge_1_8\"/><path d=\"M 40.0 170.0 V 103.71901009108313 H 66.25\" id=\"edge_1_0\"/><path d=\"M 92.5 131.2394055670778 V 103.71901009108313 H 66.25\" id=\"edge_1_6\"/><path d=\"M 70.0 170.0 V 131.2394055670778 H 92.5\" id=\"edge_1_4\"/><path d=\"M 115.0 133.39274008685192 V 131.2394055670778 H 92.5\" id=\"edge_1_5\"/><path d=\"M 100.0 170.0 V 133.39274008685192 H 115.0\" id=\"edge_1_1\"/><path d=\"M 130.0 170.0 V 133.39274008685192 H 115.0\" id=\"edge_1_3\"/><path d=\"M 160.0 170.0 V 30.0 H 113.125\" id=\"edge_1_2\"/></g><g id=\"symbols\"><g class=\"nodes\"><circle cx=\"113.125\" cy=\"30.0\" id=\"node_1_9\" r=\"3\"/><circle cx=\"66.25\" cy=\"103.71901009108313\" id=\"node_1_8\" r=\"3\"/><circle cx=\"40.0\" cy=\"170.0\" id=\"node_1_0\" r=\"3\"/><circle cx=\"92.5\" cy=\"131.2394055670778\" id=\"node_1_6\" r=\"3\"/><circle cx=\"70.0\" cy=\"170.0\" id=\"node_1_4\" r=\"3\"/><circle cx=\"115.0\" cy=\"133.39274008685192\" id=\"node_1_5\" r=\"3\"/><circle cx=\"100.0\" cy=\"170.0\" id=\"node_1_1\" r=\"3\"/><circle cx=\"130.0\" cy=\"170.0\" id=\"node_1_3\" r=\"3\"/><circle cx=\"160.0\" cy=\"170.0\" id=\"node_1_2\" r=\"3\"/></g><g class=\"mutations\" fill=\"red\"/></g><g font-size=\"14\" id=\"labels\"><g class=\"nodes\"><g text-anchor=\"start\"><text x=\"97.5\" y=\"126.2394055670778\">6</text><text x=\"120.0\" y=\"128.39274008685192\">5</text></g><g text-anchor=\"middle\"><text x=\"113.125\" y=\"25.0\">9</text><text x=\"40.0\" y=\"190.0\">0</text><text x=\"70.0\" y=\"190.0\">4</text><text x=\"100.0\" y=\"190.0\">1</text><text x=\"130.0\" y=\"190.0\">3</text><text x=\"160.0\" y=\"190.0\">2</text></g><g text-anchor=\"end\"><text x=\"61.25\" y=\"98.71901009108313\">8</text></g></g><g class=\"mutations\" font-style=\"italic\"><g text-anchor=\"start\"/><g text-anchor=\"end\"/></g></g></g></svg>"
],
"text/plain": [
"<IPython.core.display.SVG object>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ts = msprime.simulate(sample_size=5, Ne=1, recombination_rate=.2, random_seed=1)\n",
"\n",
"for tree in ts.trees():\n",
" display(SVG(tree.draw(height=200)))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Say we want to take a census of each lineage at time=0.8 (between nodes 7 and 8), and 'simplify' the tree below this point so that there is one leaf for each extant lineage.\n",
"\n",
"This is equivalent to placing a new node on each edge that spans time= 0.8, and 'simplifying' the tree below it. Each of these nodes corresponds to an present-day extant lineage."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"<svg baseProfile=\"full\" height=\"200\" version=\"1.1\" width=\"200\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:ev=\"http://www.w3.org/2001/xml-events\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><defs/><g id=\"tree_0\"><g fill=\"none\" id=\"edges\" stroke=\"black\"><path d=\"M 62.5 103.71901009108313 V 30.0 H 103.75\" id=\"edge_0_8\"/><path d=\"M 40.0 115.72893996564198 V 103.71901009108313 H 62.5\" id=\"edge_0_10\"/><path d=\"M 40.0 170.0 V 115.72893996564198 H 40.0\" id=\"edge_0_0\"/><path d=\"M 85.0 115.72893996564198 V 103.71901009108313 H 62.5\" id=\"edge_0_11\"/><path d=\"M 85.0 131.2394055670778 V 115.72893996564198 H 85.0\" id=\"edge_0_6\"/><path d=\"M 70.0 170.0 V 131.2394055670778 H 85.0\" id=\"edge_0_3\"/><path d=\"M 100.0 170.0 V 131.2394055670778 H 85.0\" id=\"edge_0_4\"/><path d=\"M 145.0 115.72893996564198 V 30.0 H 103.75\" id=\"edge_0_13\"/><path d=\"M 145.0 126.57854736965575 V 115.72893996564198 H 145.0\" id=\"edge_0_7\"/><path d=\"M 130.0 170.0 V 126.57854736965575 H 145.0\" id=\"edge_0_1\"/><path d=\"M 160.0 170.0 V 126.57854736965575 H 145.0\" id=\"edge_0_2\"/></g><g id=\"symbols\"><g class=\"nodes\"><circle cx=\"103.75\" cy=\"30.0\" id=\"node_0_9\" r=\"3\"/><circle cx=\"62.5\" cy=\"103.71901009108313\" id=\"node_0_8\" r=\"3\"/><circle cx=\"40.0\" cy=\"115.72893996564198\" id=\"node_0_10\" r=\"3\"/><circle cx=\"40.0\" cy=\"170.0\" id=\"node_0_0\" r=\"3\"/><circle cx=\"85.0\" cy=\"115.72893996564198\" id=\"node_0_11\" r=\"3\"/><circle cx=\"85.0\" cy=\"131.2394055670778\" id=\"node_0_6\" r=\"3\"/><circle cx=\"70.0\" cy=\"170.0\" id=\"node_0_3\" r=\"3\"/><circle cx=\"100.0\" cy=\"170.0\" id=\"node_0_4\" r=\"3\"/><circle cx=\"145.0\" cy=\"115.72893996564198\" id=\"node_0_13\" r=\"3\"/><circle cx=\"145.0\" cy=\"126.57854736965575\" id=\"node_0_7\" r=\"3\"/><circle cx=\"130.0\" cy=\"170.0\" id=\"node_0_1\" r=\"3\"/><circle cx=\"160.0\" cy=\"170.0\" id=\"node_0_2\" r=\"3\"/></g><g class=\"mutations\" fill=\"red\"/></g><g font-size=\"14\" id=\"labels\"><g class=\"nodes\"><g text-anchor=\"start\"><text x=\"90.0\" y=\"110.72893996564198\">11</text><text x=\"150.0\" y=\"110.72893996564198\">13</text></g><g text-anchor=\"middle\"><text x=\"103.75\" y=\"25.0\">9</text><text x=\"40.0\" y=\"190.0\">0</text><text x=\"70.0\" y=\"190.0\">3</text><text x=\"100.0\" y=\"190.0\">4</text><text x=\"130.0\" y=\"190.0\">1</text><text x=\"160.0\" y=\"190.0\">2</text></g><g text-anchor=\"end\"><text x=\"57.5\" y=\"98.71901009108313\">8</text><text x=\"35.0\" y=\"110.72893996564198\">10</text><text x=\"80.0\" y=\"126.2394055670778\">6</text><text x=\"140.0\" y=\"121.57854736965575\">7</text></g></g><g class=\"mutations\" font-style=\"italic\"><g text-anchor=\"start\"/><g text-anchor=\"end\"/></g></g></g></svg>"
],
"text/plain": [
"<IPython.core.display.SVG object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/svg+xml": [
"<svg baseProfile=\"full\" height=\"200\" version=\"1.1\" width=\"200\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:ev=\"http://www.w3.org/2001/xml-events\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><defs/><g id=\"tree_1\"><g fill=\"none\" id=\"edges\" stroke=\"black\"><path d=\"M 66.25 103.71901009108313 V 30.0 H 113.125\" id=\"edge_1_8\"/><path d=\"M 40.0 115.72893996564198 V 103.71901009108313 H 66.25\" id=\"edge_1_10\"/><path d=\"M 40.0 170.0 V 115.72893996564198 H 40.0\" id=\"edge_1_0\"/><path d=\"M 92.5 115.72893996564198 V 103.71901009108313 H 66.25\" id=\"edge_1_11\"/><path d=\"M 92.5 131.2394055670778 V 115.72893996564198 H 92.5\" id=\"edge_1_6\"/><path d=\"M 70.0 170.0 V 131.2394055670778 H 92.5\" id=\"edge_1_4\"/><path d=\"M 115.0 133.39274008685192 V 131.2394055670778 H 92.5\" id=\"edge_1_5\"/><path d=\"M 100.0 170.0 V 133.39274008685192 H 115.0\" id=\"edge_1_1\"/><path d=\"M 130.0 170.0 V 133.39274008685192 H 115.0\" id=\"edge_1_3\"/><path d=\"M 160.0 115.72893996564198 V 30.0 H 113.125\" id=\"edge_1_12\"/><path d=\"M 160.0 170.0 V 115.72893996564198 H 160.0\" id=\"edge_1_2\"/></g><g id=\"symbols\"><g class=\"nodes\"><circle cx=\"113.125\" cy=\"30.0\" id=\"node_1_9\" r=\"3\"/><circle cx=\"66.25\" cy=\"103.71901009108313\" id=\"node_1_8\" r=\"3\"/><circle cx=\"40.0\" cy=\"115.72893996564198\" id=\"node_1_10\" r=\"3\"/><circle cx=\"40.0\" cy=\"170.0\" id=\"node_1_0\" r=\"3\"/><circle cx=\"92.5\" cy=\"115.72893996564198\" id=\"node_1_11\" r=\"3\"/><circle cx=\"92.5\" cy=\"131.2394055670778\" id=\"node_1_6\" r=\"3\"/><circle cx=\"70.0\" cy=\"170.0\" id=\"node_1_4\" r=\"3\"/><circle cx=\"115.0\" cy=\"133.39274008685192\" id=\"node_1_5\" r=\"3\"/><circle cx=\"100.0\" cy=\"170.0\" id=\"node_1_1\" r=\"3\"/><circle cx=\"130.0\" cy=\"170.0\" id=\"node_1_3\" r=\"3\"/><circle cx=\"160.0\" cy=\"115.72893996564198\" id=\"node_1_12\" r=\"3\"/><circle cx=\"160.0\" cy=\"170.0\" id=\"node_1_2\" r=\"3\"/></g><g class=\"mutations\" fill=\"red\"/></g><g font-size=\"14\" id=\"labels\"><g class=\"nodes\"><g text-anchor=\"start\"><text x=\"97.5\" y=\"110.72893996564198\">11</text><text x=\"120.0\" y=\"128.39274008685192\">5</text><text x=\"165.0\" y=\"110.72893996564198\">12</text></g><g text-anchor=\"middle\"><text x=\"113.125\" y=\"25.0\">9</text><text x=\"40.0\" y=\"190.0\">0</text><text x=\"70.0\" y=\"190.0\">4</text><text x=\"100.0\" y=\"190.0\">1</text><text x=\"130.0\" y=\"190.0\">3</text><text x=\"160.0\" y=\"190.0\">2</text></g><g text-anchor=\"end\"><text x=\"61.25\" y=\"98.71901009108313\">8</text><text x=\"35.0\" y=\"110.72893996564198\">10</text><text x=\"87.5\" y=\"126.2394055670778\">6</text></g></g><g class=\"mutations\" font-style=\"italic\"><g text-anchor=\"start\"/><g text-anchor=\"end\"/></g></g></g></svg>"
],
"text/plain": [
"<IPython.core.display.SVG object>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# Census event: add nodes on each lineage at time = 0.8\n",
"new_ts = add_census(ts.tables.nodes,ts.tables.edges, census_time=.8, L=1)\n",
"\n",
"for t in new_ts.trees():\n",
" display(SVG(t.draw()))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The 'simplification' can be done with `map_ancestors`, setting our census nodes as `samples` and any older nodes as `ancestors."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"samples = [n.id for n in new_ts.nodes() if n.time == 0.8]\n",
"ancestors = [n.id for n in new_ts.nodes() if n.time > 0.8]"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"<svg baseProfile=\"full\" height=\"200\" version=\"1.1\" width=\"200\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:ev=\"http://www.w3.org/2001/xml-events\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><defs/><g id=\"tree_0\"><g fill=\"none\" id=\"edges\" stroke=\"black\"><path d=\"M 82.0 115.72893996564198 V 30.0 H 109.0\" id=\"edge_0_3\"/><path d=\"M 136.0 103.71901009108313 V 30.0 H 109.0\" id=\"edge_0_4\"/><path d=\"M 118.0 115.72893996564198 V 103.71901009108313 H 136.0\" id=\"edge_0_0\"/><path d=\"M 154.0 115.72893996564198 V 103.71901009108313 H 136.0\" id=\"edge_0_1\"/></g><g id=\"symbols\"><g class=\"nodes\"><circle cx=\"46.0\" cy=\"115.72893996564198\" id=\"node_0_2\" r=\"3\"/><circle cx=\"109.0\" cy=\"30.0\" id=\"node_0_5\" r=\"3\"/><circle cx=\"82.0\" cy=\"115.72893996564198\" id=\"node_0_3\" r=\"3\"/><circle cx=\"136.0\" cy=\"103.71901009108313\" id=\"node_0_4\" r=\"3\"/><circle cx=\"118.0\" cy=\"115.72893996564198\" id=\"node_0_0\" r=\"3\"/><circle cx=\"154.0\" cy=\"115.72893996564198\" id=\"node_0_1\" r=\"3\"/></g><g class=\"mutations\" fill=\"red\"/></g><g font-size=\"14\" id=\"labels\"><g class=\"nodes\"><g text-anchor=\"start\"><text x=\"141.0\" y=\"98.71901009108313\">4</text></g><g text-anchor=\"middle\"><text x=\"46.0\" y=\"135.72893996564198\">2</text><text x=\"109.0\" y=\"25.0\">5</text><text x=\"82.0\" y=\"135.72893996564198\">3</text><text x=\"118.0\" y=\"135.72893996564198\">0</text><text x=\"154.0\" y=\"135.72893996564198\">1</text></g><g text-anchor=\"end\"/></g><g class=\"mutations\" font-style=\"italic\"><g text-anchor=\"start\"/><g text-anchor=\"end\"/></g></g></g></svg>"
],
"text/plain": [
"<IPython.core.display.SVG object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/svg+xml": [
"<svg baseProfile=\"full\" height=\"200\" version=\"1.1\" width=\"200\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:ev=\"http://www.w3.org/2001/xml-events\" xmlns:xlink=\"http://www.w3.org/1999/xlink\"><defs/><g id=\"tree_1\"><g fill=\"none\" id=\"edges\" stroke=\"black\"><path d=\"M 100.0 103.71901009108313 V 30.0 H 127.0\" id=\"edge_1_4\"/><path d=\"M 82.0 115.72893996564198 V 103.71901009108313 H 100.0\" id=\"edge_1_0\"/><path d=\"M 118.0 115.72893996564198 V 103.71901009108313 H 100.0\" id=\"edge_1_1\"/><path d=\"M 154.0 115.72893996564198 V 30.0 H 127.0\" id=\"edge_1_2\"/></g><g id=\"symbols\"><g class=\"nodes\"><circle cx=\"46.0\" cy=\"115.72893996564198\" id=\"node_1_3\" r=\"3\"/><circle cx=\"127.0\" cy=\"30.0\" id=\"node_1_5\" r=\"3\"/><circle cx=\"100.0\" cy=\"103.71901009108313\" id=\"node_1_4\" r=\"3\"/><circle cx=\"82.0\" cy=\"115.72893996564198\" id=\"node_1_0\" r=\"3\"/><circle cx=\"118.0\" cy=\"115.72893996564198\" id=\"node_1_1\" r=\"3\"/><circle cx=\"154.0\" cy=\"115.72893996564198\" id=\"node_1_2\" r=\"3\"/></g><g class=\"mutations\" fill=\"red\"/></g><g font-size=\"14\" id=\"labels\"><g class=\"nodes\"><g text-anchor=\"start\"/><g text-anchor=\"middle\"><text x=\"46.0\" y=\"135.72893996564198\">3</text><text x=\"127.0\" y=\"25.0\">5</text><text x=\"82.0\" y=\"135.72893996564198\">0</text><text x=\"118.0\" y=\"135.72893996564198\">1</text><text x=\"154.0\" y=\"135.72893996564198\">2</text></g><g text-anchor=\"end\"><text x=\"95.0\" y=\"98.71901009108313\">4</text></g></g><g class=\"mutations\" font-style=\"italic\"><g text-anchor=\"start\"/><g text-anchor=\"end\"/></g></g></g></svg>"
],
"text/plain": [
"<IPython.core.display.SVG object>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"extant_lineages = map_ancestors(new_ts.tables.nodes, new_ts.tables.edges, \n",
" samples, ancestors, L=1)\n",
"\n",
"for t in extant_lineages.trees():\n",
" display(SVG(t.draw()))"
]
}
],
"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.3"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment