do example of retaining ancestral state #59
Description
Sometimes (eg for the MK test) it'd be good to retain substitutions. To do this, we'd need to establish a root, mark it as a sample, and keep it in through simplification. I'm recording how to do this here; perhaps it should be a vignette somewhere.
ts = msprime.simulate(5,recombination_rate=2.0)
e = msprime.EdgeTable()
n = msprime.NodeTable()
ts.dump_tables(edges=e, nodes=n)
root_time = max(n.time) + 1.0
current_root = n.num_rows
n.add_row(flags=1, population=0, time=root_time)
for t in ts.trees():
for r in t.roots:
e.add_row(left=t.interval[0], right=t.interval[1], parent=current_root, child=r)
msprime.sort_tables(edges=e, nodes=n)
msprime.simplify_tables(samples=ts.samples() + [current_root], edges=e, nodes=n)
rooted_ts = msprime.load_tables(edges=e, nodes=n)
current_root = ts.num_samples
print("These should all be", current_root)
for t in rooted_ts.trees():
print(t.roots)
... now, we append to the tables, and next time we simplify, we just need to add [current_root] to the list of samples. And repeat.
This doesn't deal with mutations, of course. If we add mutations after adding the root, then all is good. If we have already got mutations on the initial tree sequence, then we'd need to add mutations to the new edges we've just added that go back to that root.