From bcf78c87ee2a9e78f72895260604f400fbc0b759 Mon Sep 17 00:00:00 2001 From: peter Date: Tue, 15 Sep 2026 10:47:30 -0700 Subject: [PATCH 1/4] can reload after convert --- docs/tutorial.md | 12 +++++++++++- pyslim/methods.py | 4 ++-- tests/test_tree_sequence.py | 22 +--------------------- 3 files changed, 14 insertions(+), 24 deletions(-) diff --git a/docs/tutorial.md b/docs/tutorial.md index 5959dd7..9011905 100644 --- a/docs/tutorial.md +++ b/docs/tutorial.md @@ -481,7 +481,6 @@ This is recorded by setting the derived state of the tskit mutation to a comma-separated string of SLiM mutation IDs (or the empty string, to denote "no mutations"), and the `"derived_states"` entry of the tskit mutation's metadata to this same list. -(This is redundant, but there are good reasons for it.) So, each SLiM mutation can appear in more than one tskit mutation, and the [metadata](sec_metadata) about these SLiM mutations is stored in top-level metadata, rather than along with the tskit mutations. @@ -507,6 +506,17 @@ for x in mut.metadata["derived_states"]: See [](sec_tutorial_selected_mutations) for an example where a tskit mutation carries more than one (stacked) SLiM mutation. +You may have noticed that the same information goes both into the derived state entry +of a tskit mutation (e.g., `mut.derived_state`) and into the metadata of that +same mutation (e.g., `mut.metadata["derived_states"]`). The reason for this +redundancy is that tskit uses the `derived_state` entries (and the +`ancestral_state` entries of sites) to do VCF output and determine allelic +identity, and so it's helpful to be able to rewrite these +(as does {func}`.convert_alleles`; see [](sec_output)). +In fact, SLiM does not use the `derived_state` entry of tskit mutations +when it reads in a tree sequence (nor the `ancestral_state` entry of sites); +the definitive source is in metadata. + (sec_extracting_individuals)= diff --git a/pyslim/methods.py b/pyslim/methods.py index ff71edd..f6af881 100644 --- a/pyslim/methods.py +++ b/pyslim/methods.py @@ -521,8 +521,8 @@ def convert_alleles(ts): metadata. In SLiM's output the list of mutation IDs is recorded both in each mutation's - derived state and metadata, but SLiM uses the derived state for loading files, - so the resulting tree sequence will not be loadable by SLiM. + derived state and metadata, but SLiM only uses the metadata for loading files, + so the resulting tree sequence will still be loadable by SLiM. The main purpose of this method is for output: for instance, this code will produce a VCF file with nucleotide alleles: diff --git a/tests/test_tree_sequence.py b/tests/test_tree_sequence.py index 6a91ea1..26acd76 100644 --- a/tests/test_tree_sequence.py +++ b/tests/test_tree_sequence.py @@ -1148,20 +1148,6 @@ def test_convert_alleles(self, recipe): cts = pyslim.convert_alleles(ts) self.verify_converted_nucleotides(ts, cts) - def replace_derived_state(self, ts): - """ - Put the information from metadata back into the derived state column. - """ - t = ts.dump_tables() - t.mutations.clear() - for m in ts.mutations(): - ds = ",".join(map(str, m.metadata["derived_states"])) - t.mutations.append(m.replace(derived_state=ds)) - t.sites.clear() - for s in ts.sites(): - t.sites.append(s.replace(ancestral_state="")) - return t.tree_sequence() - @pytest.mark.parametrize( "recipe", ["recipe_nucleotides_WF.slim", "recipe_chromosomes_adds_muts.slim"], @@ -1172,13 +1158,7 @@ def test_reload_converted(self, recipe, helper_functions, tmp_path): for chrom, ts in recipe["ts"].items(): # convert alleles cts = pyslim.convert_alleles(ts) - # put the derived state column back (since SLiM uses that) - rcts = self.replace_derived_state(cts) - converted[chrom] = rcts - ts.tables.assert_equals(rcts.tables) - # given the tables are equal we don't really need to put them back through slim - # but keeping this for the future in which we don't need to replace_derived_state - # to read them back into slim + converted[chrom] = cts multichrom = "multichrom" in recipe if multichrom: slimfile = "restart_nucleotides_WF_chromosomes.slim" From da531a57cdfa5ca053b0dbde43fcc18a293756a7 Mon Sep 17 00:00:00 2001 From: peter Date: Tue, 15 Sep 2026 10:50:09 -0700 Subject: [PATCH 2/4] default sex ratio to 0.5; closes #339 --- CHANGELOG.rst | 4 ++++ pyslim/slim_metadata.py | 2 +- tests/test_annotation.py | 2 +- 3 files changed, 6 insertions(+), 2 deletions(-) diff --git a/CHANGELOG.rst b/CHANGELOG.rst index e1ef107..d91c058 100644 --- a/CHANGELOG.rst +++ b/CHANGELOG.rst @@ -47,10 +47,14 @@ https://tskit.dev/pyslim/docs/latest/previous_versions.html (`pyslim.INDIVIDUAL_FLAG_MIGRATED`) is now recorded, as `pyslim.INDIVIDUAL_MIGRATED`, in `individual.flags` (rather than `individual.metadata['flags']`). +- The default sex ratio for populations is now 0.5 instead of 0.0. + (:issue:`339`, :user:`petrelharp`) + **Bug fixes:** - `pyslim.annotate` now has a `num_chromosomes` argument. Previously it could not be easily used to annotate multichromosome simulations with more than 8 chromosomes. + (It also now has a `num_traits` argument.) - In some previous versions, converting files produced by a yet-older version of SLiM to the previously-current file version dropped some information from metadata: diff --git a/pyslim/slim_metadata.py b/pyslim/slim_metadata.py index 49fa203..ef7afe2 100644 --- a/pyslim/slim_metadata.py +++ b/pyslim/slim_metadata.py @@ -878,7 +878,7 @@ def default_slim_metadata(name, num_chromosomes=1, num_traits=1, **kwargs): "selfing_fraction": 0.0, "female_cloning_fraction": 0.0, "male_cloning_fraction": 0.0, - "sex_ratio": 0.0, + "sex_ratio": 0.5, "bounds_x0": 0.0, "bounds_x1": 1.0, "bounds_y0": 0.0, diff --git a/tests/test_annotation.py b/tests/test_annotation.py index 95eb0ab..ec8b7c5 100644 --- a/tests/test_annotation.py +++ b/tests/test_annotation.py @@ -159,7 +159,7 @@ def verify_defaults(self, ts): assert md["selfing_fraction"] == 0.0 assert md["female_cloning_fraction"] == 0.0 assert md["male_cloning_fraction"] == 0.0 - assert md["sex_ratio"] == 0.0 + assert md["sex_ratio"] == 0.5 assert md["bounds_x0"] == 0.0 assert md["bounds_x1"] == 1.0 assert md["bounds_y0"] == 0.0 From f226a0cccd6d1c89f8b4e79bbecf0c5fa55bb6a8 Mon Sep 17 00:00:00 2001 From: peter Date: Wed, 16 Sep 2026 15:05:17 -0700 Subject: [PATCH 3/4] windows cache bump --- .github/workflows/tests.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 9257372..30e0eab 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -84,7 +84,7 @@ jobs: uses: actions/cache@2c8a9bd7457de244a408f35966fab2fb45fda9c8 # v6.0.0 with: path: D:\a\pyslim\pyslim\SLiM - key: ${{runner.os}}-${{matrix.sys}}-${{matrix.env}}-v1-key + key: ${{runner.os}}-${{matrix.sys}}-${{matrix.env}}-v2-key - name: Build SLiM (Windows) if: matrix.os == 'windows-latest' && steps.cache-slim.outputs.cache-hit != 'true' From 42cb0aba750c0a81fbbba98c13363834efdcdf29 Mon Sep 17 00:00:00 2001 From: peter Date: Thu, 17 Sep 2026 06:06:15 -0700 Subject: [PATCH 4/4] switch from "derived_states" to "slim_ids" in mutation metadata --- .github/workflows/tests.yml | 2 +- docs/previous_versions.md | 4 +-- docs/traits.md | 6 ++-- docs/tutorial.md | 44 +++++++++++++-------------- docs/vignette_coalescent_diversity.md | 14 ++++----- pyslim/methods.py | 32 +++++++++---------- pyslim/slim_metadata.py | 12 ++++---- pyslim/slim_tree_sequence.py | 6 ++-- tests/test_annotation.py | 8 ++--- tests/test_tree_sequence.py | 34 ++++++++++----------- 10 files changed, 78 insertions(+), 84 deletions(-) diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 30e0eab..d43af26 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -84,7 +84,7 @@ jobs: uses: actions/cache@2c8a9bd7457de244a408f35966fab2fb45fda9c8 # v6.0.0 with: path: D:\a\pyslim\pyslim\SLiM - key: ${{runner.os}}-${{matrix.sys}}-${{matrix.env}}-v2-key + key: ${{runner.os}}-${{matrix.sys}}-${{matrix.env}}-v3-key - name: Build SLiM (Windows) if: matrix.os == 'windows-latest' && steps.cache-slim.outputs.cache-hit != 'true' diff --git a/docs/previous_versions.md b/docs/previous_versions.md index 4fb43f1..45e80a0 100644 --- a/docs/previous_versions.md +++ b/docs/previous_versions.md @@ -72,7 +72,7 @@ replaces the previous ``type`` argument to {class}`msprime.SLiMMutationModel`. to pull this information out of top-level metadata using the SLiM ID as a key. In brief, if `mut` is a mutation, then you should replace `mut.metadata["mutation_list"][j]` -with `mut_metadata[mut.metadata["derived_states"][j]]`, +with `mut_metadata[mut.metadata["slim_ids"][j]]`, where `mut_metadata` is the output of {func}`.mutation_metadata`. For instance, where before you might have done: ```python @@ -85,7 +85,7 @@ Now, you would do: ```{code-cell} mut_metadata = pyslim.mutation_metadata(ts) mut = ts.mutation(0) -for k in mut.metadata["derived_states"]: +for k in mut.metadata["slim_ids"]: md = mut_metadata[k] print(f"SLiM ID: {k}") print(f"Metadata: {md}") diff --git a/docs/traits.md b/docs/traits.md index 86be64c..f31ba42 100644 --- a/docs/traits.md +++ b/docs/traits.md @@ -192,8 +192,8 @@ mut util.pp(mut) ``` We can see which SLiM mutation(s) this mutation represents -by looking up those SLiM mutation IDs in the mutation's "derived state", -`mut.metadata["derived_states"]`. +by looking up those SLiM mutation IDs in the mutation's metadata: +`mut.metadata["slim_ids"]`. (As noted [previously](sec_tutorial_mutation_metadata), the same information is stored in `mut.derived_state`, but we recommend pulling this information out of metadata.) @@ -201,7 +201,7 @@ Then, we can find information about those SLiM mutations in the top-level mutation metadata (here, `mut_metadata`, obtained using {func}`.mutation_metadata`). ```{code-cell} -mut_metadata[mut.metadata["derived_states"][0]] +mut_metadata[mut.metadata["slim_ids"][0]] ``` So, we can use the `effect_size` and `dominance` to calculate genetic effects (we won't need `hemizygous dominance`, diff --git a/docs/tutorial.md b/docs/tutorial.md index 9011905..371cb94 100644 --- a/docs/tutorial.md +++ b/docs/tutorial.md @@ -437,9 +437,9 @@ to the SLiM mutations with {func}`generate_nucleotides`, and (2) move those nucleotides over into the "ancestral state" and "derived state" slots of the tree sequence with {func}`convert_alleles`. If all your mutations in SLiM were nucleotide mutations, you only need to do (2). -And, beware that (2) is an irreversible step: if you write the tree sequence -produced by {func}`convert_alleles` to a file, you can't load that file into SLiM -without putting those back somehow. So, to do this we'll do: +Since the information about SLiM mutation IDs is also stored in mutation metadata, +this entails no loss of information, and the tree sequence can still be read in by SLiM. +So, to do this we'll do: ```{code-cell} nts = pyslim.generate_nucleotides(ts) @@ -475,31 +475,29 @@ with open("example_sim2.vcf", "w") as vcffile: ## Mutation metadata Because of mutation stacking (see the SLiM manual), -each "tskit mutation" can represent a superposition of more than one -"SLiM mutation". -This is recorded by setting the derived state of the tskit mutation -to a comma-separated string of SLiM mutation IDs -(or the empty string, to denote "no mutations"), -and the `"derived_states"` entry of the tskit mutation's metadata to this same list. +each "tskit mutation" can represent a superposition of more than one "SLiM mutation". +This is recorded by recording the SLiM IDs of these mutations +in the `"slim_ids"` entry of the tskit mutation's metadata. +(Also, SLiM writes these to the derived state, as a comma-separated string.) So, each SLiM mutation can appear in more than one tskit mutation, and the [metadata](sec_metadata) about these SLiM mutations is stored in top-level metadata, rather than along with the tskit mutations. [](sec_tutorial_selected_mutations) has a more in-depth example, but here is a quick overview. To print out the information about each SLiM mutation "carried" by a given tskit mutation, whose SLiM IDs are stored in the mutation's metadata -under `"derived_states"`, we'd do: +under `"slim_ids"`, we'd do: ```{code-cell} :tags: ["remove-output"] mut_metadata = pyslim.mutation_metadata(ts) mut = ts.mutation(0) -for x in mut.metadata["derived_states"]: +for x in mut.metadata["slim_ids"]: print(f"SLiM mutation {x}:") print(mut_metadata[x]) ``` ```{code-cell} :tags: ["remove-input"] -for x in mut.metadata["derived_states"]: +for x in mut.metadata["slim_ids"]: print(f"SLiM mutation {x}:") util.pp(mut_metadata[x]) ``` @@ -508,7 +506,7 @@ carries more than one (stacked) SLiM mutation. You may have noticed that the same information goes both into the derived state entry of a tskit mutation (e.g., `mut.derived_state`) and into the metadata of that -same mutation (e.g., `mut.metadata["derived_states"]`). The reason for this +same mutation (e.g., `mut.metadata["slim_ids"]`). The reason for this redundancy is that tskit uses the `derived_state` entries (and the `ancestral_state` entries of sites) to do VCF output and determine allelic identity, and so it's helpful to be able to rewrite these @@ -1103,7 +1101,7 @@ some number of SLiM mutations, whose SLiM IDs are stored in the mutation's `meta For instance, here's which SLiM mutation(s) the first mutation in the tree sequence represents: ```{code-cell} -ds = ts.mutation(0).metadata["derived_states"] +ds = ts.mutation(0).metadata["slim_ids"] print(f"SLiM IDs: {ds}") ``` To see the information about these, we pull their information out @@ -1183,7 +1181,7 @@ Now, mutations have a ``nucleotide`` property in metadata that is not ``-1``: :tags: ["remove-output"] mut_metadata = pyslim.mutation_metadata(ts) m = ts.mutation(0) -md = [mut_metadata[k] for k in m.metadata["derived_states"]] +md = [mut_metadata[k] for k in m.metadata["slim_ids"]] print(m) for x in md: print(x) @@ -1203,7 +1201,7 @@ by indexing the {data}`.NUCLEOTIDES` object: for k in range(3): m = ts.mutation(k) print(f"Mutation {k}: position {ts.site(m.site).position}, time {m.time}") - for sid in m.metadata["derived_states"]: + for sid in m.metadata["slim_ids"]: md = mut_metadata[sid] print(f" nucleotide: {pyslim.NUCLEOTIDES[md['nucleotide']]}") ``` @@ -1255,7 +1253,7 @@ Here's the first mutation: :tags: ["remove-output"] mut_metadata = pyslim.mutation_metadata(ts) m = ts.mutation(0) -md = [mut_metadata[k] for k in m.metadata["derived_states"]] +md = [mut_metadata[k] for k in m.metadata["slim_ids"]] print(m) for x in md: print(x) @@ -1283,7 +1281,7 @@ and we can pull up information about that with the `ts.site( )` method: s = ts.site(m.site) md = [ mut_metadata[k] for m in s.mutations - for k in m.metadata["derived_states"] + for k in m.metadata["slim_ids"] ] print(s) for x in md: @@ -1299,7 +1297,7 @@ for x in md: This mutation occurred at the position along the genome shown in `site.position`, which previously had no mutations (since `site.ancestral_state` is the empty string, `''`) -and was given the SLiM mutation ID shown in `m.metadata["derived_states"]`. +and was given the SLiM mutation ID shown in `m.metadata["slim_ids"]`. The metadata (with `x` the mutation ID, `mut_metadata[x]`, a dict) tells us the mutation's selection coefficient and which population and at what SLiM time it occurred. This is not a nucleotide model, so the nucleotide entry is `-1`. @@ -1323,8 +1321,8 @@ for m in ts.mutations(): break pm = ts.mutation(m.parent) -md = [mut_metadata[k] for k in m.metadata["derived_states"]] -pmd = [mut_metadata[k] for k in pm.metadata["derived_states"]] +md = [mut_metadata[k] for k in m.metadata["slim_ids"]] +pmd = [mut_metadata[k] for k in pm.metadata["slim_ids"]] print(m) for x in md: @@ -1417,7 +1415,7 @@ mut_type = np.zeros(ts.num_sites) for j, s in enumerate(ts.sites()): mt = [] for m in s.mutations: - for sid in m.metadata["derived_states"]: + for sid in m.metadata["slim_ids"]: md = mut_metadata[sid] mt.append(md["mutation_type"]) if len(set(mt)) > 1: @@ -1454,7 +1452,7 @@ Finally, let's pull out information on the allele with the largest selection coe :tags: ["remove-output"] sel_coeffs = np.array([ sum(mut_metadata[k]["per_trait"][0]["effect_size"] - for k in m.metadata["derived_states"]) + for k in m.metadata["slim_ids"]) for m in ts.mutations() ]) which_max = np.argmax(sel_coeffs) diff --git a/docs/vignette_coalescent_diversity.md b/docs/vignette_coalescent_diversity.md index 458e10a..ee3f185 100644 --- a/docs/vignette_coalescent_diversity.md +++ b/docs/vignette_coalescent_diversity.md @@ -275,7 +275,7 @@ This runs quickly, since it's only 100 generations. First, let's look at what mutations are present. ```{code-cell} ts = tskit.load("vignette_annotated.trees") -num_stacked = np.array([len(m.metadata["derived_states"]) for m in ts.mutations()]) +num_stacked = np.array([len(m.metadata["slim_ids"]) for m in ts.mutations()]) init_time = ts.metadata['SLiM']['tick'] old_mut = np.array([m.time > init_time - 1 - 1e-12 for m in ts.mutations()]) assert sum(old_mut) == ots.num_mutations @@ -311,7 +311,7 @@ p = ts.sample_count_stat(nodes_by_time, lambda x: x/num_nodes, 2, windows='sites strict=False, span_normalise=False, polarised=True) mut_metadata = pyslim.mutation_metadata(ts) s = np.array([sum([sum([mut_metadata[k]["per_trait"][0]["effect_size"] - for k in m.metadata["derived_states"]]) + for k in m.metadata["slim_ids"]]) for m in site.mutations]) for site in ts.sites()]) ``` @@ -422,7 +422,7 @@ and print a picture of it, with mutations labeled by their type: ```{code-cell} mut_metadata = pyslim.mutation_metadata(mts) for t in mts.trees(): - mt = [max([mut_metadata[k]['mutation_type'] for k in m.metadata["derived_states"]]) for m in t.mutations()] + mt = [max([mut_metadata[k]['mutation_type'] for k in m.metadata["slim_ids"]]) for m in t.mutations()] if t.num_mutations > 12: break @@ -478,26 +478,26 @@ There are indeed: :tags: ['remove-output'] for site in mts.sites(): if len(site.mutations) > 1: - types = [set([mut_metadata[k]["mutation_type"] for k in mut.metadata["derived_states"]]) + types = [set([mut_metadata[k]["mutation_type"] for k in mut.metadata["slim_ids"]]) for mut in site.mutations] if max(map(len, types)) > 1: print(site) for mut in site.mutations: print(mut) - for k in mut.metadata["derived_states"]: + for k in mut.metadata["slim_ids"]: print(mut_metadata[k]) ``` ```{code-cell} :tags: ['remove-input'] for site in mts.sites(): if len(site.mutations) > 1: - types = [set([mut_metadata[k]["mutation_type"] for k in mut.metadata["derived_states"]]) + types = [set([mut_metadata[k]["mutation_type"] for k in mut.metadata["slim_ids"]]) for mut in site.mutations] if max(map(len, types)) > 1: util.pp(site) for mut in site.mutations: util.pp(mut) - for k in mut.metadata["derived_states"]: + for k in mut.metadata["slim_ids"]: util.pp(mut_metadata[k]) ``` diff --git a/pyslim/methods.py b/pyslim/methods.py index f6af881..d307359 100644 --- a/pyslim/methods.py +++ b/pyslim/methods.py @@ -431,10 +431,10 @@ def add_mutation_metadata(ts, mutation_type=0, remove_unused=False): :class:`msprime.SLiMv6MutationModel`. Any information about SLiM mutations already in top-level metadata will remain unchanged. - To do this, this method looks for all SLiM IDs that are found in the derived - state of some mutation but are not represented in the top-level metadata - (see :func:`.mutation_metadata`). This function then adds entries to that top-level - metadata with default values (see :func:`.default_slim_metadata`), + To do this, this method looks for all SLiM IDs that are found in the `"slim_ids"` + entry of some tskit mutation's metadata but are not represented in the top-level + metadata (see :func:`.mutation_metadata`). This function then adds entries to that + top-level metadata with default values (see :func:`.default_slim_metadata`), except that (a) the ``mutation_type`` can be specified; and (b) the ``slim_time`` is set using the ``tick`` value in top-level metadata and the ``time`` of the oldest tskit mutation in which the SLiM mutation occurs. @@ -442,7 +442,7 @@ def add_mutation_metadata(ts, mutation_type=0, remove_unused=False): :param tskit.TreeSequence ts: The tree sequence to transform. :param int mutation_type: The numeric ID of the mutation type in SLiM. :param bool remove_unused: Whether to also remove from metadata information about any - mutations not seen in the derived states of the tree sequence. + mutations not referenced by mutations in the tree sequence. :return tskit.TreeSequence: A copy of the tree sequence with mutation information in metadata. """ @@ -460,7 +460,7 @@ def add_mutation_metadata_tables(tables, mutation_type=0, remove_unused=False): :param tskit.TableCollection tables: The table collection to be modified. :param int mutation_type: The numeric ID of the mutation type in SLiM. :param bool remove_unused: Whether to also remove from metadata information about any - mutations not seen in the derived states of the tree sequence. + mutations not referenced by mutations in the tree sequence. """ ts_metadata = tables.metadata if ( @@ -475,9 +475,7 @@ def add_mutation_metadata_tables(tables, mutation_type=0, remove_unused=False): num_traits = len(ts_metadata["SLiM"]["traits"]) existing_muts = {x["mutation_id"] for x in ts_metadata["SLiM_mutation_list"]} mut_ids = [ - (int(j), mut.time) - for mut in tables.mutations - for j in mut.metadata["derived_states"] + (int(j), mut.time) for mut in tables.mutations for j in mut.metadata["slim_ids"] ] if len(mut_ids) > 0: mut_ids.sort() @@ -516,8 +514,8 @@ def convert_alleles(ts): their corresponding nucleotides. For sites, SLiM-produced tree sequences have "" (the empty string) for the ancestral state at each site; this method will replace this with the corresponding nucleotide from the reference sequence. - For mutations, SLiM records the 'derived state' as a SLiM mutation ID; this - method will replace the derived state with the nucleotide from the mutation's + For mutations, SLiM records the 'derived state' as a list of SLiM mutation IDs; + this method will replace the derived state with the nucleotide from the mutation's metadata. In SLiM's output the list of mutation IDs is recorded both in each mutation's @@ -553,7 +551,7 @@ def convert_alleles(ts): alleles = np.array([x["nucleotide"] for x in mut_metadata.values()], dtype="int") # First, do this for the unstacked mutations quickly # mut_index will map from tskit-mutations to slim-mutations - mut_index = np.array([mut.metadata["derived_states"][-1] for mut in ts.mutations()]) + mut_index = np.array([mut.metadata["slim_ids"][-1] for mut in ts.mutations()]) assert np.all(mut_index >= 0), "This should not occur: please file a bug report." nucs = alleles[np.searchsorted(mut_ids, mut_index)] # Now, update those where necessary @@ -605,7 +603,7 @@ def generate_nucleotides(ts, reference_sequence=None, keep=True, seed=None): Technical note: in the case of stacked mutations, the SLiM mutation that determines the nucleotide state of the (tskit) mutation is the last one in the list - of "derived states" in the tskit mutation metadata. This method tries to + of "slim_ids" in the tskit mutation metadata. This method tries to assign nucleotides so that each mutation differs from the previous state, but this is not always possible if some mutations already have nucleotides and others do not. @@ -659,8 +657,8 @@ def generate_nucleotides(ts, reference_sequence=None, keep=True, seed=None): pds = [] else: pa = states[mut.parent] - pds = ts.mutation(mut.parent).metadata["derived_states"] - for i in mut.metadata["derived_states"]: + pds = ts.mutation(mut.parent).metadata["slim_ids"] + for i in mut.metadata["slim_ids"]: md = mut_info[i] da = md["nucleotide"] if da == -1 or not keep: @@ -1232,7 +1230,7 @@ def next_slim_mutation_id(ts): `next_id` in your :class:`msprime.SLiMv6MutationModel` to be larger than any existing mutation IDs. Setting `next_id` equal to the output of this function will allow the mutated tree sequence to be read in by SLiM. - To do this, recall that the `derived_states` attribute of each mutation's metadata + To do this, recall that the `slim_ids` attribute of each mutation's metadata is a list of SLiM mutation IDs; this function just parses all these metadata entries and returns one larger than the largest integer found. """ @@ -1241,7 +1239,7 @@ def next_slim_mutation_id(ts): try: max_id = functools.reduce( max, - (x for mut in ts.mutations() for x in mut.metadata["derived_states"]), + (x for mut in ts.mutations() for x in mut.metadata["slim_ids"]), -1, ) except TypeError: diff --git a/pyslim/slim_metadata.py b/pyslim/slim_metadata.py index ef7afe2..4c39f7c 100644 --- a/pyslim/slim_metadata.py +++ b/pyslim/slim_metadata.py @@ -370,17 +370,17 @@ def is_vacant_num_bytes(num_chromosomes): "codec": "struct", "type": "object", "description": "SLiM schema for representing binary derived state data in mutation metadata (the metadata for each unique SLiM mutation is stored in top-level metadata).", - "examples": [{"derived_states": [0, 1, 17]}], + "examples": [{"slim_ids": [0, 1, 17]}], "properties": { - "derived_states": { + "slim_ids": { "index": 1, "type": "array", "noLengthEncodingExhaustBuffer": True, - "description": "An array of SLiM mutation IDs (int64t), representing the (stacked) mutations contained by the derived state for the mutation.", + "description": "An array of SLiM mutation IDs (int64t), representing the (stacked) mutations represented by this tskit mutation.", "items": {"binaryFormat": "q", "type": "number"}, } }, - "required": ["derived_states"], + "required": ["slim_ids"], }, "node": { "$schema": "http://json-schema.org/schema#", @@ -831,7 +831,7 @@ def default_slim_metadata(name, num_chromosomes=1, num_traits=1, **kwargs): elif name == "site": out = None elif name == "mutation": - out = {"derived_states": []} + out = {"slim_ids": []} elif name == "mutation_list_entry": out = { "mutation_id": 0, @@ -2094,7 +2094,7 @@ def update_tables(tables): muts.metadata_schema = old_schema tables.mutations.metadata_schema = slim_metadata_schemas["mutation"] for mut in muts: - md = {"derived_states": [int(j) for j in mut.derived_state.split(",")]} + md = {"slim_ids": [int(j) for j in mut.derived_state.split(",")]} tables.mutations.append(mut.replace(metadata=md)) if file_version == "0.1": diff --git a/pyslim/slim_tree_sequence.py b/pyslim/slim_tree_sequence.py index 617661f..331881d 100644 --- a/pyslim/slim_tree_sequence.py +++ b/pyslim/slim_tree_sequence.py @@ -15,7 +15,7 @@ def mutation_metadata(ts, check=True, ts_metadata=None): Returns a dictionary whose keys are the numeric SLiM IDs of mutations, and whose values are metadata entries for those mutations. These SLiM IDs are found in the metadata of tskit mutations: - for each mutation ``mut``, as ``mut.metadata["derived_states"]``. + for each mutation ``mut``, as ``mut.metadata["slim_ids"]``. This is a simple extraction function that places the list of metadata entries stored in ``ts.metadata["SLiM_mutation_list"]`` in a dictionary @@ -44,7 +44,7 @@ def mutation_metadata(ts, check=True, ts_metadata=None): ml.sort(key=lambda x: x["mutation_id"]) out = {mut["mutation_id"]: mut for mut in ml} if check: - ids = {j for mut in ts.mutations() for j in mut.metadata["derived_states"]} + ids = {j for mut in ts.mutations() for j in mut.metadata["slim_ids"]} for k in ids: if k not in out: raise ValueError( @@ -141,7 +141,7 @@ def nucleotide_at(ts, node, position, time=None, mut_metadata=None): if mut_id == tskit.NULL: out = NUCLEOTIDES.index(ts.reference_sequence.data[int(position)]) else: - md = ts.mutation(mut_id).metadata["derived_states"] + md = ts.mutation(mut_id).metadata["slim_ids"] out = -1 for k in md[::-1]: out = mut_metadata[k]["nucleotide"] diff --git a/tests/test_annotation.py b/tests/test_annotation.py index ec8b7c5..c465123 100644 --- a/tests/test_annotation.py +++ b/tests/test_annotation.py @@ -252,8 +252,8 @@ def verify_remapping(self, ts, rts, subpop_map): rmut_info = pyslim.mutation_metadata(rts) assert len(mut_info) == len(rmut_info) # mutations may have changed order - tsm = {m.derived_state: m.metadata["derived_states"] for m in ts.mutations()} - rtsm = {m.derived_state: m.metadata["derived_states"] for m in rts.mutations()} + tsm = {m.derived_state: m.metadata["slim_ids"] for m in ts.mutations()} + rtsm = {m.derived_state: m.metadata["slim_ids"] for m in rts.mutations()} for x in tsm: assert x in rtsm assert tsm[x] == rtsm[x] @@ -314,7 +314,7 @@ def test_warns_overwriting_mutations(self, helper_functions): t.mutations.metadata_schema = pyslim.slim_metadata_schemas["mutation"] t.mutations.clear() for j, mut in enumerate(ts.mutations()): - t.mutations.append(mut.replace(metadata={"derived_states": [j]})) + t.mutations.append(mut.replace(metadata={"slim_ids": [j]})) ts = t.tree_sequence() ts = pyslim.add_mutation_metadata(ts) with pytest.warns(Warning, match="already has.*metadata"): @@ -1291,7 +1291,7 @@ def test_add_mutation_metadata_removes(self): mut_info = pyslim.mutation_metadata(ts) nmut_info = pyslim.mutation_metadata(nsts) mut_ids = np.unique( - [k for mut in nsts.mutations() for k in mut.metadata["derived_states"]] + [k for mut in nsts.mutations() for k in mut.metadata["slim_ids"]] ) assert len(mut_ids) == len(nsts.metadata["SLiM_mutation_list"]) for k in mut_ids: diff --git a/tests/test_tree_sequence.py b/tests/test_tree_sequence.py index 26acd76..2ef92a9 100644 --- a/tests/test_tree_sequence.py +++ b/tests/test_tree_sequence.py @@ -56,13 +56,13 @@ def naive_mutation_at(ts, node, pos, time=None): def verify_mutation_metadata(ts): - # Verify that all derived states are properly accounted for + # Verify that all SLiM mutation IDs are properly accounted for # in mutation metadata. mdl = ts.metadata["SLiM_mutation_list"] mut_info = pyslim.mutation_metadata(ts) assert len(mut_info) == len(mdl) for mut in ts.mutations(): - for j in mut.metadata["derived_states"]: + for j in mut.metadata["slim_ids"]: assert j in mut_info @@ -120,9 +120,7 @@ def test_slim_time(self, recipe): stage = "early" if "init_mutated" in recipe else None slim_times = pyslim.slim_time(ts, ts.mutations_time, stage=stage) for t, mut in zip(slim_times, ts.mutations()): - mut_time = max( - [muts[j]["slim_time"] for j in mut.metadata["derived_states"]] - ) + mut_time = max([muts[j]["slim_time"] for j in mut.metadata["slim_ids"]]) assert mut_time == t @@ -847,7 +845,7 @@ def test_mutation_consistency(self, recipe): } mut_info = pyslim.mutation_metadata(ts) for mut in ts.mutations(): - for k in mut.metadata["derived_states"]: + for k in mut.metadata["slim_ids"]: assert k in debug_info or mut_info[k]["mutation_id"] == 2 assert k in mut_info assert debug_info[k]["chromosome_id"] == chrom_id @@ -963,12 +961,12 @@ def test_nucleotide_at(self, recipe): for k in np.where(node == ts.tables.mutations.node)[0]: mut = ts.mutation(k) if ts.site(mut.site).position == pos: - j = len(mut.metadata["derived_states"]) - 1 - k = mut.metadata["derived_states"][j] + j = len(mut.metadata["slim_ids"]) - 1 + k = mut.metadata["slim_ids"][j] b = mut_metadata[k]["nucleotide"] while j > 0 and b == -1: j -= 1 - k = mut.metadata["derived_states"][j] + k = mut.metadata["slim_ids"][j] b = mut_metadata[k]["nucleotide"] assert a == b @@ -1007,7 +1005,7 @@ def test_nucleotide_spectrum(self, recipe): pos = ts.site(mut.site).position if pos > 0 and pos < ts.sequence_length - 1: nmuts += 1 - mut_list = [mut_metadata[k] for k in mut.metadata["derived_states"]] + mut_list = [mut_metadata[k] for k in mut.metadata["slim_ids"]] k = np.argmax([u["slim_time"] for u in mut_list]) derived_nuc = mut_list[k]["nucleotide"] left_nuc = pyslim.nucleotide_at( @@ -1055,13 +1053,13 @@ def last_slim_mutations(self, ts): mut_info = pyslim.mutation_metadata(ts) for mut in ts.mutations(): slim_muts = { - k: v for k, v in mut_info.items() if k in mut.metadata["derived_states"] + k: v for k, v in mut_info.items() if k in mut.metadata["slim_ids"] } if mut.parent == tskit.NULL: parent_slim_ids = [] else: parent_mut = ts.mutation(mut.parent) - parent_slim_ids = parent_mut.metadata["derived_states"] + parent_slim_ids = parent_mut.metadata["slim_ids"] max_time = max([md["slim_time"] for md in slim_muts.values()]) any_new = any( [ @@ -1213,7 +1211,7 @@ def verify_generate_nucleotides(self, ts, check_transitions=False): } for mut in ts.mutations(): aa = ts.reference_sequence.data[int(ts.site(mut.site).position)] - for i in mut.metadata["derived_states"]: + for i in mut.metadata["slim_ids"]: md = mut_info[i] nuc = md["nucleotide"] assert nuc in [0, 1, 2, 3] @@ -1225,10 +1223,10 @@ def verify_generate_nucleotides(self, ts, check_transitions=False): assert pyslim.NUCLEOTIDES[nuc] != aa else: mp = ts.mutation(mut.parent) - if mp.metadata["derived_states"] != mut.metadata["derived_states"]: + if mp.metadata["slim_ids"] != mut.metadata["slim_ids"]: assert (ts_muts[mut.parent] != ts_muts[mut.id]) or ( - len(mut.metadata["derived_states"]) - > 1 + len(mp.metadata["derived_states"]) + len(mut.metadata["slim_ids"]) + > 1 + len(mp.metadata["slim_ids"]) ) @pytest.mark.parametrize("recipe", recipe_eq(exclude="old_mutations"), indirect=True) @@ -1289,11 +1287,11 @@ def test_generate_nucleotides_keep(self): mut_info2 = pyslim.mutation_metadata(nts2) muts1 = {} for mut in nts1.mutations(): - for i in mut.metadata["derived_states"]: + for i in mut.metadata["slim_ids"]: md = mut_info1[i] muts1[i] = md["nucleotide"] for mut in nts2.mutations(): - for i in mut.metadata["derived_states"]: + for i in mut.metadata["slim_ids"]: md = mut_info2[i] if md["mutation_type"] == 1: assert i in muts1