Re: [R-sig-phylo] estimating mutation rates
Thanks, guys. I just have a follow-up question: Does anyone know of a way to optimize BLs using optim.pml() while enforcing a molecular clock? Thanks! John On Thu, Nov 8, 2012 at 6:06 PM, Rob Lanfear rob.lanf...@gmail.com wrote: Hi John, Brian's suggestion is a good one. Another option would be to use any of a number of molecular dating programs (not strictly an R response, sorry) like R8S (for ML) or BEAST (for Bayesian). You can then specify your node dates and topology as constraints, and use the molecular data appropriately. The advantage of these approaches is that you'll have full control over the models of rate change, and over the treatment of uncertainty in the dates and rates (depending on which method you choose). However, if you really do believe your dates without question, and if your molecular dataset has sufficient power for rate estimates (i.e. at least a handful of substitutions estimated on each branch), I would just do as Brian suggested. Rob On 9 November 2012 06:35, Brian O'Meara bome...@utk.edu wrote: If you want the branch-specific substitution rate, you could estimate the branch lengths on a topology constrained to match the chronogram, then divide each of the molecular tree's branch lengths by the corresponding chronogram branch lengths to get these rates. If you want a global rate, you might constrain the molecular tree to have the same proportional branch lengths as the chronogram. Look at pml() in phangorn for calculating the likelihoods. Best, Brian ___ Brian O'Meara Assistant Professor Dept. of Ecology Evolutionary Biology U. of Tennessee, Knoxville http://www.brianomeara.info Students wanted: Applications due Dec. 15, annually Postdoc collaborators wanted: Check NIMBioS' website Calendar: http://www.brianomeara.info/calendars/omeara On Thu, Nov 8, 2012 at 12:49 PM, john d dobzhan...@gmail.com wrote: Hi there, I have a well-supported phylogeny (with branch lengths proportional to absolute time) and I'm interested in estimating the mutation rate of an independently-obtained molecular dataset. Any suggestions? thanks! John ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo [[alternative HTML version deleted]] ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo -- Rob Lanfear Research Fellow, Ecology, Evolution, and Genetics, Research School of Biology, Australian National University www.robertlanfear.com July-September National Evolutionary Synthesis Center, 2024 W. Main Street, Suite A200, Durham, NC 27705-4667, USA phone: +1 919 668 4592 ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Re: [R-sig-phylo] prop.clades() errors
Hi Jeremy, Emmanuel and I had a off-list conversation after that email, so there was follow-up, we just never told the list about it! To my knowledge, the issue I was having where not all clades shared were being detected with prop.clades was unique to 3.0-1. For example, under my 3.0-6 install of ape, the code from my old email actually does show multiple shared clades with prop.clades. Here it is again with ape 3.0-6: library(ape) library(paleotree) set.seed(3) tree1-rtree(30) tree2-degradeTree(tree1,0.5) layout(matrix(1:2,,2));plot(tree1);plot(tree2) and prop.clades(tree1,tree2) [1] 1 2 0 0 0 1 1 0 0 1 0 0 1 1 0 0 1 1 1 0 1 0 0 1 1 0 1 0 0 Now, for the sake of posterity, it turns out the '2' which threw me in that original email was a result of me not telling prop.clades that my trees were rooted (and Emmanuel has updated the man page so this is crystal clear): prop.clades(tree1,tree2,rooted=TRUE) [1] 1 1 0 0 0 1 1 1 0 0 1 0 0 1 1 0 0 1 1 1 0 1 0 0 1 1 0 1 0 Jeremy, could you make some example data where prop.clades underestimates shared clades? Also, have you also tried applying reorder.phylo or collapse.singles to your trees before trying prop.clades? -Dave On Thu, Nov 8, 2012 at 2:04 PM, Jeremy Yoder jbyo...@gmail.com wrote: All, I'm working on some analyses that require comparison of tree structures, and I've bumped into the problem David Bapst posted about back in March: http://www.mail-archive.com/r-sig-phylo@r-project.org/msg01840.html - briefly, prop.clades() is failing to count many clades that are in common. There doesn't seem to be any followup to that post; has anyone made progress figuring it out? I'm currently using ape 3.0-6 and R 2.15.1 on a MacBook running Mountain Lion. thanks, Jeremy Jeremy B. Yoder Postdoctoral Associate Department of Plant Biology University of Minnesota www.jeremybyoder.com (http://www.jeremybyoder.com) jbyo...@gmail.com Sent with Sparrow (http://www.sparrowmailapp.com/?sig) [[alternative HTML version deleted]] ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo -- David Bapst Dept of Geophysical Sciences University of Chicago 5734 S. Ellis Chicago, IL 60637 http://home.uchicago.edu/~dwbapst/ http://cran.r-project.org/web/packages/paleotree/index.html http://home.uchicago.edu/%7Edwbapst/ [[alternative HTML version deleted]] ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Re: [R-sig-phylo] estimating mutation rates
Hi John, optim.pml has an parameter optRooted and you have to set it to TRUE. If you give the function an ultrametric tree it will optimise the edge lengths. Unfortunately it is not possible to use tree rearrangements with this option so far. Regards, Klaus On 11/9/12, john d dobzhan...@gmail.com wrote: Thanks, guys. I just have a follow-up question: Does anyone know of a way to optimize BLs using optim.pml() while enforcing a molecular clock? Thanks! John On Thu, Nov 8, 2012 at 6:06 PM, Rob Lanfear rob.lanf...@gmail.com wrote: Hi John, Brian's suggestion is a good one. Another option would be to use any of a number of molecular dating programs (not strictly an R response, sorry) like R8S (for ML) or BEAST (for Bayesian). You can then specify your node dates and topology as constraints, and use the molecular data appropriately. The advantage of these approaches is that you'll have full control over the models of rate change, and over the treatment of uncertainty in the dates and rates (depending on which method you choose). However, if you really do believe your dates without question, and if your molecular dataset has sufficient power for rate estimates (i.e. at least a handful of substitutions estimated on each branch), I would just do as Brian suggested. Rob On 9 November 2012 06:35, Brian O'Meara bome...@utk.edu wrote: If you want the branch-specific substitution rate, you could estimate the branch lengths on a topology constrained to match the chronogram, then divide each of the molecular tree's branch lengths by the corresponding chronogram branch lengths to get these rates. If you want a global rate, you might constrain the molecular tree to have the same proportional branch lengths as the chronogram. Look at pml() in phangorn for calculating the likelihoods. Best, Brian ___ Brian O'Meara Assistant Professor Dept. of Ecology Evolutionary Biology U. of Tennessee, Knoxville http://www.brianomeara.info Students wanted: Applications due Dec. 15, annually Postdoc collaborators wanted: Check NIMBioS' website Calendar: http://www.brianomeara.info/calendars/omeara On Thu, Nov 8, 2012 at 12:49 PM, john d dobzhan...@gmail.com wrote: Hi there, I have a well-supported phylogeny (with branch lengths proportional to absolute time) and I'm interested in estimating the mutation rate of an independently-obtained molecular dataset. Any suggestions? thanks! John ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo [[alternative HTML version deleted]] ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo -- Rob Lanfear Research Fellow, Ecology, Evolution, and Genetics, Research School of Biology, Australian National University www.robertlanfear.com July-September National Evolutionary Synthesis Center, 2024 W. Main Street, Suite A200, Durham, NC 27705-4667, USA phone: +1 919 668 4592 ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo -- Klaus Schliep Phylogenomics Lab at the University of Vigo, Spain ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Re: [R-sig-phylo] prop.clades() errors
Thanks for getting back, Dave! Klaus Schliep had already replied off-list â he suggested that I make sure the tip labels in the trees under comparison are in the same order, using a feature in phangorn: tree2 - phangorn:::checkLabels(tree2, tree1$tip.label) prop.clades(tree1,tree2) And that seems to have fixed the issue for me. Weirdly, I was able to reproduce the error using the code from your original e-mail on my system (on ape 3.0-6); but now it's working correctly. I can break it again by deliberately mismatching the order of taxa in the two trees: library(ape) library(paleotree) set.seed(3) tree1 - rtree(30) tree2 - degradeTree(tree1,0.5) layout(matrix(1:2,,2));plot(tree1);plot(tree2) prop.clades(tree1,tree2,rooted=T) [1] 1 1 0 0 0 1 1 1 0 0 1 0 0 1 1 0 0 1 1 1 0 1 0 0 1 1 0 1 0 tree2 - phangorn:::checkLabels(tree2,tree1$tip.label[30:1]) prop.clades(tree1,tree2,rooted=T) [1] 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 So it looks to me like tip label ordering was my root problem. cheers, Jeremy Jeremy B. Yoder Postdoctoral Associate Department of Plant Biology University of Minnesota www.jeremybyoder.com (http://www.jeremybyoder.com) jbyo...@gmail.com Sent with Sparrow (http://www.sparrowmailapp.com/?sig) On Friday, 9 November, 2012 at 10:30 AM, David Bapst wrote: Hi Jeremy, Emmanuel and I had a off-list conversation after that email, so there was follow-up, we just never told the list about it! To my knowledge, the issue I was having where not all clades shared were being detected with prop.clades was unique to 3.0-1. For example, under my 3.0-6 install of ape, the code from my old email actually does show multiple shared clades with prop.clades. Here it is again with ape 3.0-6: library(ape) library(paleotree) set.seed(3) tree1-rtree(30) tree2-degradeTree(tree1,0.5) layout(matrix(1:2,,2));plot(tree1);plot(tree2) and prop.clades(tree1,tree2) [1] 1 2 0 0 0 1 1 0 0 1 0 0 1 1 0 0 1 1 1 0 1 0 0 1 1 0 1 0 0 Now, for the sake of posterity, it turns out the '2' which threw me in that original email was a result of me not telling prop.clades that my trees were rooted (and Emmanuel has updated the man page so this is crystal clear): prop.clades(tree1,tree2,rooted=TRUE) [1] 1 1 0 0 0 1 1 1 0 0 1 0 0 1 1 0 0 1 1 1 0 1 0 0 1 1 0 1 0 Jeremy, could you make some example data where prop.clades underestimates shared clades? Also, have you also tried applying reorder.phylo or collapse.singles to your trees before trying prop.clades? -Dave On Thu, Nov 8, 2012 at 2:04 PM, Jeremy Yoder jbyo...@gmail.com (mailto:jbyo...@gmail.com) wrote: All, I'm working on some analyses that require comparison of tree structures, and I've bumped into the problem David Bapst posted about back in March: http://www.mail-archive.com/r-sig-phylo@r-project.org/msg01840.html - briefly, prop.clades() is failing to count many clades that are in common. There doesn't seem to be any followup to that post; has anyone made progress figuring it out? I'm currently using ape 3.0-6 and R 2.15.1 on a MacBook running Mountain Lion. thanks, Jeremy Jeremy B. Yoder Postdoctoral Associate Department of Plant Biology University of Minnesota www.jeremybyoder.com (http://www.jeremybyoder.com) (http://www.jeremybyoder.com) jbyo...@gmail.com (mailto:jbyo...@gmail.com) Sent with Sparrow (http://www.sparrowmailapp.com/?sig) [[alternative HTML version deleted]] ___ R-sig-phylo mailing list R-sig-phylo@r-project.org (mailto:R-sig-phylo@r-project.org) https://stat.ethz.ch/mailman/listinfo/r-sig-phylo -- David Bapst Dept of Geophysical Sciences University of Chicago 5734 S. Ellis Chicago, IL 60637 http://home.uchicago.edu/~dwbapst/ (http://home.uchicago.edu/%7Edwbapst/) http://cran.r-project.org/web/packages/paleotree/index.html (http://home.uchicago.edu/%7Edwbapst/) [[alternative HTML version deleted]] ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Re: [R-sig-phylo] Ancestral Reconstruction Giving Identical Values
I would have to see your code to know 100% for certain - but your results make
intuitive sense, so let me just guess. When you have OU under a very high alpha
(which is what you are getting here), then you are modeling a process following
a random walk with an extremely high tendency to be pulled towards the OU
optimum value. In this case, characters do not retain any information about
their past states, and so the best guess for every node in the tree is the
mean of all of the characters. You can check to see if the value you are
getting back is actually the mean of all of your tips.
To me, what this means is that if you think the model is correct (lambda or OU)
then you have no information about any of the ancestral character states in the
tree.
lh
On Nov 8, 2012, at 3:14 PM, Zachary Chejanovski sliv...@life.bio.sunysb.edu
wrote:
Hi,
I am attempting to estimate ancestral states using the following methods-
APE::ace, GEIGER::getancstates, as well as this method recommend by Liam
Revell for reconstructing under an OU model:
ou-fitContinuous(tree,data,model=OU)$Trait1
anc.mins-vector()
N-length(tree$tip)
M-tree$Nnode
outree-ouTree(tree,ou$alpha)
for(i in 1:M+N){
anc.mins[i-N]-ace(data,multi2di(root(outree,node=i)),method=pic)$ace[1]
names(anc.mins)[i-N]-i
}
Here are the parameter estimates for one of my traits:
$Lambda
lnl betalambda aic aicc k
1 -121.4771 31.33674 0.4057359 248.9541 250.2173 3
$OU
lnl beta alpha convergence
message k aic aicc
1 -120.7098 281.9509 0.08288029 0 CONVERGENCE: REL_REDUCTION_OF_F
= FACTR*EPSMCH 3 247.4196 248.6828
When I do not transform the tree (i.e. Brownian Motion), the values returned
seem valid. But when I transform the tree using either lambdaTree, ouTree. or
even the exponentialTree, all the values returned are almost identical and
when I plot them on the phylogeny, it is evident that they do not make sense
in regards to the tip values.
Any help on this would be greatly appreciated.
Thank you,
Zac
[[alternative HTML version deleted]]
___
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Luke Harmon
Associate Professor
Biological Sciences
University of Idaho
208-885-0346
lu...@uidaho.edu
___
R-sig-phylo mailing list
R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
[R-sig-phylo] diversification analysis questions
Hi all I write to ask a simple question, anyone know how to plot the parameters obtained by applying the different diversification models implemented in the program laser of Rabosky ??. if anyone can guide me would be very nice cheers! -- Oscar Inostroza Michael Biólogo Laboratorio de Ecología Evolutiva y Filoinformática Facultad de Ciencias Naturales Y Oceanográficas Universidad de Concepción Barrio Universitario s/n, Casilla 160-C Concepción, CHILE Fono Oficina: (56) 41 - 220 – 7325 e-mail: oskinostr...@gmail.com oscarinostr...@udec.cl ___ R-sig-phylo mailing list R-sig-phylo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
