{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "

ipyrad-analysis toolkit: mrbayes

\n", "\n", "In these analyses our interest is primarily in inferring accurate branch lengths under a relaxed molecular clock model. This means that tips are forced to line up at the present (time) but that rates of substitutions are allowed to vary among branches to best explain the variation in the sequence data. \n", "\n", "There is a huge range of models that can be employed using mrbayes by employing different combinations of parameter settings, model definitions, and prior settings. The ipyrad-analysis tool here is intended to make it easy to run such jobs many times (e.g., distributed in parallel) once you have decided on your settings. In addition, we provide a number of pre-set models (e.g., clock_model=2) that may be useful for simple scenarios. \n", "\n", "Here we use simulations to demonstrate the accuracy of branch length estimation when sequences come from a single versus multiple distinct genealogies (e.g., gene tree vs species tree), and show an option to fix the topology to speed up analyses when your only interest is to estimate branch lengths. " ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# conda install ipyrad toytree mrbayes -c conda-forge -c bioconda" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import toytree\n", "import ipcoal\n", "import ipyrad.analysis as ipa" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Simulate a gene tree with 14 tips and MRCA of 1M generations" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
r0r1r2r3r4r5r6r705000001000000
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "TREE = toytree.rtree.bdtree(ntips=8, b=0.8, d=0.2, seed=123)\n", "TREE = TREE.mod.node_scale_root_height(1e6)\n", "TREE.draw(ts='o', layout='d', scalebar=True);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Simulate sequences on single gene tree and write to NEXUS\n", "When Ne is greater the gene tree is more likely to deviate from the species tree topology and branch lengths." ] }, { "cell_type": "code", "execution_count": 68, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "wrote concat locus (8 x 20000bp) to /tmp/mbtest-1.nex\n" ] } ], "source": [ "# init simulator\n", "model = ipcoal.Model(TREE, Ne=2e4, nsamples=2, recomb=0)\n", "\n", "# simulate sequence data on coalescent genealogies\n", "model.sim_loci(nloci=1, nsites=20000)\n", "\n", "# write results to database file\n", "model.write_concat_to_nexus(name=\"mbtest-1\", outdir='/tmp', diploid=True)" ] }, { "cell_type": "code", "execution_count": 69, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
r1-1r1-0r0-1r0-0r2-1r2-0r3-1r3-0r4-1r4-0r5-1r5-0r6-1r6-0r7-1r7-005077001015400
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# the simulated genealogy of haploid alleles\n", "gene = model.df.genealogy[0]\n", "\n", "# draw the genealogy\n", "toytree.tree(gene).draw(ts='o', layout='d', scalebar=True);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### View an example locus\n", "This shows the 2 haploid samples simulated for each tip in the species tree." ] }, { "cell_type": "code", "execution_count": 70, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/html": [ "
r0-0r0-1r1-0r1-1r2-0r2-1r3-0r3-1r4-0r4-1r5-0r5-1r6-0r6-1r7-0r7-1
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "model.draw_seqview(idx=0, start=0, end=50);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### (1) Infer a tree under a relaxed molecular clock model" ] }, { "cell_type": "code", "execution_count": 71, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Existing results loaded for run [3], see .trees attribute.\n", "brlenspr clock:uniform \n", "clockratepr normal(0.01,0.005) \n", "clockvarpr igr \n", "igrvarpr exp(10.0) \n", "nchains 4 \n", "ngen 1000000 \n", "nruns 2 \n", "samplefreq 1000 \n", "topologypr fixed(fixedtree) \n", "\n" ] } ], "source": [ "# init the mb object\n", "mb = ipa.mrbayes(\n", " data=\"/tmp/mbtest-1.nex\",\n", " name=\"itest-1\",\n", " workdir=\"/tmp\",\n", " clock_model=2,\n", " constraints=TREE,\n", " ngen=int(1e6),\n", " nruns=2,\n", ")\n", "\n", "# summary of priors/params\n", "print(mb.params)" ] }, { "cell_type": "code", "execution_count": 72, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "job itest-1 finished successfully\n" ] } ], "source": [ "# start the run\n", "mb.run(force=True)" ] }, { "cell_type": "code", "execution_count": 76, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
r1r0r2r3r4r5r7r6r0r1r2r3r4r5r6r705000001000000
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# load the inferred tree\n", "mbtre = toytree.tree(\"/tmp/itest-1.nex.con.tre\", 10)\n", "\n", "# scale root node to 1e6\n", "mbtre = mbtre.mod.node_scale_root_height(1e6)\n", "\n", "# draw inferred tree \n", "c, a, m = mbtre.draw(ts='o', layout='d', scalebar=True);\n", "\n", "# draw TRUE tree in orange on the same axes\n", "TREE.draw(\n", " axes=a, \n", " ts='o', layout='d', scalebar=True, \n", " edge_colors=\"darkorange\",\n", " node_sizes=0,\n", ");" ] }, { "cell_type": "code", "execution_count": 77, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
MeanVarianceLowerUpperMedianminESSavgESSPSRF
Parameter
TH9.748e-035.926e-078.215e-031.137e-029.745e-03338.719492.5041.002
TL3.681e-025.710e-063.254e-024.197e-023.658e-0233.052375.6491.019
r(A<->C)1.649e-012.324e-041.332e-011.865e-011.642e-017.771353.4931.046
r(A<->G)1.626e-011.629e-041.361e-011.890e-011.628e-01683.485693.2711.000
r(A<->T)1.948e-012.129e-041.699e-012.270e-011.922e-0120.844306.5301.007
r(C<->G)1.814e-011.886e-041.550e-012.091e-011.831e-0184.910297.1401.004
r(C<->T)1.312e-011.544e-041.069e-011.559e-011.292e-0155.086365.4851.015
r(G<->T)1.651e-011.697e-041.401e-011.912e-011.631e-0136.611320.2271.009
pi(A)2.488e-017.344e-062.438e-012.544e-012.485e-01387.327441.5241.000
pi(C)2.518e-011.120e-052.462e-012.572e-012.517e-016.667320.2411.064
pi(G)2.484e-019.192e-062.434e-012.547e-012.482e-0113.418239.9501.021
pi(T)2.510e-018.676e-062.451e-012.568e-012.506e-0142.040251.6881.017
alpha2.994e+001.636e+009.028e-015.494e+002.847e+00631.629691.3151.002
igrvar1.309e-043.692e-081.019e-065.165e-045.420e-056.53978.0781.060
clockrate1.064e-021.476e-053.318e-031.871e-029.703e-0386.09591.6221.000
\n", "
" ], "text/plain": [ " Mean Variance Lower Upper Median minESS \\\n", "Parameter \n", "TH 9.748e-03 5.926e-07 8.215e-03 1.137e-02 9.745e-03 338.719 \n", "TL 3.681e-02 5.710e-06 3.254e-02 4.197e-02 3.658e-02 33.052 \n", "r(A<->C) 1.649e-01 2.324e-04 1.332e-01 1.865e-01 1.642e-01 7.771 \n", "r(A<->G) 1.626e-01 1.629e-04 1.361e-01 1.890e-01 1.628e-01 683.485 \n", "r(A<->T) 1.948e-01 2.129e-04 1.699e-01 2.270e-01 1.922e-01 20.844 \n", "r(C<->G) 1.814e-01 1.886e-04 1.550e-01 2.091e-01 1.831e-01 84.910 \n", "r(C<->T) 1.312e-01 1.544e-04 1.069e-01 1.559e-01 1.292e-01 55.086 \n", "r(G<->T) 1.651e-01 1.697e-04 1.401e-01 1.912e-01 1.631e-01 36.611 \n", "pi(A) 2.488e-01 7.344e-06 2.438e-01 2.544e-01 2.485e-01 387.327 \n", "pi(C) 2.518e-01 1.120e-05 2.462e-01 2.572e-01 2.517e-01 6.667 \n", "pi(G) 2.484e-01 9.192e-06 2.434e-01 2.547e-01 2.482e-01 13.418 \n", "pi(T) 2.510e-01 8.676e-06 2.451e-01 2.568e-01 2.506e-01 42.040 \n", "alpha 2.994e+00 1.636e+00 9.028e-01 5.494e+00 2.847e+00 631.629 \n", "igrvar 1.309e-04 3.692e-08 1.019e-06 5.165e-04 5.420e-05 6.539 \n", "clockrate 1.064e-02 1.476e-05 3.318e-03 1.871e-02 9.703e-03 86.095 \n", "\n", " avgESS PSRF \n", "Parameter \n", "TH 492.504 1.002 \n", "TL 375.649 1.019 \n", "r(A<->C) 353.493 1.046 \n", "r(A<->G) 693.271 1.000 \n", "r(A<->T) 306.530 1.007 \n", "r(C<->G) 297.140 1.004 \n", "r(C<->T) 365.485 1.015 \n", "r(G<->T) 320.227 1.009 \n", "pi(A) 441.524 1.000 \n", "pi(C) 320.241 1.064 \n", "pi(G) 239.950 1.021 \n", "pi(T) 251.688 1.017 \n", "alpha 691.315 1.002 \n", "igrvar 78.078 1.060 \n", "clockrate 91.622 1.000 " ] }, "execution_count": 77, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# check convergence statistics\n", "mb.convergence_stats" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### (2) Concatenated sequences from a species tree\n", "Here we use concatenated sequence data from 100 loci where each represents one or more distinct genealogies. In addition, Ne is increased to 1e5, allowing for more genealogical variation. *We expect the accuracy of estimated edge lengths will decrease since we are now adequately modeling the genealogical variation when using concatenation.*" ] }, { "cell_type": "code", "execution_count": 78, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "wrote concat locus (8 x 20000bp) to /tmp/mbtest-2.nex\n" ] } ], "source": [ "# init simulator\n", "model = ipcoal.Model(TREE, Ne=1e5, nsamples=2, recomb=0)\n", "\n", "# simulate sequence data on coalescent genealogies\n", "model.sim_loci(nloci=100, nsites=200)\n", "\n", "# write results to database file\n", "model.write_concat_to_nexus(name=\"mbtest-2\", outdir='/tmp', diploid=True)" ] }, { "cell_type": "code", "execution_count": 87, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
r2-1r2-0r0-1r4-1r4-0r3-1r3-0r1-0r0-0r1-1r5-1r5-0r7-1r7-0r6-0r6-1r1-1r0-0r1-0r0-1r2-1r2-0r3-1r3-0r4-1r4-0r5-1r5-0r7-1r7-0r6-1r6-0r3-1r3-0r0-0r1-1r1-0r0-1r2-1r2-0r4-1r4-0r5-1r5-0r6-1r6-0r7-1r7-0r0-1r0-0r1-1r2-0r1-0r2-1r4-1r4-0r3-1r3-0r5-1r5-0r6-1r6-0r7-1r7-0
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# the simulated genealogies of haploid alleles\n", "genes = model.df.genealogy[:4]\n", "\n", "# draw the genealogies of the first four loci\n", "toytree.mtree(genes).draw_tree_grid(ts='o', layout='r');" ] }, { "cell_type": "code", "execution_count": 88, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "brlenspr clock:uniform \n", "clockratepr normal(0.01,0.005) \n", "clockvarpr igr \n", "igrvarpr exp(10.0) \n", "nchains 4 \n", "ngen 1000000 \n", "nruns 2 \n", "samplefreq 1000 \n", "topologypr fixed(fixedtree) \n", "\n", "job itest-2 finished successfully\n" ] } ], "source": [ "# init the mb object\n", "mb = ipa.mrbayes(\n", " data=\"/tmp/mbtest-2.nex\",\n", " workdir=\"/tmp\",\n", " name=\"itest-2\",\n", " clock_model=2,\n", " constraints=TREE,\n", " ngen=int(1e6),\n", " nruns=2,\n", ")\n", "\n", "# summary of priors/params\n", "print(mb.params)\n", "\n", "# start the run\n", "mb.run(force=True)" ] }, { "cell_type": "code", "execution_count": 91, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
r1r0r2r3r4r5r7r6r0r1r2r3r4r5r6r705000001000000
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# load the inferred tree\n", "mbtre = toytree.tree(\"/tmp/itest-2.nex.con.tre\", 10)\n", "\n", "# scale root node from unitless to 1e6\n", "mbtre = mbtre.mod.node_scale_root_height(1e6)\n", "\n", "# draw inferred tree \n", "c, a, m = mbtre.draw(ts='o', layout='d', scalebar=True);\n", "\n", "# draw true tree in orange on the same axes\n", "TREE.draw(\n", " axes=a, \n", " ts='o', layout='d', scalebar=True, \n", " edge_colors=\"darkorange\",\n", " node_sizes=0,\n", ");" ] }, { "cell_type": "code", "execution_count": 92, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
MeanVarianceLowerUpperMedianminESSavgESSPSRF
Parameter
TH1.438e-028.827e-071.262e-021.641e-021.437e-02561.836656.4181.001
TL6.213e-021.072e-055.609e-026.894e-026.199e-02528.868598.3301.000
r(A<->C)1.307e-011.363e-041.094e-011.541e-011.309e-01656.427670.7891.000
r(A<->G)1.323e-011.300e-041.110e-011.544e-011.324e-01576.671663.8361.000
r(A<->T)2.159e-012.071e-041.889e-012.456e-012.159e-01589.471617.0540.999
r(C<->G)2.303e-012.222e-041.996e-012.583e-012.304e-01552.482599.8170.999
r(C<->T)1.479e-011.360e-041.258e-011.703e-011.473e-01593.189604.3391.000
r(G<->T)1.428e-011.458e-041.210e-011.671e-011.421e-01457.278527.6561.004
pi(A)2.479e-018.765e-062.416e-012.532e-012.479e-01515.518517.6031.000
pi(C)2.490e-018.842e-062.436e-012.552e-012.489e-01360.866439.7211.000
pi(G)2.509e-018.810e-062.454e-012.570e-012.510e-01388.195427.8271.000
pi(T)2.522e-019.324e-062.465e-012.582e-012.523e-01378.545414.8680.999
alpha4.138e-027.639e-046.544e-058.931e-023.943e-02544.383647.6921.000
igrvar1.155e-043.021e-081.036e-063.720e-046.472e-05402.092483.5180.999
clockrate1.194e-021.689e-054.819e-031.972e-021.162e-0274.12384.9881.012
\n", "
" ], "text/plain": [ " Mean Variance Lower Upper Median minESS \\\n", "Parameter \n", "TH 1.438e-02 8.827e-07 1.262e-02 1.641e-02 1.437e-02 561.836 \n", "TL 6.213e-02 1.072e-05 5.609e-02 6.894e-02 6.199e-02 528.868 \n", "r(A<->C) 1.307e-01 1.363e-04 1.094e-01 1.541e-01 1.309e-01 656.427 \n", "r(A<->G) 1.323e-01 1.300e-04 1.110e-01 1.544e-01 1.324e-01 576.671 \n", "r(A<->T) 2.159e-01 2.071e-04 1.889e-01 2.456e-01 2.159e-01 589.471 \n", "r(C<->G) 2.303e-01 2.222e-04 1.996e-01 2.583e-01 2.304e-01 552.482 \n", "r(C<->T) 1.479e-01 1.360e-04 1.258e-01 1.703e-01 1.473e-01 593.189 \n", "r(G<->T) 1.428e-01 1.458e-04 1.210e-01 1.671e-01 1.421e-01 457.278 \n", "pi(A) 2.479e-01 8.765e-06 2.416e-01 2.532e-01 2.479e-01 515.518 \n", "pi(C) 2.490e-01 8.842e-06 2.436e-01 2.552e-01 2.489e-01 360.866 \n", "pi(G) 2.509e-01 8.810e-06 2.454e-01 2.570e-01 2.510e-01 388.195 \n", "pi(T) 2.522e-01 9.324e-06 2.465e-01 2.582e-01 2.523e-01 378.545 \n", "alpha 4.138e-02 7.639e-04 6.544e-05 8.931e-02 3.943e-02 544.383 \n", "igrvar 1.155e-04 3.021e-08 1.036e-06 3.720e-04 6.472e-05 402.092 \n", "clockrate 1.194e-02 1.689e-05 4.819e-03 1.972e-02 1.162e-02 74.123 \n", "\n", " avgESS PSRF \n", "Parameter \n", "TH 656.418 1.001 \n", "TL 598.330 1.000 \n", "r(A<->C) 670.789 1.000 \n", "r(A<->G) 663.836 1.000 \n", "r(A<->T) 617.054 0.999 \n", "r(C<->G) 599.817 0.999 \n", "r(C<->T) 604.339 1.000 \n", "r(G<->T) 527.656 1.004 \n", "pi(A) 517.603 1.000 \n", "pi(C) 439.721 1.000 \n", "pi(G) 427.827 1.000 \n", "pi(T) 414.868 0.999 \n", "alpha 647.692 1.000 \n", "igrvar 483.518 0.999 \n", "clockrate 84.988 1.012 " ] }, "execution_count": 92, "metadata": {}, "output_type": "execute_result" } ], "source": [ "mb.convergence_stats" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### To see the NEXUS file (data, parameters, priors):" ] }, { "cell_type": "code", "execution_count": 93, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "#NEXUS\n", "\n", "[log block]\n", "log start filename=/tmp/itest-2.nex.log replace;\n", "\n", "[data block]\n", "execute /tmp/mbtest-2.nex;\n", "\n", "[tree block]\n", "begin trees;\n", " tree fixedtree = ((r7,r6),(r5,(r4,(r3,(r2,(r1,r0))))));\n", "end;\n", "\n", "[mb block]\n", "begin mrbayes;\n", " set autoclose=yes nowarn=yes;\n", "\n", " lset nst=6 rates=gamma;\n", "\n", " prset brlenspr=clock:uniform;\n", " prset clockvarpr=igr;\n", " prset igrvarpr=exp(10.0);\n", " prset clockratepr=normal(0.01,0.005);\n", " prset topologypr=fixed(fixedtree);\n", "\n", " \n", "\n", " mcmcp ngen=1000000 nrun=2 nchains=4;\n", " mcmcp relburnin=yes burninfrac=0.25;\n", " mcmcp samplefreq=1000;\n", " mcmcp printfreq=10000 diagnfr=5000;\n", " mcmcp filename=/tmp/itest-2.nex;\n", " mcmc;\n", "\n", " sump filename=/tmp/itest-2.nex;\n", " sumt filename=/tmp/itest-2.nex contype=allcompat;\n", "end;\n", "\n", "[log block]\n", "log stop filename=/tmp/itest-2.nex.log append;\n", "\n" ] } ], "source": [ "mb.print_nexus_string()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### (3) Tree inference (not fixed topology) and plotting support values\n", "Here we will try to infer the topology from a concatenated data set (i.e., not set a constraint on the topology). I increased the ngen setting since the MCMC chain takes longer to converge when searching over topology space. " ] }, { "cell_type": "code", "execution_count": 95, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "brlenspr clock:uniform \n", "clockratepr normal(0.01,0.005) \n", "clockvarpr igr \n", "igrvarpr exp(10.0) \n", "nchains 4 \n", "ngen 2000000 \n", "nruns 2 \n", "samplefreq 1000 \n", "topologypr uniform \n", "\n", "job itest-3 finished successfully\n" ] } ], "source": [ "# init the mb object\n", "mb = ipa.mrbayes(\n", " data=\"/tmp/mbtest-2.nex\",\n", " name=\"itest-3\",\n", " workdir=\"/tmp\",\n", " clock_model=2,\n", " ngen=int(2e6),\n", " nruns=2,\n", ")\n", "\n", "# summary of priors/params\n", "print(mb.params)\n", "\n", "# start run\n", "mb.run(force=True)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# load the inferred tree\n", "mbtre = toytree.tree(\"/tmp/itest-3.nex.con.tre\", 10)\n", "\n", "# scale root node from unitless to 1e6\n", "mbtre = mbtre.mod.node_scale_root_height(1e6)\n", "\n", "# draw inferred tree \n", "c, a, m = mbtre.draw(\n", " #ts='s', \n", " layout='d',\n", " scalebar=True, \n", " node_sizes=18, \n", " node_labels=\"prob{percent}\",\n", ");" ] } ], "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.6.7" } }, "nbformat": 4, "nbformat_minor": 4 }