################################################################################
#
# RevBayes: Discrete morphology tutorials -- among-character state-frequency
#           variation using a discretized Beta distribution
##
################################################################################

# Specify the number of categories we would like to use to describe our data.
# For simplicity, we will use 5.
num_cats = 5

# Specify the two parameters to the Beta distribution, and the moves on these
# parameters.
beta_scale ~ dnLognormal( 0.0, sd=2*0.587405 )
beta_scale.setValue( 1.0 )
moves.append( mvScale(beta_scale, lambda=1, weight=5.0 ) )

# Create the Beta distribution, according to the two parameters and the number
# of categories.
cats := fnDiscretizeBeta(beta_scale, beta_scale, num_cats)

# Create the Q matrices. These data are binary, so we initialize the
# Jukes-Cantor matrix with two states. The matrix is initialized from the
# median values of our categories from the Beta distribution.
for (i in 1:cats.size()) {
    Q[i] := fnF81(simplex(abs(1-cats[i]), cats[i]))
}

# Tell the model what the probability of a character going into any particular
# category is. This prior says that a character is equally likely to be put
# into any category.
matrix_probs <- simplex( rep(1, num_cats) )

# Set up Gamma-distributed rate variation.
alpha_morpho ~ dnUniform( 0.0, 1E8 )
alpha_morpho.setValue( 1.0 )
rates_morpho := fnDiscretizeGamma( alpha_morpho, alpha_morpho, 4 )
moves.append( mvScale(alpha_morpho, lambda=1, weight=2.0) )

# Combine all of our elements into a CTMC. Because we have not observed any
# invariant sites, we specify the coding as 'variable'. Note that we are now
# feeding the site matrices to the CTMC.
phyMorpho ~ dnPhyloCTMC(tree=phylogeny, siteRates=rates_morpho, Q=Q, type="Standard", coding="variable", siteMatrices=matrix_probs)
phyMorpho.clamp(morpho)