Hi Liam and all users
My morphological data varies from binary and multisate characters,
including polymorphism and non-applicable coding for many of those
characters. I've made a lot of progress after Liam's advices, but lately I
have been getting an error from make.simmap.
By trying to run a simulation with make.simmap under ARD model on
characters with states that the probability is zero on occuring in a given
terminal, after calling summary(), this error is displayed:
Error in colMeans(x$count) : 'x' must be an array of at least two dimensions
When changing the input of the same data in a two column array, the error
still persists. The same data was analyzed with make.simmap under "ER" and
"SYM" models, and they run just fine. I am also comparing the results using
fitDiscrete (geiger) and ace (ape) functions and they display results for
"ARD" model without a glitch. According to output from fitDiscrete
function, ARD model is indicated as the "best" model comparing AIC
likelihoods.
I've noticed that in this particular character, the "Q" matrix extracted
from make.simmap simmulation under "ARD" is the following:
Q =
0 1
0 -1.916782 1.916782
1 0.000000 0.000000
Are those "zero" probabilities the cause of the error? And by that I mean,
should I take from this evidence that "ARD" is not an appropriate model to
this character?
Also, results from characters with "non-applicable" coding (NA in matrix)
for some terminal are coming out very different from what I would expect.
Is this coding allowed for ancestral reconstruction?
I am sorry for the many questions, and I would be very glad if any
suggestions come to mind.
Best wishes,
Carolina Vieira.
---------------------------------------------
This is the error rundown:
> anc.tree_ard<-make.simmap(tree,ch434,model="ARD")
make.simmap is sampling character histories conditioned on the transition
matrix
Q =
0 1
0 -1.916782 1.916782
1 0.000000 0.000000
(estimated using likelihood);
and (mean) root node prior probabilities
pi =
0 1
0.5 0.5
Done.
> plot(anc.tree_ard,node.numbers=T,cores,fsize=0.5,ftype="i")
> add.simmap.legend(colors=cores,prompt=FALSE,x=0.5*par()$usr[1],
+ y=-30*par()$usr[3],fsize=0.8)
> title(main="ACSR \"ARD\" model (make.simmap)")
> anc.tree_ard<-make.simmap(tree,ch434,model="ARD",nsim=1000)
make.simmap is sampling character histories conditioned on the transition
matrix
Q =
0 1
0 -1.916782 1.916782
1 0.000000 0.000000
(estimated using likelihood);
and (mean) root node prior probabilities
pi =
0 1
0.5 0.5
Done.
> anc.tree_ard
1000 phylogenetic trees with mapped discrete characters
> simmap_ard<-summary(anc.tree_ard)
> simmap_ard
1000 trees with a mapped discrete character with states:
29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46,
47, 48, 49, 50, 51, 52, 53, 54, 55
Error in colMeans(x$count) : 'x' must be an array of at least two
dimensions
Carolina Vieira
Bióloga - Mestre em Ecologia e Conservação
Doutoranda em Biologia Animal (PPGBAN - UFRGS)
<https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail>
Livre
de vírus. www.avast.com
<https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail>.
<#DAB4FAD8-2DD7-40BB-A1B8-4E2AA1F9FDF2>
Em qua., 26 de fev. de 2020 às 23:54, Liam J. Revell <liam.rev...@umb.edu>
escreveu:
> Dear Carolina.
>
> After running make.simmap with nsim=100 or more, just run summary( ) on
> your "multiSimmap" object and then plot the result. For instance:
>
> map.trees<-make.simmap(tree,x,model="SYM",nsim=200)
> map.summary<-summary(map.trees)
> print(map.summary)
> plot(map.summary,ftype="off",lwd=1,cex=c(0.7,0.4))
> map.summary$ace ## probabilities at nodes
>
> All the best, Liam
>
> Liam J. Revell
> Associate Professor, University of Massachusetts Boston
> Profesor Asistente, Universidad Católica de la Ssma Concepción
> web: http://faculty.umb.edu/liam.revell/, http://www.phytools.org
>
> Academic Director UMass Boston Chile Abroad:
> https://www.umb.edu/academics/caps/international/biology_chile
>
> On 2/24/2020 2:33 PM, Carolina Santos Vieira wrote:
> > **[EXTERNAL SENDER]
> >
> > Hi Liam,
> >
> > Thank you very much for clearing that up!
> > I was able to run the make.simmap script, but in the first plot instead
> > of getting one symbol of the corresponding state at the tips, they are
> > displayed in proportional sizes (I've attached the plot).
> >
> > I have also looked into how to get the likelihoods from each ancestral
> > state from this analysis (like we get from "lik.anc" in ace function),
> > but didn't find it. Is it possible to get it?
> >
> > Thank you again for the help!
> > Best regards,
> > Carolina Vieira
> > Bióloga - Mestre em Ecologia e Conservação
> > Doutoranda em Biologia Animal (PPGBAN - UFRGS)
> >
> >
> > Em sex., 21 de fev. de 2020 às 22:53, Liam J. Revell
> > <liam.rev...@umb.edu <mailto:liam.rev...@umb.edu>> escreveu:
> >
> > Dear Carolina.
> >
> > I don't know how rayDISC works (but the authors of corHMM are on this
> > list, so they will probably respond); however, for
> > phytools::make.simmap
> > & phytools::rerootingMethod, which can also handle missing data,
> > missingness is coded by supplying the input data in the form of an N
> > x k
> > matrix for N species & k character states. In this case, we would
> > put 0s
> > and 1s in any row in which the character state was known; and a
> series
> > of values of 1/k for any completely unknown state.
> >
> > For instance:
> > 0 1 2
> > sp A 0 0 1
> > sp B 1 0 0
> > sp C 1/3 1/3 1/3
> > sp D 0 1/2 1/2
> >
> > corresponds to species A in state '2', species B in state '0', and
> > species C unknown or missing. Species D is in either state '1' or
> '2',
> > but not in state '0'. (Instead of 1/k you can also set them all to
> > be 1.
> > The likelihood changes, but ancestral states are the same.)
> >
> > I hope this is helpful. All the best, Liam
> >
> > Liam J. Revell
> > Associate Professor, University of Massachusetts Boston
> > Profesor Asistente, Universidad Católica de la Ssma Concepción
> > web: http://faculty.umb.edu/liam.revell/, http://www.phytools.org
> > <
> https://nam01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.phytools.org%2F&data=02%7C01%7Cliam.revell%40umb.edu%7Cbeaccbb894534e0785d708d7b9607ccd%7Cb97188711ee94425953c1ace1373eb38%7C0%7C0%7C637181696380978269&sdata=%2Fialzss4jrVm8Rpt6jDUXwkdGqQZMMFA4N84HgVhWps%3D&reserved=0
> >
> >
> > Academic Director UMass Boston Chile Abroad:
> > https://www.umb.edu/academics/caps/international/biology_chile
> >
> > On 2/21/2020 8:22 PM, Carolina Santos Vieira wrote:
> > > [EXTERNAL SENDER]
> > >
> > > Dear all,
> > >
> > > I am working with a morphological matrix which contains multistate
> > > characters (including polymorphism) and some missing data entries
> > (coded as
> > > NA in data file for input). There are no unvariable characters in
> > coding.
> > > For the reconstruction of ancestral characters, I came across
> > with rayDISC
> > > function from corHMM package (since fitDiscrete doesn't seem to
> > account for
> > > it).
> > >
> > > However, when running the script two kinds of errors are
> displayed:
> > >
> > > - Calling the full data file with the rayDISC function, no matter
> > which
> > > character is selected for the analysis, it returns with this
> error:
> > >
> > > rayDISC(phy,data, ntraits=1, charnum=1, rate.mat=NULL,
> model=c("ER"),
> > > node.states=c("marginal"), state.recon=c("subsequently"),
> > > lewis.asc.bias=FALSE, p=NULL, root.p=NULL, ip=NULL, lb=0, ub=100,
> > verbose=TRUE,
> > > diagn=FALSE)
> > >
> > > $diagnostic
> > > [1] "Character 1 is invariant. Analysis stopped."
> > >
> > >
> > > - I have tried clipping the data creating a new object with just
> one
> > > of the characters, but a new error appears:
> > >
> > > Error in data[, 1] : incorrect number of dimensions
> > >
> > > I also looked into running make.simmap in parallel
> > >
> > (
> https://nam01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fblog.phytools.org%2F2017%2F11%2Frunning-makesimmap-in-parallel.html&data=02%7C01%7Cliam.revell%40umb.edu%7Ce210401fdd8c483d2bde08d7b735d3f4%7Cb97188711ee94425953c1ace1373eb38%7C0%7C1%7C637179314141242835&sdata=%2FHwtbpjRxREIS7T2PPIQZ%2B5iyRpVbGaHMKWX6e8%2Fb5U%3D&reserved=0
> > <
> https://nam01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fblog.phytools.org%2F2017%2F11%2Frunning-makesimmap-in-parallel.html&data=02%7C01%7Cliam.revell%40umb.edu%7Cbeaccbb894534e0785d708d7b9607ccd%7Cb97188711ee94425953c1ace1373eb38%7C0%7C0%7C637181696380988260&sdata=Xbyv%2F2R6lsf7ArCi%2F5SmiV2r2nE57ZMsAsYrsL5uhhI%3D&reserved=0
> >),
> > > and this happens: " Error in if (any(x < 0)) x <- x - min(x) :
> > > missing Value when TRUE/FALSE needed ". Trying to replace NA as
> > > 'unknown' to be considered in the analysis by changing its
> > "status" to
> > > FALSE, it returns in the matrix coded as "0", which is not
> properly
> > > valid since this implies coding it as a known state (and not
> missing
> > > data, as it truly is).
> > >
> > >
> > > I appreciate if anyone could shine a light on this matter.
> > > With regards,
> > > Carolina Vieira
> > > Bióloga - Mestre em Ecologia e Conservação
> > > Doutoranda em Biologia Animal (PPGBAN - UFRGS)
> > >
> > > [[alternative HTML version deleted]]
> > >
> > > _______________________________________________
> > > R-sig-phylo mailing list - R-sig-phylo@r-project.org
> > <mailto:R-sig-phylo@r-project.org>
> > >
> >
> https://nam01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fstat.ethz.ch%2Fmailman%2Flistinfo%2Fr-sig-phylo&data=02%7C01%7Cliam.revell%40umb.edu%7Ce210401fdd8c483d2bde08d7b735d3f4%7Cb97188711ee94425953c1ace1373eb38%7C0%7C1%7C637179314141252833&sdata=feVEKUpekTBV5FVzfZDdDU7sJGi6x%2FGtumn3TpmhRuk%3D&reserved=0
> > <
> https://nam01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fstat.ethz.ch%2Fmailman%2Flistinfo%2Fr-sig-phylo&data=02%7C01%7Cliam.revell%40umb.edu%7Cbeaccbb894534e0785d708d7b9607ccd%7Cb97188711ee94425953c1ace1373eb38%7C0%7C0%7C637181696380988260&sdata=3OcfvhEJpXJEhCeN8M4LnvU1Kn385TZUY7mewlTtnB0%3D&reserved=0
> >
> > > Searchable archive at
> >
> https://nam01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.mail-archive.com%2Fr-sig-phylo%40r-project.org%2F&data=02%7C01%7Cliam.revell%40umb.edu%7Ce210401fdd8c483d2bde08d7b735d3f4%7Cb97188711ee94425953c1ace1373eb38%7C0%7C1%7C637179314141252833&sdata=8otz0Vyk02Uo2uGw4fWbOqfuA4EKp8GQ71cZufu3XFg%3D&reserved=0
> > <
> https://nam01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.mail-archive.com%2Fr-sig-phylo%40r-project.org%2F&data=02%7C01%7Cliam.revell%40umb.edu%7Cbeaccbb894534e0785d708d7b9607ccd%7Cb97188711ee94425953c1ace1373eb38%7C0%7C0%7C637181696380998258&sdata=RhuYW6qY%2ByiZBEBZT6J0buhER4Ks9s74JFKr%2BV%2Bnsfw%3D&reserved=0
> >
> > >
> >
>
<https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail>
Livre
de vírus. www.avast.com
<https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail>.
<#DAB4FAD8-2DD7-40BB-A1B8-4E2AA1F9FDF2>
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/
