Skip to content

Commit 43d5e27

Browse files
committed
fix #529, extend calcFST() to support multiple chromosomes
1 parent d3130a1 commit 43d5e27

8 files changed

Lines changed: 71 additions & 48 deletions

File tree

QtSLiM/help/SLiMHelpClasses.html

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1064,8 +1064,8 @@
10641064
<p class="p6">Declare the current simulation finished.<span class="Apple-converted-space">  </span>This method is equivalent to the <span class="s1">Community</span> method <span class="s1">simulationFinished()</span>, except that this method is only legal to call in single-species models (to provide backward compatibility).<span class="Apple-converted-space">  </span>It is recommended that new code should call the <span class="s1">Community</span> method; this method may be deprecated in the future.</p>
10651065
<p class="p3">– (void)skipTick(void)</p>
10661066
<p class="p6">Deactivate the target species for the current tick.<span class="Apple-converted-space">  </span>This sets the <span class="s1">active</span> property of the species to <span class="s1">F</span>; it also set the <span class="s1">active</span> property of all callbacks that belong to the species (with the species as their <span class="s1">species</span> specifier) to <span class="s1">F</span>, and sets the active property of all events that are synchronized with the species (with the species as their <span class="s1">ticks</span> specifier) to <span class="s1">F</span>.<span class="Apple-converted-space">  </span>The cycle counter for the species will not be incremented at the end of the tick.<span class="Apple-converted-space">  </span>This method may only be called in <span class="s1">first()</span> events, to ensure that species are either active or inactive throughout a given tick.</p>
1067-
<p class="p5">– (object&lt;Mutation&gt;)subsetMutations([No&lt;Mutation&gt;$ exclude = NULL], [Nio&lt;MutationType&gt;$ mutType = NULL], [Ni$ position = NULL], [Nis$ nucleotide = NULL], [Ni$ tag = NULL], [Ni$ id = NULL], [Niso&lt;Chromosome&gt;$ chromosome = NULL])</p>
1068-
<p class="p6">Returns a vector of mutations subset from the list of all active mutations in the species (as would be provided by the <span class="s1">mutations</span> property).<span class="Apple-converted-space">  </span>The parameters specify constraints upon the subset of mutations that will be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">exclude</span>, if non-<span class="s1">NULL</span>, may specify a specific mutation that should not be included (typically the focal mutation in some operation).<span class="Apple-converted-space">  </span>Parameter <span class="s1">mutType</span>, if non-<span class="s1">NULL</span>, may specify a mutation type for the mutations to be returned (as either a <span class="s1">MutationType</span> object or an <span class="s1">integer</span> identifier).<span class="Apple-converted-space">  </span>Parameter <span class="s1">position</span>, if non-<span class="s1">NULL</span>, may specify a base position for the mutations to be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">nucleotide</span>, if non-<span class="s1">NULL</span>, may specify a nucleotide for the mutations to be returned (either as a string, <span class="s1">"A"</span> / <span class="s1">"C"</span> / <span class="s1">"G"</span> / <span class="s1">"T"</span>, or as an integer, <span class="s1">0</span> / <span class="s1">1</span> / <span class="s1">2</span> / <span class="s1">3</span> respectively).<span class="Apple-converted-space">  </span>Parameter <span class="s1">tag</span>, if non-<span class="s1">NULL</span>, may specify a tag value for the mutations to be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">id</span>, if non-<span class="s1">NULL</span>, may specify a required value for the <span class="s1">id</span> property of the mutations to be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">chromosome</span>, if non-<span class="s1">NULL</span>, may specify a chromosome with which the mutations returned must be associated (as either an <span class="s1">integer</span> id, a <span class="s1">string</span> symbol, or a <span class="s1">Chromosome</span> object).</p>
1067+
<p class="p5">– (object&lt;Mutation&gt;)subsetMutations([No&lt;Mutation&gt;$ exclude = NULL], [Nio&lt;MutationType&gt;$ mutType = NULL], [Ni$ position = NULL], [Nis$ nucleotide = NULL], [Ni$ tag = NULL], [Ni$ id = NULL], [Niso&lt;Chromosome&gt; chromosome = NULL])</p>
1068+
<p class="p6">Returns a vector of mutations subset from the list of all active mutations in the species (as would be provided by the <span class="s1">mutations</span> property).<span class="Apple-converted-space">  </span>The parameters specify constraints upon the subset of mutations that will be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">exclude</span>, if non-<span class="s1">NULL</span>, may specify a specific mutation that should not be included (typically the focal mutation in some operation).<span class="Apple-converted-space">  </span>Parameter <span class="s1">mutType</span>, if non-<span class="s1">NULL</span>, may specify a mutation type for the mutations to be returned (as either a <span class="s1">MutationType</span> object or an <span class="s1">integer</span> identifier).<span class="Apple-converted-space">  </span>Parameter <span class="s1">position</span>, if non-<span class="s1">NULL</span>, may specify a base position for the mutations to be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">nucleotide</span>, if non-<span class="s1">NULL</span>, may specify a nucleotide for the mutations to be returned (either as a string, <span class="s1">"A"</span> / <span class="s1">"C"</span> / <span class="s1">"G"</span> / <span class="s1">"T"</span>, or as an integer, <span class="s1">0</span> / <span class="s1">1</span> / <span class="s1">2</span> / <span class="s1">3</span> respectively).<span class="Apple-converted-space">  </span>Parameter <span class="s1">tag</span>, if non-<span class="s1">NULL</span>, may specify a tag value for the mutations to be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">id</span>, if non-<span class="s1">NULL</span>, may specify a required value for the <span class="s1">id</span> property of the mutations to be returned.<span class="Apple-converted-space">  </span>Parameter <span class="s1">chromosome</span>, if non-<span class="s1">NULL</span>, may specify a chromosome or chromosomes with which the mutations returned must be associated (as either <span class="s1">integer</span> ids, <span class="s1">string</span> symbols, or <span class="s1">Chromosome</span> objects).</p>
10691069
<p class="p6">This method is shorthand for getting the <span class="s1">mutations</span> property of the subpopulation, and then using operator <span class="s1">[]</span> to select only mutations with the desired properties; besides being much simpler than the equivalent Eidos code, it is also much faster.<span class="Apple-converted-space">  </span>Note that if you only need to select on mutation type, the <span class="s1">mutationsOfType()</span> method will be even faster.</p>
10701070
<p class="p5">– (object&lt;Substitution&gt;)substitutionsOfType(io&lt;MutationType&gt;$ mutType)</p>
10711071
<p class="p6">Returns an <span class="s1">object</span> vector of all the substitutions that are of the type specified by <span class="s1">mutType</span>, out of all of the substitutions that are currently present in the species.<span class="Apple-converted-space">  </span>This method is provided for speed; it is much faster than the corresponding Eidos code.<span class="Apple-converted-space">  </span>See also <span class="s1">mutationsOfType()</span>.</p>

