Question: Stochastic mapping & discrete characters
0
NK0 wrote:

Hi all,

I am trying to examine a correlation between two discrete traits, using stochastic mapping.

I tried two options and ran into problems with both:

1. Using stochastic mapping of Phytools –

make.simmap - The problem is that I need to use constraints - for one character

I need to set the state of the root to 0 (I figured how to do that), while for the second

character I need to set the reverse rate to 0. Couldn’t figure out this second constraints.

Tried this

mtrees<-make.simmap(tree,mtree\$states,nsim=10, model=matrix(c(0,1,0,0),2))

But model=matrix(c(0,1,0,0)) leads to the following error:

BLAS/LAPACK routine 'DGEBAL' gave error code -3

In matrix(XX\$rates[II], m, m, dimnames = list(lvls, lvls)) :

data length  is not a sub-multiple or multiple of the number of rows 

The same command works if I use , model=matrix(c(0,1,0,0),2))

In addition, I wonder if anyone has examined the distribution of two characters -
i.e. how can I calculate the percentage of time spent in (0,0), (0,1), etc. for the

classes generated by phytools.

2. Using Diversitree -

I used asr.stoch with a likelihood function and parameters. I managed to make constraints in the

two needed ways (i.e. on the root, and on the rates)

pars <- c(.1, .1, .03, .03, .03, .06)

set.seed(1)

phy <- trees(pars, "bisse", max.taxa=4, max.t=Inf, x0=0)[]

h <- history.from.sim.discrete(phy, 0:1)

plot(h, phy, main="True history")

lik <- make.mk2(phy, phy\$tip.state)

#lik<-constrain(q13~q12,q21~q12,q23~q12,q31~q12,q32~q12)

argnames(lik)

lik2<-constrain(lik, q10~0)

argnames(lik2)

fit = find.mle(lik2,c(.1))

st.s1 <- asr.stoch(lik2, fit\$par)

plot(st.s1, phy)

st.s2 <- asr.stoch(lik2, fit\$par)

What I do not understand is how to simulate stochastic mappings on my trees randomly – I want to

use the same rate parameters, but with no constraints on the tips. These simulations will help me

Determine whether the statistic I am examining is significantly large for my data.

In addition, I wonder if anyone has examined the distribution of two characters -

i.e. how can I calculate the percentage of time spent in (0,0), (0,1), etc. for the

classes generated by Diversitree.

Thanks a lot for your help!

software error genome • 1.8k views
modified 5.8 years ago • written 5.9 years ago by NK0
1
Asaf7.6k wrote:

I used the tool discrete of BayesTraits with good success for this kind of problem.

0
NK0 wrote:

Thanks Asaf, I am also using Bayestraits, but I would like to compare its results to another tool/method.