[R-sig-phylo] axisChrono and axisPhylo
Greetings R-sig-phylo, This is my first post but happy to join this community I have been having some trouble with the axisChrono of the phyloch package and axisPhylo functions of the ape package that I was hoping someone might be able to help me with. When using either function, the spacing between ticks on the scale is too large for the figure I’m trying to generate. I have tried passing variations of at=c(1,2,3,4,5), or at=seq(0,10,1) etc, but for either function I get the following error: Error in axis(side = side, at = c(maxi - x), labels = abs(x * fact), ...) : formal argument “at matched by multiple actual arguments Does anyone know how I might be able to specify the tick marks that I want for these particular functions? Best, Jacob Berv Ph.D. Student Lovette Lab Cornell Laboratory of Ornithology ___ 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/
Re: [R-sig-phylo] phytools - evaluating significance of pgls.Ives
Thank you all very much for the comments! They are so helpful. Yes, I do only have 8 species, and 3 replicates each. It is not ideal, but it's what we have and we have a phylogeny, so I'd like to try some tests incorporating phylogeny. I probably should have added that I plan on running OLS to test each regression without the tree. This should give me an idea of the general relationship. Regarding Liam's comment in pgls.Ives: Is the lower=c(1e-8,1e-8) call of the pgls.Ives call the part where it constrains the slope to (almost)zero? I'll definitely give the MERegPHYSIGv2.m method a try as well, with Tony's diagnostic suggestions. I will also go the LRT route with the data. Thanks very much again for your help! Andrea ~~ Andrea Berardi, PhD Postdoctoral Researcher, Smith Lab EBIO, University of Colorado-Boulder andrea.bera...@colorado.edu On Mar 1, 2015, at 8:42 PM, Liam J. Revell liam.rev...@umb.edu wrote: Hi Andrea. This is not presently implemented, but since this is a likelihood method it would be straightforward to constrain to a slope of zero and then do a LR test. This would be probably be the easiest way to test a hypothesis about the regression. That being said, as noted in the function documentation, some problems have been reported with the optimization algorithm for this model, which is simple and thus may fail to find the ML solution. Consequently, I would encourage you to look for other implementations of the method so that you can be confident in your result. I'm not aware of one in R at this time. All the best, Liam Liam J. Revell, Assistant Professor of Biology University of Massachusetts Boston web: http://faculty.umb.edu/liam.revell/ email: liam.rev...@umb.edu blog: http://blog.phytools.org On 3/1/2015 10:31 PM, Andrea Berardi wrote: Hi all, I'm just learning how to do PGLS analyses, and I'm looking for advice on how to evaluate the significance of the regression fit using pgls.Ives in the phytools package. I'm using this function because it incorporates sampling error of species means, and my data has about 3 individuals per species, with 8 species. My goal is to test whether a flower trait predicts the leaf trait, while controlling for shared ancestry. Here is the output from pgls.Ives: fit - pgls.Ives(Tree, Flower_trait, Leaf_trait) fit $beta [1] 96.3963098 0.1292656 $sig2x [1] 22218901073 $sig2y [1] 23027587 $a [1] -10063.150 -1204.422 $logL [1] -158.2337 $convergence [1] 0 $message [1] CONVERGENCE: REL_REDUCTION_OF_F = FACTR*EPSMCH I am also running pgls on species averages for the traits using the gls function in nlme and the corBrownian and corMartins functions in ape. But, we are interested in incorporating the within-species variation in our small dataset. Any suggestions would be welcome! Thanks for your help, Andrea ~~ Andrea Berardi, PhD Postdoctoral Researcher, Smith Lab EBIO, University of Colorado-Boulder ___ 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/ ___ 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/
Re: [R-sig-phylo] phytools - evaluating significance of pgls.Ives
Andrea, remember that you can and should also do the OLS models (i.e., assuming a star phylogeny) with measurement error considered. That's in the programs that accompany Ives, Midford, and Garland (2007, Syst. Biol. 56:252–270), and were in the batch I just sent you. Cheers, Ted Theodore Garland, Jr., Professor Department of Biology University of California, Riverside Riverside, CA 92521 Office Phone: (951) 827-3524 Facsimile: (951) 827-4286 (not confidential) Email: tgarl...@ucr.edu http://www.biology.ucr.edu/people/faculty/Garland.html http://scholar.google.com/citations?hl=enuser=iSSbrhwJ Director, UCR Institute for the Development of Educational Applications Editor in Chief, Physiological and Biochemical Zoology Fail Lab: Episode One http://testtube.com/faillab/zoochosis-episode-one-evolution http://www.youtube.com/watch?v=c0msBWyTzU0 From: R-sig-phylo [r-sig-phylo-boun...@r-project.org] on behalf of Andrea Berardi [andrea.bera...@colorado.edu] Sent: Monday, March 02, 2015 3:57 PM To: r-sig-phylo@r-project.org Cc: Anthony R Ives; Peter Smits Subject: Re: [R-sig-phylo] phytools - evaluating significance of pgls.Ives Thank you all very much for the comments! They are so helpful. Yes, I do only have 8 species, and 3 replicates each. It is not ideal, but it's what we have and we have a phylogeny, so I'd like to try some tests incorporating phylogeny. I probably should have added that I plan on running OLS to test each regression without the tree. This should give me an idea of the general relationship. Regarding Liam's comment in pgls.Ives: Is the lower=c(1e-8,1e-8) call of the pgls.Ives call the part where it constrains the slope to (almost)zero? I'll definitely give the MERegPHYSIGv2.m method a try as well, with Tony's diagnostic suggestions. I will also go the LRT route with the data. Thanks very much again for your help! Andrea ~~ Andrea Berardi, PhD Postdoctoral Researcher, Smith Lab EBIO, University of Colorado-Boulder andrea.bera...@colorado.edu On Mar 1, 2015, at 8:42 PM, Liam J. Revell liam.rev...@umb.edu wrote: Hi Andrea. This is not presently implemented, but since this is a likelihood method it would be straightforward to constrain to a slope of zero and then do a LR test. This would be probably be the easiest way to test a hypothesis about the regression. That being said, as noted in the function documentation, some problems have been reported with the optimization algorithm for this model, which is simple and thus may fail to find the ML solution. Consequently, I would encourage you to look for other implementations of the method so that you can be confident in your result. I'm not aware of one in R at this time. All the best, Liam Liam J. Revell, Assistant Professor of Biology University of Massachusetts Boston web: http://faculty.umb.edu/liam.revell/ email: liam.rev...@umb.edu blog: http://blog.phytools.org On 3/1/2015 10:31 PM, Andrea Berardi wrote: Hi all, I'm just learning how to do PGLS analyses, and I'm looking for advice on how to evaluate the significance of the regression fit using pgls.Ives in the phytools package. I'm using this function because it incorporates sampling error of species means, and my data has about 3 individuals per species, with 8 species. My goal is to test whether a flower trait predicts the leaf trait, while controlling for shared ancestry. Here is the output from pgls.Ives: fit - pgls.Ives(Tree, Flower_trait, Leaf_trait) fit $beta [1] 96.3963098 0.1292656 $sig2x [1] 22218901073 $sig2y [1] 23027587 $a [1] -10063.150 -1204.422 $logL [1] -158.2337 $convergence [1] 0 $message [1] CONVERGENCE: REL_REDUCTION_OF_F = FACTR*EPSMCH I am also running pgls on species averages for the traits using the gls function in nlme and the corBrownian and corMartins functions in ape. But, we are interested in incorporating the within-species variation in our small dataset. Any suggestions would be welcome! Thanks for your help, Andrea ~~ Andrea Berardi, PhD Postdoctoral Researcher, Smith Lab EBIO, University of Colorado-Boulder ___ 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/ ___ 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/ ___ 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/
Re: [R-sig-phylo] axisChrono and axisPhylo
Thanks Franz,
I tried modifying the axisChrono function, but when I change the “at=“ argument
to equal something like “at= c(1,2,3,4,5)” I can get it to list out specific
tick marks but they are always in increasing order, not reverse order as you
want with a right facing cladogram. I’m not sure how to reverse the sequence
within the function.
Jake
Here is the function.
function (side = 1, unit = Ma, fact = 1, ...)
{
lastPP - get(last_plot.phylo, envir = .PlotPhyloEnv)
if (lastPP$type %in% c(phylogram, cladogram)) {
if (lastPP$direction %in% c(rightwards, leftwards)) {
x - pretty(lastPP$xx)
if (lastPP$direction == rightwards)
maxi - max(lastPP$xx)
else {
maxi - min(lastPP$xx)
x - -x
}
}
else {
x - pretty(lastPP$yy)
if (lastPP$direction == upwards)
maxi - max(lastPP$yy)
else {
maxi - min(lastPP$yy)
x - -x
}
}
}
axis(side = side, at = c(maxi - x), labels = abs(x * fact),
...)
mtext(text = unit, side = side, at = 1.07 * maxi, line = 1,
...)
}
On Mar 2, 2015, at 5:03 PM, Franz Krah f.k...@mailbox.org wrote:
Hi Jacob,
as far as I can see 'axisPhylo' calculates the at argument of 'axis'
itself...
I think that produces the error if you just put the argument into the
'axisPhylo' function.
So you have to alter the 'axisPhylo' function so you can define your own
'axis' at arguments.
Hope that helps and is correct...
Best Franz
p.s. that's my first answer here but I feel today is a god day to start
(thanks Dave Bapst for encouraging ;-))
Jacob Berv jakeberv.r.sig.ph...@gmail.com hat am 2. März 2015 um 22:38
geschrieben:
Greetings R-sig-phylo,
This is my first post but happy to join this community
I have been having some trouble with the axisChrono of the phyloch package
and
axisPhylo functions of the ape package that I was hoping someone might be
able
to help me with.
When using either function, the spacing between ticks on the scale is too
large for the figure I’m trying to generate. I have tried passing variations
of at=c(1,2,3,4,5), or at=seq(0,10,1) etc, but for either function I get the
following error:
Error in axis(side = side, at = c(maxi - x), labels = abs(x * fact), ...) :
formal argument “at matched by multiple actual arguments
Does anyone know how I might be able to specify the tick marks that I want
for
these particular functions?
Best,
Jacob Berv
Ph.D. Student
Lovette Lab
Cornell Laboratory of Ornithology
___
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/
MSc. Franz Krah
Mycology, Phylogenetics, PCM
http://www.biodiv.wzw.tum.de/index.php?id=18
Jacob Berv
Ph.D. Student
Lovette Lab
Cornell Laboratory of Ornithology
___
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/
Re: [R-sig-phylo] The Value of Responding on R-sig-phylo, was: phytools - evaluating significance of pgls.Ives
Hi Everyone, I agree with what has been said so far. Tony, don't take my comments to heart. Your contribution is most welcome as your opinion is highly respected (by me and hopefully by others). Same goes for all the other authorities on the list e.g Liam, Joe, Ben, Jarrod, Ted, Marguerite and others. I don't contribute much these days due to time constraints, which is a shame because I think that lists like this do help people, even if particular answers to particular questions are not quite right or are incomplete in some sense for the original poster. And I like helping people too! Questions can get answers that evolve into interesting threads in their own right (like this one, although a bit off-topic). Mailing lists like this are more of a group conversation rather than a simple do this or don't do that type of thing. I've been on lists since the '90s (does anyone remember sci.bio.evolution on USENET?). This is one of the best lists I've been on, largely because I know that some of the best researchers in the field read the list and sometimes respond. Tony, I know that statistical advice is easier to give when the person and data are live right in front of you. That's why I invariably only see people by appointment when I give stats advice. I don't like giving advice over email and I don't give advice over the phone. (Actually, I never answer my phone but that's another story.) It's a bit different on the list because other authorities can read the advice an chime in if they disagree with anything I say. Hopefully we can all benefit from that. Keep on posting everybody! Simon. On 03/03/15 03:43, Hilmar Lapp wrote: +100 !! -hilmar On 3/2/15, 11:07 AM, David Bapst wrote: Off-topic, but I wanted to comment on Anthony's wariness about commenting on R-sig lists... Yes, it can certainly be very difficult to give advice that is tailored specifically to the problem of a particular worker on R-sig lists, particularly as one can't just tell the other person to open their data files and show them to you. However, answering a question on R-sig lists (regardless of whether it is the right answer, but what answer remains right forever anyway?) records that answer virtually forever, for future posterity. This means the answer can be found by any future individuals who encounter this problem, via a simple google search. Since I first found R-sig-phylo in early 2010 (five years ago!), I have made almost constant use of the knowledge base contained with R-sig-phylo's archive. It doesn't mean the information is always right, or always the best answers for the actual person who asked originally, but it certainly helps point out the right line of thinking or the right literature to investigate. StackExchange discussion archives have become just as valuable (but there appears to be very little R phylo discussion over there). It feels like the amount of discussion on the list has dropped off slowly over the last year, and as phylogenetics in R seems as widely used as ever, I have often wondered if this reflects that most issues people may run into are now more easily solved with google searches leading to the R-sig-phylo archives. Anyway, I just hope that this element is not forgotten about why I think replying to question on R-sig-phylo remains a valuable contribution to our community! Cheers, -Dave Bapst On Mon, Mar 2, 2015 at 7:58 AM, Anthony Ives ari...@wisc.edu wrote: Simon and Ben, Of course, sample size of 8 is going to be an issue in almost any analysis. But sometimes that is all the data there are. Incidentally, this exchange reminded me that I’m still wary of making comments on r-sig. If somebody comes into my office, I have the time to discuss with them their data, so I can learn more about it. Then I feel I can at least make informed recommendations for analyses — they might still be badly wrong recommendations, but at least they are informed. I’m still uncomfortable about making suggestions on r-sig, when I don’t really have full information, or the time to think. Therefore, the few comments I’ve made have been very general about methods, rather than specific about data sets. I think this is just a matter of me waking up to the 21st century. I do like the idea of crowdsourcing; I just need to get comfortable with it. Cheers, Tony Anthony Ives Department of Zoology 459 Birge Hall (4th floor, E end of bldg) UW-Madison Madison, WI 53706 608-262-1519 On Mar 1, 2015, at 10:53 PM, Simon Blomberg s.blombe...@uq.edu.au wrote: Hi Ben, Yes, you would have to assume constant variance across species to use N=24. I think that is the only option. But given that biological data often has a positive mean-variance relationship, again I'm dubious about the exercise. YMMV, however! Cheers, Simon. Simon Blomberg, BSc (Hons), PhD, MAppStat, AStat. Senior Lecturer and Consultant Statistician School of Biological Sciences The University of Queensland St. Lucia Queensland
Re: [R-sig-phylo] axisChrono and axisPhylo
Hi Jacob, What I do is to use the basic axis() function. The trick is to transform the scale by subtracting the depth of your tree like this: axis(1, at=max(branching.times(tree))-0:10, labels=0:10) Best, Santiago On Mar 2, 2015, at 4:38 PM, Jacob Berv jakeberv.r.sig.ph...@gmail.com wrote: Greetings R-sig-phylo, This is my first post but happy to join this community I have been having some trouble with the axisChrono of the phyloch package and axisPhylo functions of the ape package that I was hoping someone might be able to help me with. When using either function, the spacing between ticks on the scale is too large for the figure I’m trying to generate. I have tried passing variations of at=c(1,2,3,4,5), or at=seq(0,10,1) etc, but for either function I get the following error: Error in axis(side = side, at = c(maxi - x), labels = abs(x * fact), ...) : formal argument “at matched by multiple actual arguments Does anyone know how I might be able to specify the tick marks that I want for these particular functions? Best, Jacob Berv Ph.D. Student Lovette Lab Cornell Laboratory of Ornithology ___ 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/ Santiago Claramunt, Ph.D. Research Associate Department of Ornithology American Museum of Natural History sclaram...@amnh.org https://sites.google.com/site/sclaramuntuy/ ___ 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/
Re: [R-sig-phylo] axisChrono and axisPhylo
Excellent! Thanks Santiago. In case anyone is interested, I ended up using the following three lines to generate an axis with labeled and thick tick marks every 10 and thin (unlabeled) ticks every 2. There is probably a more elegant way to do this but this works: axis(3, at=max(branching.times(timetree))-seq(0,70,10), labels=seq(0,70,10),cex.axis=0.5, line=-0.95, lwd=2) axis(3, at=max(branching.times(timetree))-seq(0,70,2),cex.axis=0.25, line=-0.95, lwd=0.75, labels=FALSE, tick=TRUE) mtext(text = Ma, side = 3, at = 1.04 * max(branching.times(timetree)), line = -1) On Mar 2, 2015, at 11:10 PM, Santiago Claramunt sclaram...@amnh.org wrote: axis(1, at=max(branching.times(tree))-0:10, labels=0:10) Jacob Berv Ph.D. Student Lovette Lab Cornell Laboratory of Ornithology [[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/
Re: [R-sig-phylo] phytools - evaluating significance of pgls.Ives
Simon and Ben, Of course, sample size of 8 is going to be an issue in almost any analysis. But sometimes that is all the data there are. Incidentally, this exchange reminded me that I’m still wary of making comments on r-sig. If somebody comes into my office, I have the time to discuss with them their data, so I can learn more about it. Then I feel I can at least make informed recommendations for analyses — they might still be badly wrong recommendations, but at least they are informed. I’m still uncomfortable about making suggestions on r-sig, when I don’t really have full information, or the time to think. Therefore, the few comments I’ve made have been very general about methods, rather than specific about data sets. I think this is just a matter of me waking up to the 21st century. I do like the idea of crowdsourcing; I just need to get comfortable with it. Cheers, Tony Anthony Ives Department of Zoology 459 Birge Hall (4th floor, E end of bldg) UW-Madison Madison, WI 53706 608-262-1519 On Mar 1, 2015, at 10:53 PM, Simon Blomberg s.blombe...@uq.edu.au wrote: Hi Ben, Yes, you would have to assume constant variance across species to use N=24. I think that is the only option. But given that biological data often has a positive mean-variance relationship, again I'm dubious about the exercise. YMMV, however! Cheers, Simon. Simon Blomberg, BSc (Hons), PhD, MAppStat, AStat. Senior Lecturer and Consultant Statistician School of Biological Sciences The University of Queensland St. Lucia Queensland 4072 Australia T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au http://www.evolutionarystatistics.org Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem Basically, I'm not interested in doing research and I never have been. I'm interested in understanding, which is quite a different thing. - David Blackwell From: R-sig-phylo [r-sig-phylo-boun...@r-project.org] on behalf of Ben Bolker [bbol...@gmail.com] Sent: Monday, March 02, 2015 2:49 PM To: r-sig-phylo@r-project.org Subject: Re: [R-sig-phylo] phytools - evaluating significance of pgls.Ives -BEGIN PGP SIGNED MESSAGE- Hash: SHA1 On 15-03-01 11:40 PM, Simon Blomberg wrote: Am I missing something? The OP only has 8 species in the data set. I wouldn't put much store in fancy PCM modelling based on such a small data set. And 3 individuals per species is not enough for a good estimate of the within-species variance. Simon. Agree wholeheartedly with the first point -- but for the second, isn't 24 rather than 8 the relevant number for estimating within-species variance (since presumably we are assuming the same variance within every species, thus we can effectively pool within-species variation \across species for this purpose) ? Simon Blomberg, BSc (Hons), PhD, MAppStat, AStat. Senior Lecturer and Consultant Statistician School of Biological Sciences The University of Queensland St. Lucia Queensland 4072 Australia T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au http://www.evolutionarystatistics.org Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem Basically, I'm not interested in doing research and I never have been. I'm interested in understanding, which is quite a different thing. - David Blackwell From: R-sig-phylo [r-sig-phylo-boun...@r-project.org] on behalf of Anthony R Ives [ari...@wisc.edu] Sent: Monday, March 02, 2015 2:14 PM To: Andrea Berardi Cc: r-sig-phylo@r-project.org Subject: Re: [R-sig-phylo] phytools - evaluating significance of pgls.Ives Andrea, I second Liam’s recommendation to use a LRT. For measurement error, the latest code I have in matlab is MERegPHYSIGv2.m, which does both measurement error and an OU or Pagel-lambda transform (see Johnson, M. T. J., A. R. Ives, J. Ahern, and J. P. Salminen. 2014. Macroevolution of plant defenses against herbivores in the evening primroses. New Phytologist 203:267-279). Measurement-error models are always going to have difficulties at parameter boundaries; for example, if the assumed measurement error is large, it can exceed the observed variation in the data, which of course causes problems (statistical and logical). In MERegPHYSIGv2.m, I did a round or two of simulated annealing first, before polishing the results with a Nelder-Mead optimizer. It seems like you could do the same with Liam’s code pretty easily by changing the method of optimization (using edit()). Before doing this, thought, I would take a careful look at your data and your estimates of measurement error. An easy diagnostic is to start with 10% of your estimated measurement standard errors and then increase slowly to 100%. When I have done this, I’ve been able to see problems when parameter values go awry. It is not a fail-safe
