TMRCA analysis of BEAST EBOV trees from Dudas et al.

ebov_tmrca.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Loading required package: ade4\n",
      "\n",
      "Attaching package: ‘adephylo’\n",
      "\n",
      "The following object is masked from ‘package:ade4’:\n",
      "\n",
      "    orthogram\n",
      "\n",
      "\n",
      "Attaching package: ‘lubridate’\n",
      "\n",
      "The following object is masked from ‘package:base’:\n",
      "\n",
      "    date\n",
      "\n"
     ]
    }
   ],
   "source": [
    "library(ape)\n",
    "library(adephylo)\n",
    "library(lubridate)\n",
    "library(ggplot2)\n",
    "library(coda)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "trees <- read.nexus(\"Makona_1610_cds_ig.1K.trees\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "tmrca <- function(tr){\n",
    "  decimal_date(ymd(\"2015-10-24\"))-unname(distRoot(tr,\"'EBOV|EM_FORE_2015_2878||GIN|Forecariah|2015-10-24'\"))\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "trees.tmrca <- unlist(lapply(trees,tmrca))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {},
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAC31BMVEUAAAABAQEDAwMEBAQF\nBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8QEBARERESEhITExMUFBQVFRUWFhYX\nFxcYGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIjIyMkJCQlJSUnJycoKCgpKSkq\nKiorKyssLCwtLS0uLi4vLy8wMDAxMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6Ojo8PDw+\nPj4/Pz9AQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhJSUlKSkpLS0tNTU1OTk5PT09QUFBR\nUVFSUlJTU1NUVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1eXl5fX19gYGBhYWFiYmJj\nY2NkZGRlZWVmZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29wcHBxcXFycnJzc3N0dHR1\ndXV2dnZ3d3d4eHh5eXl6enp7e3t8fHx9fX1+fn5/f3+AgICBgYGCgoKDg4OEhISFhYWGhoaH\nh4eIiIiJiYmKioqLi4uMjIyNjY2Ojo6Pj4+QkJCRkZGTk5OUlJSVlZWWlpaXl5eYmJiZmZma\nmpqbm5ucnJydnZ2enp6fn5+goKChoaGioqKkpKSmpqaoqKipqamrq6usrKytra2urq6vr6+w\nsLCxsbGysrKzs7O0tLS1tbW2tra3t7e4uLi5ubm6urq7u7u8vLy9vb2+vr6/v7/AwMDBwcHC\nwsLDw8PExMTFxcXGxsbHx8fIyMjJycnKysrLy8vMzMzNzc3Ozs7Pz8/R0dHS0tLT09PU1NTV\n1dXW1tbX19fY2NjZ2dna2trb29vc3Nzd3d3e3t7f39/g4ODh4eHi4uLj4+Pk5OTl5eXm5ubn\n5+fo6Ojp6enq6urr6+vs7Ozt7e3u7u7v7+/w8PDx8fHy8vLz8/P09PT19fX29vb39/f4+Pj5\n+fn6+vr7+/v8/Pz9/f3+/v7///9DuzCnAAAACXBIWXMAABJ0AAASdAHeZh94AAAgAElEQVR4\nnO3d/YNc1WGf8QHS2m7axqUkMXWaNk6bFLtxq9ROm/Qlbu20Y2FgUYQXVpYqsNM4rQIF2zhW\nHCUE2yGqX2JLqpMISOIaioFQyxQjDEbI1BbIMrWwsUsCsuTV2672bXbPH9CZHTH73ZnZ7557\nde7emdnn+WEk7c6ec+6Z89HO7C6iEojonKuUvQCiQQhIRAkCElGCgESUICARJQhIRAkCElGC\nYiDVrhqt357++Mi6Dz4fwtztmzbsrBW9LqK+KgLS1O5qA9KHNu8/+IGRsXDn+kf3XbOj9d4z\no1Edn52Mu2POTk0UOvzY7Fih40+cKnT4ydnjhY4/U+joo7NThQ5/LmcnA6R73lZtQBqrPlFX\nc8VDteH7Q9h75eRL7x4/GtWxMBV3x5ydnCh0+NPhdKHjT5wsdPipcKzQ8WcLHf1YmC50/BOT\n+T82A6QTz+1tQPrullP1p3Vv/+zh6pG6nuoBIKUMSKYBgRTCN6svfQZ7pHpof7Xx+mhob/1m\n9vl634v7FHg8TOf//BnR6WKfOY6H8ULHnzxd6PDTodindrOFjn48FPvU8VT+Z47fzwWpdtfa\nj4WHLmv8duS++s3omno73EcSDXQLX3TLAOm5dw/dMxfOfkZ6uH4zdlO9BybjCrORd8zXdK3Q\n4WfCTKHj16YLHX42FDr85Fyxwxd9dvIPP5ED0oHLtjV+OVytPy+c4DVS2niNZBqs10jTV39q\nrvGH2vCeEB4falEEUoqAZBosSPurDz5R78Vwx8jBQ5t3td4JpBQByTRYkO6uzndvmNu9ccOu\nhRdZQEoRkEwDA8kGpBQByQQkCUg2IJmAJAHJBiQTkCQg2YBkApIEJBuQTECSgGQDkglIEpBs\nQDIBSQKSDUgmIElAsgHJBCQJSDYgmYAkAckGJBOQJCDZgGQCkgQkG5BMQJKAZAOSCUgSkGxA\nMgFJApINSCYgSUCyAckEJAlINiCZgCQByQYkE5AkINmAZAKSBCQbkExAkoBkA5IJSBKQbEAy\nAUkCkg1IJiBJQLIByQQkCUg2IJmAJAHJBqRure+okGmAJAHJBSQXkCQguYDkApIEJBeQXECS\ngOQCkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA\n5AKSC0gSkFxAcgFJApILSC4gSUByAckFJAlILiC5gCQByQUkF5AkILmA5AKSBCQXkFxAkoDk\nApILSBKQXEByAUkCkgtILiBJQHIByQUkCUguILmAJAHJBSQXkCQguYDkApIEJBeQXECSgOQC\nkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA5AKS\nC0gSkFxAcgFJApILSC4gSUByAckFJAlILiC5gCQByQUkF5AkILmA5AKSBCQXkFxAkoDkApIL\nSBKQXEByAUkCkgtILiBJQHIByQUkCUguILmAJAHJBSQXkCQguYDkApIEJBeQXECSgOQCkgtI\nEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcvQHpzImoToaZ\nuDvmbGyq0OHPxF5nzqbGCh1+JpwsdPzZYobthFTINGPT+T82HaTJiagmw2zcHXM2VSt0+Okw\nXej4talCh5+NfZhyNlfMsJ2QCpnmHM7OeDpIPLVLEU/turWqntoBKUVA6haQOgOSDUjdAlJn\nQLIBqVtA6gxINiB1C0idAckGpG4BqTMg2YDULSB1BiQbkLoFpM6AZANSt4DUGZBsQOoWkDoD\nkg1I3QJSZ0CyAalbQOoMSDYgdQtInQHJBqRuAakzINmA1C0gdQYkG5C6BaTOgGQDUreA1BmQ\nbEDqFpA6A5INSF3UdCnBPJ0BSQKSC0guIElAcgHJBSQJSC4guYAkAckFJBeQJCC5gOQCkgQk\nV9mQYkQAKSog2YAEpLgJgWQDEpDiJgSSDUhAipsQSDYgASluQiDZgASkuAmBZAMSkOImBJIN\nSECKmxBINiABKW5CINmABKS4CYFkAxKQ4iYEkg1IQIqbEEg2IAEpbkIg2YAEpLgJgWQDEpDi\nJgSSDUhAipsQSDYgASluQiDZgASkuAmBZAMSkOImBJINSECKmxBINiABKW5CINmABKS4CYFk\nAxKQ4iYEkg1IQIqbEEg2IAEpbkIg2YAEpLgJgWQDEpDiJgSSDUhAipsQSDYgASluQiDZgASk\nuAmBZAMSkOImBJINSECKmxBINiABKW5CINmABKS4CYFkAxKQ4iYEkg1IQIqbEEg2IAEpbkIg\n2YAEpLgJgWQDEpDiJgSSDUhAipsQSDYgASluQiDZgASkuAmBZAMSkOImBJINSECKmxBINiAB\nKW5CINmABKS4CYFkAxKQ4iYEkg1IQIqbEEg2IAEpbkIg2YAEpLgJgWQDEpDiJgSSDUhAipsQ\nSDYgASluQiDZgASkuAmBZAMSkOImBJINSECKmxBINiABKW5CINmABKS4CYFkAxKQ4iYEkg1I\nQIqbEEg2IAEpbkIg2YC0eiCd+ujV7/jEVAhzt2/asLMGpKQBadVAmrtpy9f2v/MjIdy5/tF9\n1+wAUtKAtGogfb16JIRn1n6/Nnx/CHuvnARSyoC0aiD9r+HGp6VLHz/cADVePQCklAFp1UB6\nfO2ZEI5W799fbbw+Gtpbv5naU+/QqahOh5m4O+bszHShw0+EiULHnz5T6PAz4XSh488t8/5O\nETH36VIhqx/Pf3ZOZoc0Nrz91NGt1Xseuqzxp5H76jeja+rtWObjiEKniJj7dGmlF75cC190\ni/+q3dObqpfuvupLZz8jPVy/mfjjevvHohoPtbg75mxyptjhw2Sh488UO3wtjBc6/twy7+8U\nEXOfLhWy+on8Z+d0Dkhh7ujU1FufOVytPy+c4DVS2niNtGpeI5343e+E8MWNs7XhPfUXTEMT\nQEoZkFYNpLn/fP2TD6+/P4Q7Rg4e2ryr9XYgpQhIqwZSOLJ16N33N0Tt3rhhFz/ZkDYgrR5I\nSwSkFAEJSHETAskGJCDFTQgkG5CAFDchkGxAAlLchECyAQlIcRMCyQYkIMVNCCQbkIAUNyGQ\nbEACUtyEQLIBCUhxEwLJBiQgxU0IJBuQgBQ3IZBsQAJS3IRAsgEJSHETAskGJCDFTQgkG5CA\nFDchkGxAAlLchECyAQlIcRMCyQYkIMVNCCQbkIAUNyGQbEACUtyEQLIBCUhxEwLJBiQgxU0I\nJBuQgBQ3IZBsQAJS3IRAsgEJSHETAsm26iDlrJDVA0kCkgtILiBJQHIByQUkCUguILmAJAHJ\nBSQXkCQguYDkApIEJBeQXECSgOQCkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckF\nJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxAcgFJApILSC4gSUByAckFJAlILiC5gCQByQUk\nF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEByAUkCkmtgIBVCC0gSkFxAcgFJApILSC4gSUBy\nAckFJAlILiC5gCQByQUkF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEByAUkCkgtILiBJQHIB\nyQUkCUguILmAJAHJBSQXkCQguYDkApIEJBeQXECSgOQCkgtIEpBcQHIBSQKSC0guIElAcgHJ\nBSQJSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxAcgFJApILSC4gSUByAckF\nJAlILiC5gCQByQUkF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEByAUkCkgtILiBJQHIBydUb\nkCbGohoPtbg75mxiptDhJ8NkoePPRG5jzmphvNDx55Z5f4GQEqz+HM7OaSBlC0g2IPHULi6e\n2tl4agekuIBkAxKQ4gKSDUhAigtINiABKS4g2YAEpLiAZAMSkOICkg1IQIoLSDYgASkuINmA\nBKS4gGQDEpDiApINSECKC0g2IAEpLiDZgASkuIBkAxKQ4gKSDUhAigtINiABKS4g2YAEpLiA\nZAMSkOICkg1IQIoLSDYgASkuINmABKS4gGQDEpDiApINSECKC0g2IAEpLiDZgASkuIBkAxKQ\n4gKSDUhAigtINiABKS4g2YAEpLiAZAMSkOICkg1IQIoLSDYgASkuINmABKS4gGQDEpDiApIN\nSECKC0g2IAEpLiDZgASkuIBkAxKQ4gKSDUhAigtINiABKS4g2YAEpLiAZAMSkOICkg1IQIoL\nSDYgASkuINmABKS4gGQDEpDiApINSECKC0g2IAEpLiA1W+IgAwlIcQGpGZDaAlK2gNQMSG0B\nKVtAagaktoCULSA1A1JbQMoWkJoBqS0gZQtIzXIebSABqRmQmgGpLSBlC0jNgNQWkLIFpGZA\nagtI2QJSMyC1BaRsAakZkNoCUraA1AxIbQEpW0BqBqS2gJQtIDUDUltAyhaQmgGpLSBlC0jN\ngNQWkLIFpGZAagtI2QJSMyC1BaRsAalZgSJyluCigCQByQUkF5AkILmA5AKSBCQXkFxAkoDk\nApILSBKQXEByAUkCkgtILiBJQHIByQUkCUguILmAJAHJBSQXkCQguYDkApIEJBeQXECSgOQC\nkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckFJFcBkNYdbP764LVAWhSQmpXNprME\nF5UaUp1E5e759x65/hVdyHzx14fe/3wIc7dv2rCzBqSkASlvCS4qNaStlYXe1OnoS5c/8NX3\nvHM23Ln+0X3X7ABS0oCUtwQXlRrSl7dvr/yn7fP9wfOdkLb8WQh/cePzteH7Q9h75SSQUgak\nvCW4qAJeI/3Lr3Z5StfsxeoL878erh6p66keAFLKgJS3BBe1sl+1e6r65S3rf+M7YX+18fpo\naG/9ZuKP6+0fi2o81OLumLPJmWKHD5OFjj9T7PC1MJ5moLLZdJbgoibyn53TS0A6/o4fu7BZ\nB6S91esefXrb1WMPXdb408h99ZvRNfV2dNyTBrey2XRW7n4sfNFtMaSNlTduvm6+jg95rPpU\nCGeuePDsZ6SH6zdTe+odOhXV6TATd8ecnZkudPiJMFHo+NNnCh1+JpxOM1DZbDpLcFHj+c/O\nySUgXbRubil7h6qn6rfv+pPD1frzwgleI6WN10h5S3BRBbxGevmuJT+JnRl6MoSxyx+pDe8J\n4fGhCSClDEh5S3BRBUD6xXctCSns3PjYwfdfNxXuGDl4aPMCOCClCEh5S3BRBUD6+o9+fHop\nSLVPb1p/S/0D53Zv3LCLn2xIG5DyluCiCoC09g2VV1zS+FrcmqU/M7UFpBQBKW8JLqoASG9u\nBaRFAalZ2Ww6S3BR/GcUEpBcQHIBSQKSC0iuAiBd0gpIiwJSs7LZdJbgoor4YkOjN/9E5V9d\nD6RFAalZ2Ww6S3BRhT21m7v7h/cDaVFAalY2m84SXFSBr5G2vAVIiwJSs7LZdJbgogqE9Ikf\nAtKigNSsbDadJbio4iBN/eKrgbQoIDUrm01nCS6qsG/I/ttXV24A0qKA1KxsNp0luKgCIP1M\ns3/+gSkgLQpIzcpm01mCi+IbshKQXEByFQJp9lt7/vxwLcQHpBQBKW8JLqoISA+8tvGv2v2j\nB4C0OCA1K5tNZwkuqgBI+y64+Lf+x12/9aMX8A3ZxQGpWdlsOktwUQVA+nc/Pv+eFy/mG7KL\nA1Kzstl0luCiCoD0w+9t/nr9jwBpUUBqVjabzhJcVAGQLjoL6QYgLQ5Izcpm01mCiyruqd33\nXs1Tu8UBqVnZbDpLcFHFfLHht++667cvvuAJIC0KSM3KZtNZgosq5MvflzS+/P3Tfx7tCEhJ\nAlLeElxUId+QrR1+4IH/yzdk2wNSs7LZdJbgooqANLrz8yH80bZjQFockJqVzaazBBdVAKTD\nf7fy4RC2VS7+NpAWBaRmZbPpLMFFFQDpir/x+ca/ov+VC9cDaVFAalY2m84SXFQR35A9+4+e\n3PB3gLQoIDUrm01nCS6qAEiv3Nr89X2vBNKigNSsbDadJbioIv5vFK8dm8fx0/8aSIsCUrOy\n2XSW4KIKgPTI+a/774/t2/36874ApEUBqVnZbDpLcFFFfPn77p9ofEP2VXdEOwJSkoCUtwQX\nVcg3ZKcev+MP947HOwJSkoCUtwQXxb/ZIAHJBSQXkCQguYDkApIEJBeQXECSgOQCkgtIEpBc\nQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxA\ncgFJApILSC4gSUByAckFJAlILiC5gCQByQUkF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEBy\nAUkCkgtILiBJQHIByQUkCUguILmAJAHJBSQXkCQguYDkApIEJBeQXECSgOQCkgtIEpBcQHIB\nSQKSC0guIElAcgHJBSQJSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFz5IJVt\nJKoEuwMkCUguILmAJAHJBSQXkCQguYDkApIEJBeQXECSgOQCkgtIEpBcQHIBSQKSC0guIElA\ncgHJBSQJSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA5BpgSF3Kep1AkoDkApILSBKQXEBy\nAUkCkgtILiBJQHIBydUbkCYnopoMs3F3zNlUrdDhp8N0oePXpgodfjb2YVpU2SDylvU6z+Hs\njKeDdOZEVCfDTNwdczY2VejwZ2KvM2dTY4UOPxNO5vioskHkLet1jk3n2JyzpYPEU7sU8dQu\nYVmvszee2gEpRUBKWNbrBJIEJBeQXECSgOQCkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4g\nuYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxAcgFJApILSC4gSUByAckFJAlILiC5\ngCQByQUkF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEByAUkCkgtILiBJQHIByQUkCUguILmA\nJAHJBSQXkCQguYDkApIEJBeQXECSgOQCkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAk\nAckFJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxAcgFJApILSC4gSUByAckFJAlILiC5gCQB\nyQUkF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEByAUkCkgtILiBJQHIByQUkCUguILmAJAHJ\nBSQXkCQguYDkApIEJBeQXECSgOQCkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckF\nJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxAcgFJApILSC4gSUByAckFJAlILiC5gCQByQUk\nF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEByAUkCkgtILiBJQHIByQUkCUguILmAJAHJBSQX\nkCQguYDkApIEJBeQXECSgOQCkgtIEpBcQHIBSQKSC0guIElAcgHJBSQJSC4guYAkAckFJBeQ\nJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxAcgFJApIrBlLZxz9dWXdnhSGN3nrV1R89FcLc7Zs2\n7KwBKWlASljW3VlZSHPvueGpJ6+9NYQ71z+675odQEoakBKWdXdWFtIL1WdDeOTSWm34/hD2\nXjkJpJQBKWFZd2dlIR1+70wIB9ZOHK4eqeupHgBSyoCUsKy7s+JfbJg7fuvWsL/aeH00tLd+\nM3ZTvQcm4wqzkXfM13St0OFnwkyh49emCx1+Nix/n7KPf7qy7s50/qM5kQvSzdXhE+Ghyxq/\nHbmvfjO6pt6OZT6I+qSyj3+6VnDTFr7olgXS0Wdue8f42c9ID9dvZp+v973RqI6H6bg75uz0\nZKHDj4fxQsefPF3o8NPh+LL3Kfv4pyvr7pyayrWpjb6fHdKRZxsAr3jscLX+vHCC10hp4zVS\nwrLuzsq+RtozXP9MNHHpvtrwnhAeH2o9OQRSioCUsKy7s7KQRtdt/8bTv7lpPNwxcvDQ5l2t\ntwMpRUBKWNbdWeGv2h26cWj4lr8MYW73xg27+MmGtAEpYVl3h5+1k4Dk6oRU9mEvsqy7AyQJ\nSC4guYAkAckFJBeQJCC5gOQCkgQkF5BcQJKA5AKSC0gSkFxAcgFJApILSC4gSUByAckFJAlI\nLiC5gCQByQUkF5AkILmA5AKSBCQXkFxAkoDkApILSBKQXEByAUkCkgtILiBJQHIByQUkCUgu\nILmAJAHJBSQXkCQguVYXpM787gBJApILSC4gSUByAckFJAlILiC5gCQByQUkF5AkILmA5AKS\nBCQXkFxAkoDkApILSBKQXEByAUkCkgtILiBJQHIByQUkCUguILmAJAHJtdohdaZ7ASQJSC4g\ntad7ASQJSC4gtad7ASQJSC4gtad7ASQJSC4gtad7ASQJSC4gtad7ASQJSK26nBogtad7ASQJ\nSK26nBogtad7ASQJSK26nBogtad7ASQJSK26nBogtad7ASQJSK26nBogtad7ASQJSK26nBog\ntad7ASQJSK26nBogtad7ASQJSK26nBogtad7ASQJSK26nBogtad7ASQJSK26nBogtad7ASQJ\nSK26nBogtad7ASQJSK26nBogtad7ASQJSK26nBogtad7ASQJSK3KPqN9kW4YkCQgtSr7jPZF\numFAkoDUquwz2hfphgFJAlKrss9oX6QbBiQJSK3KPqN9kW4YkCQgtSr7jPZFumFAkoDUquwz\n2hfphgFJAlKrss9oX6QbBiQJSK3KPqN9kW4YkCQgtSr7jPZFumFAkoDUquwz2hfphgFJAlKr\nss9oX6QbBiQJSK3KPqN9kW4YkCQgtSr7jPZFumFAkoDUquwz2hfphgFJAlKrss9oX6QbBiQJ\nSK3KPqN9kW4YkCQgtSr7jPZFumFAkoDUquwz2hfphgFJAlKrss9oX6QbBiQJSK3KPqN9kW4Y\nkCQgtSr7jPZFumFAkoDUquwz2hfphgFJAlKrss9oX6QbBiQJSK3KPqN9kW4YkCQgtSr7jPZF\numFAkoDUquwz2hfphgFJAlKrss9oX6QbBiQJSK3KPqN9kW4YkCQgtSr7jPZFumFAkoDUquwz\n2hfphgFJAlKrss9oX6QbBiQJSK3KPqN9kW4YkCQgtSr7jPZFumFAkoDUquwz2hfphgFJAlKr\nss9oX6Qb1huQpmbiCnORd8xXbbbQ4WdDwePX0o1V9hnti3TDavmP5lQ6SGdGozoepuPumLPT\nk4UOPx7GCx1/8nS6sco+o32Rbtipqdx7/f10kHhqlyKe2q1wumG98dQOSCkC0gqnGwYkCUit\nyj6jfZFuGJAkILUq+4z2RbphQJKA1KrsM9oX6YYBSQJSq7LPaF+kGwYkCUityj6jfZFuGJAk\nILUq+4z2RbphQJKA1KrsM9oX6YYBSQJSq7LPaF+kGwYkCUityj6jfZFuGJCkVQup7BPZp+kW\nAkkCEmVJtxBIEpAoS7qFQJKARFnSLQSSBCTKkm4hkCQgUZZ0C4EkrQ5IZR+/wUl3FUgSkChL\nuqtAkoBEWdJdBZIEJMqS7iqQJCBRlnRXgSQBibKkuwokqf8h+Qe72cofuEFNdxVIEpAoS7qr\nQJKARFnSXQWSBCTKku4qkCQgUZZ0V4EkAYmypLsKJAlIlCXdVSBJQKIs6a4CSQISZUl3FUgS\nkChLuqtAkoBEWdJdBZIEJMqS7iqQJCBRlnRXgSQBibKkuwokCUiUJd1VIElAoizprgJJAhJl\nSXcVSBKQKEu6q0CSgERZ0l0FkgQkypLuKpAkIFGWdFeBJAGJsqS7CiQJSJQl3VUgSUCiLOmu\nAkkCEmVJdxVIEpAoS7qrQJKARFnSXQWSBCTKku4qkCQgUZZ0V4EkAYmypLsKJAlIlCXdVSBJ\nQKJzLNcDB6RsAWnwy/XAASlbQBr8cj1wQMoWkAa/XA8ckLIFpMEv1wMHpGwBafDL9cABKVup\nIeV8ZAs+Squ7XA8kkLIFpMEv1wMJpGwBafDL9UACKVtAGvxyPZBAyhaQBr9cDySQsgWkwS/X\nA7m6IJ37pgFp8Mv1QAIpW0Aa/HI9kEDKFpAGv1wPJJCyBaTBL9cDCaRsAWnwy/VAAilbQBr8\ncj2QQMoWkAa/XA8kkLKVBVKBj2yyoamzrGdiPiBlC0iDX9YzMR+QsgWkwS/rmZgPSNkC0uCX\n9UzMB6RsAWnwy3om5gNStoA0+GU9E/MNMKR025j0YaJeL9e5A1LENiZ9mKjXy3XugBSxjUkf\nJur1cp07IEVsY9KHiXq9XOcOSBHbmPRhol4v17kDUsQ2Jn2YqNfLde5WFFKqRbuARMmLOHfn\nAGnu9k0bdtaARANfxLk7B0h3rn903zU7gEQDX8S5yw+pNnx/CHuvnAQSDXoR5y4/pMPVI3U9\n1QNAokEv4tzlh7S/2nh9NLS3fnNyuN5nZpatyxKX/6CM1WbNbDHFrZtWU8sfu6n8kB66rHE7\ncl/9ZnRNvR3L3J9ocFv4olvez0gPv/Tnnvs+UiGtwP9orMimwrFCx58tdPRjYbrQ8U9M5v/Y\n/JAOV+sfPJHpNVIjINmAZBpMSLXhPSE8PjQBpJQByTSYkMIdIwcPbd7V+iOQUgQk04BCmtu9\nccOuTD/Z0AhINiCZBhRSW0BKEZBMQJKAZAOSCUgSkGxAMgFJApINSCYgSUCyAckEJAlINiCZ\ngCQByQYkE5AkINmAZAKSBCQbkExAkoBkA5IJSBKQbEAyAUkCkg1IJiBJQLIByQQkCUg2IJmA\nJAHJBiQTkCQg2YBkApIEJBuQTECSgGQDkglIEpBsQDIBSQKSDUgmIElAsgHJBCQJSDYgmYAk\nAckGJBOQJCDZgGQCkgQkG5BMQJKAZAOSCUgSkGxAMgFJApINSCYgSUCyAckEJOmFXZ/Lv9qI\nRos9iU/uerLQ8U+OFjr853a9UOj4Y4WO/uKuuwsdf/RU/o9NB/8cc3MAAAZwSURBVCmyk2t+\nbYVmKqS71txV9hLOpV9bc7LsJZxDU2uuLXsJywekqIBUYkCSgFRmQCo8IEUFpBIDkjR3cnyF\nZiqkqZNTZS/hXBo/OVf2Es6hvjg7KwWJaKADElGCgESUoEIgnf74yLoPPl9/cnv7pg07a423\n1K4ard+O3nrV1R891brbobeOFjH7ORe1/C/++tD7ny9vjUsXs/ozn7xm+L+dMoOU11LLD3Je\nFt7XOxUC6UOb9x/8wMhYuHP9o/uu2VF/w9Tuan0P5t5zw1NPXnvrS/ea2FztTUgxy//S5Q98\n9T3vnC11nd2LWf32zU9+7Ve3lbrMpVpi+UHPS+t9PVQRkMaqT9T/0rviodrw/SHsvXIy3PO2\namMPXqg+G8Ijl770V8nHf7U3IUUtf8ufhfAXN/bgp6SY1dfe9sUQvlI9U/Zau7TU8oOcl9b7\neqkiIH13S/1Zw9zbP3u4eiSE8eqBcOK5vY09OPzemRAOrJ1o3mvfhq/2JqSY5b9YfaHsZS5R\nzOqnLq2f1kPVXnxut9Ty9by03tdLFfbFhkeqh/ZXG3//De2t33zz7B7MHb91a/P9J4af/GZv\nQppvmeU/Vf3ylvW/8Z3y1udbbvNved/J8W03l7a85eq2fDkv8r7eqSBItbvWfiw8dFnjtyP3\nhYXH8ubq8In538xt+4PQu5CWXf7e6nWPPr3t6rHSVuhadvVh7O3V6roTZa1vmbouX8+LvK93\nKgbSc+8eumfupb85Hg4Lj+XRZ257x/y3qb+weaJ3IS2//MeqTzWeyj9Y4iKXbPnVT/6XWw49\nc9uv9OJTu6WWr+dF3tc7FQLpwGXbGpd8uHo0hIn557Lze3Ck/nI31K54rHGXT1TXrn1r9a0f\nLWL6cy1i+c0XGO/6kzKXuUQRq398qP56afaaPaWuc4mWWL6eF3lf71QEpOmrPzX/o1214T2N\nB63x+nZ+M/YM1/8mmbh0X+N9x557rv4q8sBRO1A5xSz/zNCT9SdIlz9S6kK7FrP6h4em63e4\n+t5SF9q9pZav50Xe1zsVAWl/9cEn6r0Y7hg5eGjzrsab5jdjdN32bzz9m5vGwwP/M7Te2HtF\nLX/nxscOvv+6HvxJ1pjVn7xm28Gvf3jd98pea5eWWn7rN43NX3hf71QEpLur890b5nZv3LBr\n/jsXzc04dOPQ8C1/GcLWLWHhjT1X1PJrn960/pZe/Hwatfrvbrvql7ceLnmlXVty+S/9prH8\nhff1TvysHVGCgESUICARJQhIRAkCElGCgESUICARJQhIRAkCElGCgNTjbaq81I+HkcrFZ//j\n9l+p/LUQPtd8+4+8qflf5oz9zut/8MKf+9TM/B+eP6/ye2UteVUGpB7vszfffPNI5Rfqtx+p\nQ6rM//R2mH1VE9Jb6m++aegHzvtS/W3Pvabyj6/b+JrKm+d/dmZ7pfLGMpe96gJSH/TlSvMf\nKhk5729df/YNF81D+uT8n75Q+aUQ5t54wR/Wfz+9uXJb420/97JfqDxXzmpXZ0Dqg1qQzt/w\n9+b/I4P/euFagRQu+rEQ7q00/zPyyVe+tn77XGXtJysfKWWxqzQg9UELkD5X+Ur917mf3Hj5\nIkivCeFNLz/W/MOfbpsI4cOV2/9f5Q2lLHaVBqQ+aAHS5F9/b/3Xpyr3KqT/XdkSwqveoB+x\n5oLR8E8q317ZZa7qgNQHLUAK6/9B/bndB39w4vKFLzas+yu/NB7GKuvkAw5X3hzC1sqHy1nu\nqgxIfZBA+kzlQAiv/+Vw+cKXvyvnfyyEY5Wr5QN+p7IzhK9UfraU1a7OgNQHCaSxl28N36r8\naVh4ajf3nWrl82Hub/782Tsf3X80vK6y79vf/taFlW+VtuRVF5D6IIEU3nZJuO2vngr6Gumb\nlRtC+Gc/dPYf2bux8ujXX/oe7ofKWvHqC0h9kEK6o/KNn6+GRZBmKiMhfLLS/Afyp//hK6a2\nVt55V73fr6wpbcmrLiD1QQrpxA+8+7w/WgypVvkPdUA/+bLbG7+/sXL93E+dP/8vk8/+/cqz\npa15tQWkPkghhX9/3vnHFkMKr37dXAj/529XfuY/XntJ5WfPfK1Sbb59a+XW7gNS8oDUBy2C\n9OnKvwltkN5S+Uz99uhNP/Wyi/7F9unwvspnm29/tvJPy1juqgxIRAkCElGCgESUICARJQhI\nRAkCElGCgESUICARJQhIRAn6/6D8WzGabfYYAAAAAElFTkSuQmCC",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ggplot(data.frame(TMRCA=trees.tmrca))+geom_histogram(aes(x=TMRCA),binwidth = 1./52)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "                     2.5%                       50%                     97.5% \n",
       "\"2013-09-13 19:34:44 UTC\" \"2013-12-13 07:04:00 UTC\" \"2014-01-26 16:06:57 UTC\" "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "date_decimal(quantile(trees.tmrca,probs=c(0.025,0.5,0.975)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[1] \"2013-09-30 06:03:56 UTC\" \"2014-01-31 08:59:19 UTC\""
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "date_decimal(HPDinterval(as.mcmc(trees.tmrca)))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "3.3.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}

can you try removing, say, 10% burn-in and see if the tMRCA estimates change?