QtSLiM/help/SLiMHelpFunctions.html

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -209,9 +209,9 @@
209209
<p class="p3">All haplosomes and mutations must be associated with the same chromosome.<span class="Apple-converted-space">  </span>If <span class="s3">muts</span> is <span class="s3">NULL</span> (the default), all mutations in the population associated with the same chromosome as the given haplosomes will be used.</p>
210210
<p class="p3">This function was written by Vitor Sudbrack (currently affiliated with University of Lausanne).</p>
211211
<p class="p4">(float$)calcFST(object&lt;Haplosome&gt; haplosomes1, object&lt;Haplosome&gt; haplosomes2, [No&lt;Mutation&gt; muts = NULL], [Ni$ start = NULL], [Ni$ end = NULL])</p>
212-
<p class="p3">Calculates the <i>F</i><span class="s4"><sub>ST</sub></span> between two <span class="s3">Haplosome</span> vectors – typically, but not necessarily, the haplosomes that constitute two different subpopulations (which we will assume for the purposes of this discussion).<span class="Apple-converted-space">  </span>In general, higher <i>F</i><span class="s4"><sub>ST</sub></span> indicates greater genetic divergence between subpopulations.</p>
213-
<p class="p3">The calculation is done using only the mutations in <span class="s3">muts</span>; if <span class="s3">muts</span> is <span class="s3">NULL</span>, all mutations are used.<span class="Apple-converted-space">  </span>The <span class="s3">muts</span> parameter can therefore be used to calculate the <i>F</i><span class="s4"><sub>ST</sub></span> only for a particular mutation type (by passing only mutations of that type).</p>
214-
<p class="p3">The calculation can be narrowed to apply to only a window – a subrange of the full chromosome – by passing the interval bounds [<span class="s3">start</span>, <span class="s3">end</span>] for the desired window.<span class="Apple-converted-space">  </span>In this case, the vector of mutations used for the calculation will be subset to include only mutations within the specified window.<span class="Apple-converted-space">  </span>The default behavior, with <span class="s3">start</span> and <span class="s3">end</span> of <span class="s3">NULL</span>, provides the haplosome-wide <i>F</i><span class="s4"><sub>ST</sub></span>, which is often used to assess the overall level of genetic divergence between sister species or allopatric subpopulations.</p>
212+
<p class="p3">Calculates the <i>F</i><span class="s4"><sub>ST</sub></span> between two <span class="s3">Haplosome</span> vectors – typically, but not necessarily, the haplosomes that constitute two different subpopulations (which we will assume for the purposes of this discussion).<span class="Apple-converted-space">  </span>In general, higher <i>F</i><span class="s4"><sub>ST</sub></span> indicates greater genetic divergence between subpopulations.<span class="Apple-converted-space">  </span>The haplosomes may be associated with more than one chromosome, in a multi-chromosome model; if so, <span class="s3">haplosomes1</span> and <span class="s3">haplosomes2</span> must be associated with the same set of chromosomes, defining the focal set of chromosomes for the calculation.</p>
213+
<p class="p3">The calculation is done using only the mutations in <span class="s3">muts</span>; if <span class="s3">muts</span> is <span class="s3">NULL</span>, all mutations associated with the focal chromosomes are used.<span class="Apple-converted-space">  </span>The <span class="s3">muts</span> parameter can be used to calculate the <i>F</i><span class="s4"><sub>ST</sub></span> only for a particular mutation type (by passing only mutations of that type), for example; it can focus the calculation on particular mutations of interest.<span class="Apple-converted-space">  </span>The mutations in <span class="s3">muts</span> must always be associated with the focal chromosomes.</p>
214+
<p class="p3">If there is a single focal chromosome, the calculation can be narrowed to apply to only a window – a subrange of the focal chromosome – by passing the interval bounds [<span class="s3">start</span>, <span class="s3">end</span>] for the desired window.<span class="Apple-converted-space">  </span>In this case, the vector of mutations used for the calculation will be subset to include only mutations within the specified window.<span class="Apple-converted-space">  </span>The default behavior, with <span class="s3">start</span> and <span class="s3">end</span> of <span class="s3">NULL</span>, provides the chromosome-wide <i>F</i><span class="s4"><sub>ST</sub></span>, which is often used to assess the overall level of genetic divergence between sister species or allopatric subpopulations.</p>
215215
<p class="p3">The code for <span class="s3">calcFST()</span> is, roughly, an Eidos implementation of Wright’s definition of <i>F</i><span class="s4"><sub>ST</sub></span> (but see below for further discussion and clarification):</p>
216216
<p class="p3"><i>F</i><span class="s4"><sub>ST</sub></span><span class="s1"> = 1 - <i>H</i></span><span class="s4"><sub>S</sub></span><span class="s1"> / <i>H</i></span><span class="s4"><sub>T</sub></span></p>
217217
<p class="p3">where <i>H</i><span class="s4"><i><sub>S</sub></i></span> is the average heterozygosity in the two subpopulations, and <i>H</i><span class="s4"><i><sub>T </sub></i></span>is the total heterozygosity when both subpopulations are combined.<span class="Apple-converted-space">  </span>In this implementation, the two haplosome vectors are weighted equally, not weighted by their size.<span class="Apple-converted-space">  </span>In SLiM 3, the implementation followed Wright’s definition closely, and returned the <i>average of ratios</i>: <span class="s3">mean(1.0 - H_s/H_t)</span>, in the Eidos code.<span class="Apple-converted-space">  </span>In SLiM 4, it returns the <i>ratio of averages</i> instead: <span class="s3">1.0 - mean(H_s)/mean(H_t)</span>.<span class="Apple-converted-space">  </span>In other words, the <i>F</i><span class="s4"><sub>ST</sub></span> value reported by SLiM 4 is an average across the specified mutations in the two sets of haplosomes, where <span class="s3">H_s</span> and <span class="s3">H_t</span> are first averaged across all specified mutations prior to taking the ratio of the two.<span class="Apple-converted-space">  </span>This ratio of averages is less biased than the average of ratios, and and is generally considered to be best practice (see, e.g., Bhatia et al., 2013).<span class="Apple-converted-space">  </span>This means that the behavior of <span class="s3">calcFST()</span> differs between SLiM 3 and SLiM 4.</p>

0 commit comments

Comments
 (0)