Accessing mutations in the genome of one individual?

[elissasoroj]

Hello!

I am trying to figure out how to look up which mutations occurred in the genome of an individual.

#for example, looking up individual 10:
postslim.tscoal.individual(10)

I get a bunch of information, but nothing that screams "use this to look up mutations". The only information seemingly shared between individuals and mutations is nodes, but it seems like the nodes assigned to individuals do not overlap the nodes assigned to mutations. What am I missing?

Apologies for the naive question and thanks for any help!
~Elissa

[petrelharp]

Hm, that's a good question - that ought to be an easy thing to do. This reason you haven't found it is because in the tree sequence it ends up being more natural to look up things by site instead of by genome - e.g., to iterate across trees.

Now, a question: do you want to find all the mutations that either genome of a given individual has inherited? Or, all the alleles? The distinction is that a mutation object is unique, but there may be more than one unique mutation at a site that's produce the same allele (IBS but not IBD). The best answer is also going to depend on what you want to do with them.

1. Edit: my first answer was wrong; here's one based on @jeromekelleher's in this issue. If you want to get the mutations then we need to iterate over the trees and ask which nodes have inherited which mutations without intervening mutations. I'll then tabulate a vector of length equal to the number of mutations with the number of copies of each mutation that the individual carries, using some methods of the Tree object:

ts = msprime.sim_mutations(msprime.sim_ancestry(10, sequence_length=40, random_seed=6), rate=0.01, random_seed=6)
num_copies = np.zeros(ts.num_mutations, dtype='int')
ind_nodes = ts.individual(0).nodes  # the list of nodes that we're counting copies in
def inherited_mutations(ts, node):
   """
   Returns the list of mutation IDs directly inherited by the specified node. A
   mutation is inherited by a node if it is ancestral  to the node in any of
   the trees along the sequence, and there are no "intervening" mutations.
   """
   inherited = []
   for tree in ts.trees():
       for site in tree.sites():
           mutation_map = {mutation.node: mutation for mutation in site.mutations}
           u = node
           while u != tskit.NULL:
               if u in mutation_map:
                   inherited.append(mutation_map[u].id)
                   break
               u = tree.parent(u)
   return inherited

for node in ind_nodes:
    inherited = inherited_mutations(ts, node=node)
    for mut in inherited:
        num_copies[mut] += 1

Then you could iterate over these mutations, for instance, like:

for j in np.where(num_copies > 0)[0]:
  mut = ts.mutation(j)
  print(f"Has {num_copies[j]} copies of mutation:")
  print(mut)

2. If you want to get the alleles then here's one way to do it, using the ts.variants() iterator.

for v in ts.variants(samples=ind_nodes):
   print(f"The individual's genotype here is {'/'.join([v.alleles[g] for g in v.genotypes])}")