Meta tags:
Headings (most frequently used words):
model, 2026, character, multi, july, discrete, dependent, ou, phytools, for, regime, in, about, blog, joint, under, the, using, on, phylogenetic, tools, comparative, biology, pages, sunday, august, wednesday, 29, tuesday, 14, monday, 13, saturday, 11, this, me, links, archive, small, experiment, but, strong, case, estimation, of, continuous, trait, evolution, stochastic, mapping, threshold, fitthresh, comparing, to, ouwie, visualizing, process, tree, more, optimum, null, ornstein, uhlenbeck,
Text of the page (most frequently used words):
the (244), model (132), this (93), our (92), #discrete (91), theta (83), character (82), alpha (70), for (65), that (65), and (61), with (55), low (54), high (53), phytools (49), tree (48), phy (47), trait (46), log (45), set (45), letters (42), sigsq (41), using (39), history (38), let (37), multi (37), levels (36), ouwie (34), fit (33), setnames (32), but (31), fitmultiou (31), continuous (31), true (30), sim_tree (30), likelihood (29), 100 (29), regime (29), fitted (28), data (28), threshold (27), 200 (26), now (25), object (24), fitmk (24), root (23), simulate (23), plot (23), matrix (22), colors (22), smp (22), stochastic (21), july (20), going (20), just (20), may (19), some (19), can (19), share (18), mk_fit (18), parameters (18), tips (18), lwd (18), not (18), 267 (18), phylogenetic (17), parameter (17), estimates (17), generating (17), hcl (17), includes (17), process (17), june (16), joint (16), levs (16), values (16), done (16), null (16), same (16), which (16), value (16), head (16), from (16), will (16), cols (16), december (15), january (15), march (15), august (15), september (15), fit_null (15), sig2 (15), has (15), solution (15), t49 (15), off (15), nodes (15), tip (15), dependent (15), init (15), april (14), november (14), more (14), class (14), ncores (14), sigma (14), simmap (14), october (13), should (13), estimated (13), prior (13), optimization (13), fitcontinuous (13), state (13), under (13), each (13), would (13), have (13), very (13), all (13), see (13), february (12), new (12), pretty (12), these (12), thinks (12), warning (12), geiger (12), multiou (12), simulation (12), use (12), ftype (12), internal (12), follows (12), max (12), merged (12), zero (12), t38 (12), rates (11), based (11), found (11), mle (11), get (11), states (11), upwards (11), col (11), labels (11), rooted (11), branch (11), length (11), because (11), fit_mou (11), fitthresh (11), t16 (11), t21 (11), t123 (11), t124 (11), loglik (10), compare (10), ou_fit (10), though (10), t19 (10), need (10), next (10), two (10), direction (10), also (10), transition (10), library (10), out (10), mean (10), map (10), results (10), ouwie_results (10), smps (10), revell (9), package (9), posted (9), here (9), optimizing (9), was (9), summary (9), iterations (9), function (9), start (9), t53 (9), pull (9), getstates (9), different (9), are (9), make (9), fixed (9), rstate (9), sum (9), probabilities (9), slightly (9), setting (9), other (9), liam (8), about (8), email (8), cool (8), null_logl (8), again (8), list (8), iter (8), parallel (8), trace (8), then (8), method (8), only (8), rep (8), par (8), legend (8), bottomleft (8), bty (8), mapped (8), load (8), nodeheights (8), across (8), maps (8), rate (8), models (8), case (7), 2026 (7), blog (7), expect (7), attr (7), even (7), closely (7), discretization (7), post (7), used (7), aic (7), number (7), might (7), sim (7), approximation (7), 283 (7), optimum (7), between (7), regimes (7), numeric (7), ms_thresh (7), liability (7), foreach (7), pinterest (6), facebook (6), blogthis (6), comments (6), 300 (6), comparative (6), best (6), factor (6), step (6), conditions (6), since (6), cex (6), lengths (6), anc (6), dimnames (6), pbtree (6), scale (6), time (6), how (6), close (6), 2094242 (6), 500 (6), 255589 (6), what (6), 2554668 (6), when (6), generate (6), changes (6), quite (6), genus_species (6), reg (6), t35 (6), frame (6), diag (6), t28 (6), t31 (6), t32 (6), t23 (6), t24 (6), mergemappedstates (6), liabilities (6), following (6), 472 (6), ouwie_data (6), small (5), posts (5), similar (5), 1883095 (5), well (5), equal (5), due (5), treating (5), nlminb (5), aicc (5), search (5), add (5), lend (5), estimation (5), trees (5), been (5), evolution (5), single (5), one (5), diffusion (5), print (5), 287 (5), seems (5), 400 (5), maxit (5), 2000 (5), run (5), first (5), already (5), false (5), course (5), fitting (5), t89 (5), t90 (5), vary (5), mapping (5), thresh_fit (5), sapply (5), 350 (5), histories (5), 2016 (4), you (4), before (4), version (4), both (4), likelihoods (4), increase (4), far (4), finally (4), once (4), opt (4), comparison (4), match (4), those (4), null_model (4), while (4), 188436 (4), free (4), t10 (4), t54 (4), simulated (4), specify (4), good (4), nrow (4), probably (4), note (4), them (4), few (4), pre (4), there (4), 7812 (4), better (4), than (4), 0693 (4), 600 (4), var (4), bit (4), starting (4), t44 (4), t57 (4), t136 (4), t177 (4), analysis (4), into (4), loading (4), thought (4), hidden (4), characters (4), original (4), infrequently (4), oum (4), finished (4), station (4), 02713072 (4), 142 (4), long (4), 028064 (4), 02921177 (4), they (4), finite (4), space (4), mkn (4), fitzjohn (4), nsim (4), mk_model (4), mkm (4), fun (4), bias (4), 7823842 (4), options (4), scipen (4), 355 (4), foo (4), ouwie_fits (4), type (4), 458785 (4), 250 (4), 03312624 (4), 241 (4), 036172 (4), 2013 (3), 2014 (3), experiment (3), page (3), tools (3), analyses (3), recommend (3), cran (3), home (3), 5433 (3), converge (3), obj (3), don (3), ready (3), takes (3), earlier (3), flat (3), convergence (3), diagnostics (3), failed (3), frequency (3), lik (3), bnd (3), bounds (3), res (3), iteration (3), maximum (3), still (3), root_state (3), code (3), actually (3), ahead (3), decided (3), reasonable (3), obtain (3), estimate (3), another (3), ornstein (3), uhlenbeck (3), boucher (3), démery (3), interested (3), biorxiv (3), weight (3), 0493 (3), know (3), much (3), 7939 (3), 0513 (3), 0689 (3), 4188483 (3), 700 (3), 06931472 (3), sensible (3), 284 (3), obtained (3), 511178 (3), 333333 (3), confirm (3), 0829 (3), 5109335 (3), 199 (3), want (3), work (3), recent (3), recently (3), readers (3), line (3), stat_theta (3), seq (3), min (3), packages (3), usually (3), way (3), begin (3), known (3), indeed (3), 8523 (3), 3246 (3), 4163 (3), 0271 (3), matches (3), ou1 (3), result (3), names (3), goes (3), great (3), rowsums (3), jointly (3), along (3), something (3), methods (3), biased (3), 0253 (3), 1673 (3), 352 (3), 354 (3), likely (3), 249 (3), 9258 (3), 159 (2), 2017 (2), 2025 (2), strong (2), lab (2), webpage (2), latest (2), similarity (2), 58921 (2), ignore (2), unlist (2), colnames (2), quick (2), too (2), 5912 (2), 5499 (2), 1496 (2), 430 (2), 8024 (2), evolutionary (2), traits (2), however (2), including (2), who (2), simulating (2), isn (2), demo (2), practice (2), packageversion (2), check (2), corresponding (2), trying (2), dependence (2), hope (2), separately (2), identifiable (2), without (2), ago (2), details (2), general (2), compared (2), anova (2), naturally (2), whether (2), were (2), 1335 (2), 1852 (2), 0707 (2), 2202 (2), 800 (2), 2003 (2), 900 (2), 2156 (2), 7829 (2), 2094 (2), converged (2), try (2), 1910 (2), 0528 (2), 2053 (2), 0527 (2), 9134 (2), fact (2), caution (2), 185 (2), 5743 (2), 0141 (2), unknown (2), taxa (2), none (2), hand (2), required (2), consider (2), future (2), morning (2), layout (2), mar (2), xlim (2), mtext (2), adj (2), transparent (2), predicted (2), distribution (2), stat (2), norm (2), singleton (2), later (2), illustration (2), look (2), fairly (2), keeping (2), mind (2), lnl (2), bic (2), ntax (2), 3007500 (2), 1198701 (2), optima (2), half (2), life (2), reporting (2), 304729 (2), arrived (2), reliable (2), fitoum (2), initializing (2), thorough (2), summarizing (2), algorithm (2), specified (2), defaulting (2), computing (2), determinant (2), inversion (2), vcv (2), 4519 (2), 1141 (2), 4520 (2), matter (2), got (2), fitou (2), give (2), simply (2), 3416784 (2), 6347880 (2), 3128115 (2), 5001694 (2), 7747887 (2), 5183583 (2), special (2), format (2), higher (2), fix (2), 6338 (2), 1838 (2), 175 (2), towards (2), infty (2), squared (2), sample (2), over (2), thesis (2), interesting (2), uses (2), review (2), comparing (2), thresholds (2), does (2), lowermost (2), through (2), complicated (2), many (2), back (2), makes (2), sense (2), standard (2), clearly (2), side (2), big (2), discretized (2), totally (2), times (2), t25 (2), t26 (2), t27 (2), t33 (2), thresh (2), underlying (2), functions (2), making (2), mvmorph (2), lot (2), days (2), away (2), seem (2), traditional (2), perhaps (2), difference (2), underestimate (2), 4625 (2), 4682 (2), 2082 (2), 473 (2), colmeans (2), average (2), 360 (2), 356 (2), 358 (2), printing (2), round (2), digits (2), rownames (2), makecluster (2), psock (2), registerdoparallel (2), optimize (2), dopar (2), stopcluster (2), parallelize (2), 142639 (2), 634035 (2), 209871 (2), 537554 (2), 871243 (2), 072081 (2), specific (2), mfrow (2), nulo (2), nuisance (2), differ (2), mathbf (2), ylim (2), 1009 (2), 9043 (2), 2719 (2), 0331 (2), 227 (2), 193 (2), cluster (2), doparallel (2), 6548297 (2), 7918578 (2), 1488759 (2), 4502061 (2), 3959087 (2), 2302126 (2), sunday (2), workflow (2), biology (2), powered, blogger, 2010, 162, 2011, 2012, 123, 2015, 132, 154, 2018, 2019, 2020, 2021, 2022, 2023, 2024, archive, umass, boston, faculty, phylogenetics, links, view, complete, profile, web, chronicles, development, install, attempting, replicate, any, site, beta, release, subscribe, atom, older, soon, come, micro, experimentation, confirms, limit, 5498946, 5498956, 1521741, 1495582, 5928578, 5912444, 1884358, cbind, throw, 0991, 0802, 0040, 0662, 5097, 6013, 3368, 0960, 1504, 6002, 5436, 1250, 2330, 1476, 5919, 5487, 1494, 1888, 5434, 331, 1883, having, effect, mode, had, realized, matched, ran, 29465, 549895, 152174, 592858, 294564, 589128, 017699, quickly, 5425203, 2690402, 4526387, 2366645, 1525741, 9399830, people, like, adapt, necessary, anyway, visualize, could, take, larger, started, colleague, prototype, consequently, write, today, demonstrating, illustrate, extent, numerical, precision, update, inadvertently, x_0, geq, github, weeks, technical, approach, saturday, preliminary, tells, basically, evidence, falls, 582, 0985, 166496e, 547, 5624, 000000e, explains, independent, easy, remarkably, remember, definitely, 2324, 6486, 9372, 0436, 2556, 315, 4681, 5480, 6211, 3112, 1711, 0816, 1765, 270, 5138, 9475, 9025, 0608, 1907, 0725, 2221, 268, 0952, 8971, 7990, 0222, 1960, 2176, 9648, 8941, 8285, 0684, 1935, 0701, 9573, 8831, 8318, 0694, 1936, 2166, 9496, 8828, 8724, 0711, 1903, 0708, 2171, 9394, 8728, 8737, 0625, 1894, 0710, 2152, 9342, 8513, 8355, 9719, 0724, 2116, 9037, 8589, 8295, 9641, 0720, 2139, 9002, 1000, 8613, 8330, 9623, 0719, 2149, 8995, 1100, 8608, 8220, 9628, 0713, 2154, 8972, 1200, 8600, 8276, 9672, 2008, 0712, 8957, 1300, 8269, 0429, 0894, 1860, 0700, 2134, 8096, 1400, 8193, 0484, 1319, 1836, 0691, 2161, 7998, 1500, 8143, 0503, 1171, 1844, 0696, 2132, 7953, 1600, 7951, 0509, 1326, 1842, 0690, 2107, 1700, 1340, 1851, 2095, 1723, 23239962, 64856597, 93721136, 04357107, 25558924, difficult, suspicious, 3575, 3680, 8215, 2403, 0115, 2066, 309, 6850, 3023, 0148, 6594, 2779, 0615, 1947, 8126, 2191, 1288, 9057, 2837, 0611, 1861, 2140, 1372, 9328, 2794, 0554, 1914, 0423, 2355, 1808, 9946, 2713, 0499, 1976, 9893, 2229, 1901, 0192, 2716, 0514, 2043, 9538, 2210, 9799, 2748, 0519, 2054, 9183, 2172, 2083, 9723, 2747, 2056, 9137, 2164, 2088, 9738, 2751, 9135, 9724, 2753, 2052, 913, 2157, 2097, 9725, 2752, farther, apart, comfortable, undoubtedly, substantially, 285, 5048, 195, 760975, 012919, 081833, 655514, 743851, 487702, 610151, hard, 0601, 2682, 0177, 409, 2704, 4702, 0182, 0846, 2622, 2005, 5735, 0150, 2569, 0941, 279, 2555, awesome, vector, 9418267, 5788125, 2039455, 7955580, 6064443, 9103976, nearly, allows, multiple, assumes, constant, continous, empirical, shall, aside, move, forward, readily, ape, writing, versions, checking, utmost, nonetheless, works, simpler, alternative, covered, monday, comment, nothing, heights, margins, second, subplot, phenogram, axis, las, equilibrium, dnorm, sqrt, polygon, border, unbranching, era, limits, 1001, chunk, objects, become, useful, contrast, three, phytools_liam, pic, twitter, com, jrmalmqpbh, https, lplhskjoo6, tweet, included, derived, visualizing, tuesday, folks, hold, whenever, unambiguous, diminish, increases, implements, 47594, 95187, 37292, 107, 3726, 75307314, 5184579, 09854245, 1132707, forget, almost, never, genuinely, knew, always, 2784, 6184, 1283, 0281, 1417, 4558, 8575, 3260, 1157, 0268, 4242, 8537, 3254, 1143, 0267, 4172, 4531, 8535, 3252, 1145, 0269, 4170, 1142, 404, 27841916, 61842864, 12828760, 02806402, identify, blogging, remind, reader, hypothesizing, hopefully, 25048, 501, 751, 150, 3165, 524044e, 510672e, 6250224, 3348563, 45480781, why, compile, put, bring, picture, updated, closer, currently, implemented, really, 18378, 41855, compute, assure, ourselves, converging, 168076, 000000, 085106, 625022, 250476, 500953, 750953, 086, 9456, 0721, 1329, 0914, 6004, 6244, 0116, 0897, 3007, 6123, 0053, 0892, 0283, 2222, 271, 0004, 0860, 0292, degree, exactly, assumptions, twice, completeness, point, generator, according, dependentn, hypothesis, explore, integrates, uncertainty, highly, implement, popular, brian, meara, jeremy, beaulieu, described, wow, beautiful, repeat, estimable, otherwise, non, normally, centers, whole, instead, unitless, scaleless, doesn, relative, positions, annoying, interpret, lower, 7568, 341537, 847661, sequence, threshstate, reuse, things, little, instance, crossing, involves, switches, forth, perfect, capturing, form, seeing, looks, running, turn, merge, contains, finely, effectively, typically, accomplished, modifying, iterating, subsequent, steps, message, elements, safely, believe, element, feel, undertake, implicit, within, 04421, 08842, 1050167, 90152, 80305, 8949833, symmetric, distinguishable, 901524, proceed, t37, t79, t80, t29, t30, t46, t47, else, 9810162, 1667724, 9394853, 5019343, 2609373, 1475927, 4670823, 2442964, 7697432, 0546414, fastbm, simplest, example, call, yet, software, straightforward, automated, show, inasmuch, crossed, fractal, moves, didn, trust, technically, enable, sorts, actual, weird, question, regarding, your, bins, test, multivariate, shape, reconstructing, ancestral, reconstruction, sort, convert, full, transitions, recommendation, play, nice, etc, answer, right, worries, couldn, sus, own, mvbm, unfortunately, responsive, inquiries, sorry, yesterday, user, contacted, subject, text, getting, style, wednesday, imo, admittedly, fine, generally, affirms, ought, cited, becomes, tends, affect, measurement, underestimated, overestimated, astonishing, convincing, 9220, 2096, 1930, 4588, 622, 5731, 1219, 2495, 4753, 6045, 5548, 490, 4077, 6367, 8417, 4181, 2777, 6418, 6047, 4362, 4288, 2951, 6848, 2549, 5098, 0866, 4541, 1942, 7713, 2389, 4734, 0244, 4641, 2102, 7708, 1732, 4698, 0304, 4678, 2093, 7854, 1688, 4622, 4686, 2081, 732, 7824, particular, downwards, 0776599, 7045410, 2024869, 7556554, 7203710, 349, 353, 342, 359, 341, summarize, error, estimating, consequence, evident, think, arisen, widely, did, 499556, 500444, 945824, expected, change, gave, accurate, particularly, regard, 632, 1410, 4069, 5262, 0369, 378, 6136, 3327, 7524, 3626, 0282, 258, 4337, 0864, 8329, 5574, 2065, 0202, 247, 0498, 1137, 9956, 5351, 2351, 0203, 245, 0731, 0581, 8858, 6185, 2743, 0387, 242, 3865, 0999, 9000, 6364, 2739, 0330, 9284, 1005, 9045, 6323, 2720, 666, 6320, 14104895, 40690078, 52619893, 03687111, predict, 06588211, 91247554, 49325269, 40862951, 208, 88237689, shows, 182, 194, 224, 206, 209, 216, 214, 192, 186, 232, 178, 213, 180, 228, 210, 239, 207, 196, created, easier, pattern, socket, host, localhost, detectcores, speed, notice, immediately, resemble, 32705, 67295, 278772, imagine, unobserved, terrific, assume, where, comes, theta_a, theta_b, conditioning, depends, total, depth, intend, engineer, simple, circumstances, words, uncertain, among, magnitude, typical, huelsenbeck, pointed, heterogeneous, tend, knowledge, studied, explicitly, talking, documented, couple, entry, princeton, pages,
Text of the page (random words):
a b a 0 03312624 0 03312624 b 0 03312624 0 03312624 log likelihood 241 9258 r thinks it has found the ml solution even though both analyses gave as parameter estimates quite close to the generating values joint estimation seems to be even more accurate particularly with regard to alpha and sigma 2 that s cool now let s consider the case of a high rate of transition for our discrete character q 0 8 k 2 high q matrix q k k dimnames list letters 1 k letters 1 k diag high q 0 diag high q rowsums high q high q a b a 0 8 0 8 b 0 8 0 8 now the discrete character is expected to change a lot in the true history high sim_tree sim history phy high q anc y0 done simulation s high sim_tree phylogenetic tree with 250 tips and 249 internal nodes tip labels t16 t21 t123 t124 t38 t49 the tree includes a mapped 2 state discrete character with states a b rooted includes branch lengths plot this tree just for fun cols setnames hcl colors n k letters 1 k plot high sim_tree cols ftype off lwd 1 direction upwards ylim c 1 10 par lend 1 legend bottomleft letters 1 k lwd 4 col hcl colors n k cex 0 7 bty n let s generate continuous trait data under the same alpha sigma 2 and mathbf theta 0 5 as before but with our new high q history simulate continuous trait high x multiou high sim_tree alpha sig2 theta a0 theta y0 head high x t16 t21 t123 t124 t38 t49 1 142639 1 634035 1 209871 1 537554 2 871243 1 072081 of course now we might expect our stochastic character histories to differ more from each other and from the true history high y as factor getstates high sim_tree type tips head high y t16 t21 t123 t124 t38 t49 a a a a a a levels a b high mk_fit fitmk phy high y model er pi fitzjohn high mk_fit object of class fitmk fitted or set value of q a b a 0 458785 0 458785 b 0 458785 0 458785 fitted or set value of pi a b 0 499556 0 500444 due to treating the root prior as a nuisance log likelihood 159 945824 optimization method used was nlminb r thinks it has found the ml solution high smps simmap high mk_fit high smps 100 phylogenetic trees with mapped discrete characters much as we did earlier let s plot a few of these just to see what we ve got par mfrow c 5 5 nulo sapply high smps 1 25 plot lwd 1 col cols ftype off direction upwards it should be pretty evident i think that even though they seem to have arisen under the same process in this case the specific details of our discrete character histories vary quite widely one from the other as well as from our generating history i expect that the consequence of this will be that each one of our stochastic character histories is likely to add error and perhaps bias to the estimating of the continuous trait multi regime ou process let s see if that s true make our ouwie data frame high ouwie_data data frame genus_species names high x reg high y x high x head high ouwie_data genus_species reg x t16 t16 a 1 142639 t21 t21 a 1 634035 t123 t123 a 1 209871 t124 t124 a 1 537554 t38 t38 a 2 871243 t49 t49 a 1 072081 once again we ll parallelize across maps using foreach foreach mc makecluster ncores type psock registerdoparallel cl mc optimize across all 100 stochastic maps high ouwie_fits foreach i 1 length high smps dopar ouwie ouwie high smps i high ouwie_data model oum simmap tree true root station false stopcluster mc let s summarize our results foo function x setnames c x theta 1 x solution 1 x loglik c the a the b alpha sigsq log l high ouwie_results t sapply high ouwie_fits foo rownames high ouwie_results 1 nrow high ouwie_results here s just the first 25 again options scipen 5 just for printing round high ouwie_results 1 25 digits 2 the a the b alpha sigsq log l 1 0 67 3 66 0 20 0 73 350 92 2 0 64 3 55 0 19 0 74 352 68 3 0 55 3 73 0 19 0 71 349 64 4 1 11 3 40 0 19 0 75 355 78 5 0 59 3 63 0 19 0 72 350 10 6 1 24 3 54 0 22 0 80 354 80 7 1 51 3 03 0 19 0 78 360 86 8 1 02 3 40 0 21 0 78 354 09 9 1 37 4 09 0 21 0 76 350 52 10 0 94 3 56 0 20 0 75 353 15 11 0 94 4 29 0 21 0 71 342 30 12 1 14 3 25 0 22 0 81 356 49 13 1 66 3 62 0 17 0 74 359 53 14 1 46 3 64 0 20 0 77 356 67 15 1 02 3 13 0 16 0 72 358 30 16 0 31 3 84 0 22 0 72 341 27 17 1 54 3 59 0 19 0 77 358 10 18 1 43 4 06 0 18 0 73 355 63 19 1 53 3 73 0 21 0 78 355 44 20 1 35 4 16 0 18 0 72 352 88 21 1 53 3 08 0 18 0 77 360 42 22 0 42 3 49 0 19 0 73 350 55 23 1 42 3 85 0 21 0 79 354 72 24 0 23 3 70 0 17 0 68 350 52 25 1 02 3 36 0 20 0 77 355 54 options scipen 0 now let s get an average across trees colmeans high ouwie_results the a the b alpha sigsq log l 1 0776599 3 7045410 0 2024869 0 7556554 352 7203710 as we predicted the estimates of theta in particular are biased towards each other but we also see that alpha is biased downwards and sigma 2 upwards finally let s compare this to phytools fitmultiou init setnames c mean high x high y levels high y 1 mean high x high y levels high y 2 log 2 max nodeheights phy var high x max nodeheights phy fitmk phy high y model er rates c theta a theta b alpha sigsq q 1 high fit_mou fitmultiou phy high x high y model er levs 100 parallel true ncores ncores root mle trace 1 maxit 2000 init init iter the a the b alpha sigsq q 1 log l 0 1 9220 3 2096 0 0693 0 1930 0 4588 622 5731 100 1 1219 4 2495 0 4753 0 6045 0 5548 490 4077 200 0 6367 4 8417 0 4181 0 2777 0 6418 473 6047 300 0 4362 4 8024 0 4288 0 2951 0 6848 473 2549 400 0 5098 5 0866 0 4541 0 1942 0 7713 472 2389 500 0 4734 5 0244 0 4641 0 2102 0 7708 472 1732 600 0 4698 5 0304 0 4678 0 2093 0 7854 472 1688 700 0 4622 5 0253 0 4686 0 2081 0 7829 472 1673 732 0 4625 5 0253 0 4682 0 2082 0 7824 472 1673 done optimizing high fit_mou object of class fitmultiou based on a discretization with k 100 levels fitted multi theta ou model parameters levels a b theta 0 4625 5 0253 alpha 0 4682 sigsq 0 2082 estimated q matrix a b a 0 7823842 0 7823842 b 0 7823842 0 7823842 log likelihood 472 1673 r thinks it has found the ml solution this is pretty astonishing seems to make quite a convincing case for joint estimation admittedly this is a very small experiment of course however it does seem to confirm the thesis that if our discrete character changes infrequently then the traditional two step process of generating stochastic character maps then fitting a fixed regime model to each map should work just fine on the other hand it also generally affirms revell 2013 which perhaps ought to be a bit better cited that as the rate of transition in our character goes up the two step process becomes biased and tends to underestimate the difference between the regimes this also may affect the measurement of alpha which is underestimated and sigma 2 which is overestimated this bias goes away when we jointly model x and y using our new finite space approximation which is very cool imo ok more on this later posted by liam revell at 3 42 pm no comments email this blogthis share to x share to facebook share to pinterest wednesday july 29 2026 stochastic character mapping under the threshold model using fitthresh in phytools unfortunately i m not usually so responsive to email inquiries these days sorry but just yesterday a phytools user contacted me with the subject line getting a simmap style tree from a threshold model result and the following email text i have a weird question regarding your threshold model functions i have a continuous trait i am making discrete bins out of to test whether multivariate rates of shape evolution vary with the levels of this now discrete trait using mvmorph fitthresh makes a lot of sense for reconstructing this character and results in ancestral values that look different more sensible than make simmap based on a reconstruction of the underlying continuous character but i need the threshold model results in some sort of format that will work with mvbm e g simmap is there a not too complicated way to convert the full q matrix into transitions between the original discrete character or do you have a better recommendation for making fitthresh play nice with the trait fitting functions in phytools mvmorph ouwie etc if there is no good answer right now no worries just couldn t sus it out on my own indeed this is something that i d thought of already that is that the discrete diffusion finite space approximation that we use in fitthresh and other new methods of phytools would enable new sorts of stochastic character mapping including of the threshold model and indeed of actual continuous traits before i go on to show how to do this i will note that this is an approximation of the corresponding stochastic character map for the threshold process only inasmuch as technically the threshold is crossed infty times by the underlying fractal bm process when it moves from one side of the threshold to the other if you didn t already know this just trust me it s true so how do we do this well so far i have not yet automated this in software but it s pretty straightforward if we go through it step by step i m going to start with the simplest example which is the threshold model with two states let s call them a and b load phytools library phytools first let s simulate some data under this process simulate tree phy pbtree n 80 scale 1 phy phylogenetic tree with 80 tips and 79 internal nodes tip labels t28 t31 t32 t23 t24 t19 rooted includes branch length s simulate liabilities liability fastbm phy a 0 5 head liability 10 t28 t31 t32 t23 t24 t19 t25 t26 0 9810162 1 1667724 0 9394853 0 5019343 0 2609373 0 1475927 0 4670823 0 2442964 t27 t33 0 7697432 1 0546414 threshold liabilities thresh function x if x 0 a else b x sapply liability thresh head x 20 t28 t31 t32 t23 t24 t19 t25 t26 t27 t33 t35 t37 t38 t79 t80 t9 t29 t30 t46 t47 b b b b b a a b a a a a a a a b a a a a now let s proceed to fit the threshold model to these data using fitthresh as follows thresh_fit fitthresh phy x thresh_fit object of class fitthresh set value of sigsq of the liability 1 0 set or estimated threshold s 0 log likelihood 33 901524 lowermost threshold is fixed totally just for fun let s fit compare a symmetric m k model to see if they re distinguishable based on these data mk_model fitmk phy x model er anova mk_model thresh_fit log l d f aic weight mk_model 36 04421 1 74 08842 0 1050167 thresh_fit 33 90152 1 69 80305 0 8949833 now to undertake stochastic mapping using our threshold model we need to pull out the implicit m k model that s hidden within our fitthresh object mkm thresh_fit mk_fit i m not going to print it because the q matrix is 200 times 200 but readers following along should feel free to go for it next because this hidden object is not quite the same as a standard fitmk object let s add the element root prior as follows mkm root prior fitzjohn now believe it or not we re totally ready for stochastic mapping when we run the following code we re probably going to see the message warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero this is just because our matrix is so big some of the elements are close to zero so we can safely ignore it i m only going to do one stochastic character map here typically we should do many of these of course this would be accomplished by modifying nsim and then iterating all subsequent steps over each stochastic map tree smp simmap mkm nsim 1 pi mle warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero once again this is a big object because it contains all 200 of the finely discretized levels of our diffusion approximation effectively a stochastic map of liabilities to turn this into a stochastic character map of our original discrete trait we need to simply merge liabilities on each side of the single threshold as follows ss colnames mkm data merged smp mergemappedstates smp ss which as numeric ss 0 a merged smp mergemappedstates merged smp ss which as numeric ss 0 b merged smp phylogenetic tree with 80 tips and 79 internal nodes tip labels t28 t31 t32 t23 t24 t19 the tree includes a mapped 2 state discrete character with states a b rooted includes branch lengths great let s plot it plot merged smp direction upwards ftype off lwd 3 colors setnames hcl colors n 2 c a b legend bottomleft c a b col hcl colors n 2 lwd 2 bty n we can see that this is clearly capturing a form of the threshold model history of our threshold character by seeing just how different it looks compared to running the same analysis with a standard m k model mk smp simmap mk_model nsim 1 plot mk smp direction upwards ftype off lwd 3 colors setnames hcl colors n 2 c a b legend bottomleft c a b col hcl colors n 2 lwd 2 bty n we see for instance that threshold crossing usually involves many switches back forth under the threshold model and none in the m k model which makes perfect sense things get only a little more complicated when i want to fit the multi state model let s see if i can run through that first i ll reuse my original tree liabilities to simulate some tip data as follows y threshstate liability setnames c 0 1 2 letters 1 3 head y t28 t31 t32 t23 t24 t19 b c b b b a next fit the multi state threshold model as follows ms_thresh fitthresh phy y sequence letters 1 3 ms_thresh object of class fitthresh set value of sigsq of the liability 1 0 set or estimated threshold s 0 7568 0 341537 log likelihood 51 847661 lowermost threshold is fixed note that the lower estimated threshold is not a separately estimable parameter otherwise the model would be non identifiable it s normally set to zero but fitthresh does something different that i ll probably fix in the future it centers the whole liability distribution on zero instead since the liabilities are unitless and scaleless this doesn t matter at all for the model fit or relative positions of the thresholds but our results are a bit more annoying to interpret pull out our hidden m k object do stochastic mapping as before mkn ms_thresh mk_fit mkn root prior fitzjohn smp simmap mkn nsim 1 warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero warning in rstate p sum p some probabilities slightly 0 setting p 0 to zero now repeat the mergemappedstates step but this time using the estimated thresholds ss mkn states merged smp 2 mergemappedstates smp ss which as numeric ss ms_thresh threshold 1 a merged smp 2 mergemappedstates merged smp 2 ss which as numeric ss ms_thresh threshold 1 as numeric ss ms_thresh threshold 2 b merged smp 2 mergemappedstates merged smp 2 ss which as numeric ss ms_thresh threshold 2 c merged smp 2 phylogenetic tree with 80 tips and 79 internal nodes tip labels t28 t31 t32 t23 t24 t19 the tree includes a mapped 3 state discrete character with states a b c rooted includes branch le...
|