Skip to content

Commit 02597d1

Browse files
committed
Minor changes to script comments. Started fHiSSE section.
1 parent 72b1bfa commit 02597d1

File tree

2 files changed

+60
-4
lines changed

2 files changed

+60
-4
lines changed

tutorials/sse/fossils.md

Lines changed: 58 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -16,11 +16,13 @@ include_files:
1616
- data/canidae_tree.nex
1717
- scripts/mcmc_fBiSSE.Rev
1818
- scripts/mcmc_fBiSSE_tp.Rev
19+
- scripts/mcmc_fHiSSE.Rev
1920
index: true
2021
---
2122

2223
{% assign bisse_script = "mcmc_fBiSSE.Rev" %}
2324
{% assign bisse_tp_script = "mcmc_fBiSSE_tp.Rev" %}
25+
{% assign hisse_script = "mcmc_fHiSSE.Rev" %}
2426

2527
{% section Introduction | introduction %}
2628

@@ -54,9 +56,9 @@ We provide the data files which we will use in this tutorial:
5456
of the tree, and our main interest on the effects of hypercarnivory specifically, we simplified
5557
the dataset.
5658

57-
> Create a new directory on your computer called `RB_fmusse_tutorial`.
59+
> Create a new directory on your computer called `RB_fsse_tutorial`.
5860
>
59-
> Within the `RB_fmusse_tutorial` directory, create a subdirectory called `data`.
61+
> Within the `RB_fsse_tutorial` directory, create a subdirectory called `data`.
6062
> Then, download the provided files and place them in the `data` folder.
6163
> You can then start a `scripts` folder to hold the script for this analysis, or use RevBayes interactively by running `rb` and adding the code described below one block at a time.
6264
{:.instruction}
@@ -311,3 +313,57 @@ to change the output, checkpointing, and ancestral state filenames if you run bo
311313
{% endaside %}
312314

313315
{% section Setting up a serially-sampled HiSSE model | sec_fHiSSE %}
316+
317+
The hidden SSE (HiSSE) model {% cite Beaulieu2016 %} was introduced to explicitly accommodate the role of
318+
unobserved effects on diversification rates.
319+
This alleviated some of the more complicated issues in BiSSE, as it was shown to falsely identify trait-associated
320+
hetereogeneity when the source of heterogeneity was not included in the model {% cite Rabosky2015 %}.
321+
See [the HiSSE tutorial]({{ base.url }}/tutorials/sse/hisse) for more information on the model.
322+
323+
We only need to make minor modifications to our BiSSE script to make it into a HiSSE model.
324+
The easiest way to proceed would therefore be to make a copy of your BiSSE script, renaming it as appropriate,
325+
and make the changes we mention below.
326+
If a part of the script is not mentioned, assume it remains unchanged (e.g. the lines to read `tree` and `data`
327+
remain the same, since we are using the same fixed tree and same observed characters).
328+
329+
First, in addition to our `num_states` variables as before, we must create new variables to
330+
hold the number of hidden states, and the total number of states.
331+
332+
{{ hisse_script | snippet:"block#", "4-5" }}
333+
334+
We are using only two hidden states, which we call A and B.
335+
While it is possible to include more hidden states, since these are unobserved states for which we have no data,
336+
one should proceed with caution including more parameters.
337+
338+
Then, we must expand our `data` object to include the hidden states
339+
340+
{{ hisse_script | snippet:"block#", "8" }}
341+
342+
Now the character in `data_exp` has four states.
343+
By default, RevBayes sets up state numbers to follow the order of observed, then hidden.
344+
So here, states 1 and 2 correspond to 0A and 1A, and states 3 and 4 to 0B and 1B.
345+
346+
To set up our speciation, extinction, and fossil sampling rates, we follow the exact same procedure as in the
347+
BiSSE case, but with `num_rates` rates instead.
348+
349+
{{ hisse_script | snippet:"block#", "9-13" }}
350+
351+
Note that this sets each rate to be independent.
352+
While that is the most flexible set up, some analyses might run into "label-swapping" issues.
353+
Label-swapping occurs because the likelihood of the model does not change if the states A and B are switched,
354+
since there is no data associated with those states.
355+
As such, the assignment of each state might switch throughout the MCMC, making it so that each rate parameter's
356+
posterior distribution is actually a mixture between the corresponding rates.
357+
If you encounter such problems in your analyses, see [the HiSSE tutorial]({{ base.url }}/tutorials/sse/hisse) for
358+
a way to set rates that guarantees there will be no label-swapping (sacrificing a small amount of flexibility in
359+
the process).
360+
361+
To set up the rate matrix, we can create variables for `q_obs` and `q_hidden`, corresponding to $q_{01}$ and \
362+
$q_{10}$, and $q_{AB}$ and $q_{BA}$ respectively.
363+
We can then use `fnHiddenStateRateMatrix` to build the full rate matrix.
364+
365+
{{ hisse_script | snippet:"block#", "14-18" }}
366+
367+
Note that this set up assumes that transitions between rate categories with the same observed state occur at the
368+
same rate (i.e. $q_{0A->0B} = q_{1A->1B} = q_{AB}$), and similarly for categories with the same hidden state.
369+
While this is often a useful assumption (and one that was

tutorials/sse/scripts/mcmc_fHiSSE.Rev

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
################################################################################
22
#
3-
# BiSSE analysis of the influence of diet on the diversification
3+
# HiSSE analysis of the influence of diet on the diversification
44
# rates of extant and extinct canids
55
#
66
# author: Bruno do Rosario Petrucci
@@ -12,7 +12,7 @@
1212

1313
# read the phylogenetic tree
1414
tree <- readTrees("data/canidae_tree.nex")[1]
15-
# this includes 94 canid species, with 3/36 extant species, and 89 extinct species
15+
# this includes 94 canid species, with 4/36 extant species, and 89 extinct species
1616
# representing all three canid subfamilies: caninae, borophaginae, and hesperocyoninae
1717

1818
# read the character data

0 commit comments

Comments
 (0)