Re: [R] TukeyHSD for multiple response
Michael, Thank you very much for your answer. I finally tried lsmeans to compare what I wanted. I'll follow your advice and explore the CDA. It's probably a better solution to assess what I want. Best, Sérgio. - Mensagem original - > De: "Michael Friendly" > Para: "Sergio Ferreira Cardoso" , > "R-help list" > Enviadas: Sábado, 26 De Maio de 2018 17:28:20 > Assunto: Re: TukeyHSD for multiple response > Hi Sergio > > Doing those tests 30 times is going to give you a huge Type I error > rate, even if there was a function that did that. There is a reason > why TukeyHSD doesn't make it easy. > > In general, if there are useful comparisons among the species, you are > better off setting up and testing contrasts than doing all-pairwise > Tukey tests. > > Also, the PCA scores are ordered in terms of variance acct'd for, so > maybe only the first few are important. > > Finally, you might be better off using Canonical Discriminant analysis > than PCA followed by MANOVA. The candisc package is well suited to this > task. It can give you HE plots in the space that best discriminates > among the levels of an effect, and show how the original variables > relate to (project into) that space. > CDA is sort of like PCA, but the goal is to account for maximum > differences among groups rather than maximum total variance. > > For proper partial Type III tests, use car::Manova rather than stats::manova > which only gives sequential, Type I tests > > HTH > -Michael > > On 5/25/2018 9:11 AM, Sergio Ferreira Cardoso wrote: >> Dear all, >> >> I'm testing the effect of species and sex in my sample by using the principal >> component scores of a PCA analysis. >> I have 30 PCs and I tried to see if there is any significant difference from >> males to females, given that there is a significant effect of phylogeny >> (factor >> with several species). >> I didi it like this: >> >> Y<-PCA$pc.scores[,1:30] >> fit <- manova(Y ~ sp*sex) >> summary(fit, test="Wilks") >> >> I get a barely significant p-value for the effect of sex and I'd like to know >> for which of the species there is a difference between males and females. >> I tried TukeyHSD(fit) but I get the following error: >> >> Error in model.tables.aov(x, "means") : >> 'model.tables' is not implemented for multiple responses >> >> So this has to do with the fact that I have a multivariate independent >> variable. >> Is there an alternative function to this? >> >> Thanks in advance, >> Sérgio. >> > > > -- > Michael Friendly Email: friendly AT yorku DOT ca > Professor, Psychology Dept. & Chair, ASA Statistical Graphics Section > York University Voice: 416 736-2100 x66249 Fax: 416 736-5814 > 4700 Keele StreetWeb: http://www.datavis.ca > Toronto, ONT M3J 1P3 CANADA __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] TukeyHSD for multiple response
Dear all, I'm testing the effect of species and sex in my sample by using the principal component scores of a PCA analysis. I have 30 PCs and I tried to see if there is any significant difference from males to females, given that there is a significant effect of phylogeny (factor with several species). I didi it like this: Y<-PCA$pc.scores[,1:30] fit <- manova(Y ~ sp*sex) summary(fit, test="Wilks") I get a barely significant p-value for the effect of sex and I'd like to know for which of the species there is a difference between males and females. I tried TukeyHSD(fit) but I get the following error: Error in model.tables.aov(x, "means") : 'model.tables' is not implemented for multiple responses So this has to do with the fact that I have a multivariate independent variable. Is there an alternative function to this? Thanks in advance, Sérgio. -- Institut des Sciences de l'Evolution UMR5554, CNRS, IRD, EPHE Université de Montpellier Place Eugène Bataillon 34095 Montpellier Cedex 05 France Email: sergio.ferreira-card...@umontpellier.fr Tel: +33 (4 ) 67 14 46 52 [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] phylo.pca
Dear all, I'm trying to use phylo.pca function from phytools for the first time. I'm using an ultrametric tree with 167 tips and all branch lengths transformed to 1. I lunched the pPCA like this: pPCA<-phyl.pca(tree,data,method = "lambda") For some reason it takes for ever and never reaches the end of the process. I tried with method="BM" and the process runs normally. Has this happened to anyone here? How can I go around this? Thanks in advance, Sérgio. -- Institut des Sciences de l'Evolution UMR5554, CNRS, IRD, EPHE Université de Montpellier Place Eugène Bataillon 34095 Montpellier Cedex 05 France Email: sergio.ferreira-card...@umontpellier.fr Tel: +33 (4 ) 67 14 46 52 [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Mean correlation within cluster
Hello all, I'd like to calculate the mean correlation within a cluster and understand if it's significantly >0. I'm using packages 'geomorph' and 'paleomorph'. #Simulate an array A <- array ( rep ( 1 : 36 , by = 4 ), dim = c ( 12 , 3 , 4 )) #Load 'geomorph' package and superimpose coordinates test.gpa <- gpagen ( A , print.progress = FALSE ) #Load 'paleomorph' and generate covariance and correlation matrices cvmatrix <- dotcvm ( test.gpa $ coords ) corrmatrix <- dotcorr ( test.gpa $ coords ) Then I do a clustering with Ward method and euclidean distance, using the cvmatrix and I get a dendrogram. This part is not the problem, so I'll go directly to what I want. I would like to calculate the mean correlation between the elements of each cluster. Since clustering methods will mandatorily produce clusters, I'd like to know if the elements of my clusters are correlated (I mean, if the clusters are valid). I believe this might not be very complicated, given that I have all the elements. I just don't know how to do it on R. I tried to do the clustering with p-values for each cluster, with pvclust(), but I'm following a paper in which the authors test the significancy of the clusters by assessing the mean correlation of the elements within each cluster. I'd like to compare this method with the one from pvclust(). Thank you in advance. Best regards, Sérgio. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Chi-square test
Dear David and John, Thank you for your replies. Indeed I'm using ape and nlme packages. Here it is: > fit<-gls(fcl~mass+activity+agility,correlation=corBrownian(phy=tree),data=df,method="ML",weights=varFixed(~vf)) > Anova(fit) Analysis of Deviance Table (Type II tests) Response: fcl Df Chisq Pr(>Chisq) mass 1 0.1756 0.6752 activity 2 0.5549 0.7577 agility 4 3.2903 0.5105 Anyway, I have the help I was looking for. Thank you vey much. Best regards, Sérgio. - Mensagem original - > De: "Fox, John" > Para: "Sergio Ferreira Cardoso" > Cc: "R-help list" > Enviadas: Sábado, 21 De Janeiro de 2017 6:09:22 > Assunto: Re: [R] Chi-square test > Dear Sergio, > > You appear to have asked this question twice on r-help. > > Anova() has no specific method for “gls” models (I assume, though you don’t > say > so, that the model is fit by gls() in the nlme package), but the default > method > works and provides Wald chi-square tests for terms in the model. I don’t > understand the model formula x ~ 1 + 2 + 3 + x, however, and so I have no idea > what gls() would do with this model, other than report an error. Perhaps you > can show us the output — or, better yet, provide a reproducible example. > > As a general matter, for 1-df terms in an additive model, the 1-df chi-square > values reported by Anova() will simply be the squares of the corresponding > Wald > statistics (labelled “t” I believe) reported in the summary of the model. > Although the p-value is from the upper tail of the chi-square distribution, > the > test is inherently two-sided. > > Best, > John > > - > John Fox, Professor > McMaster University > Hamilton, Ontario, Canada > Web: http::/socserv.mcmaster.ca/jfox > >> On Jan 20, 2017, at 8:36 AM, Sergio Ferreira Cardoso >> wrote: >> >> Dear all, >> >> Anova() for .car package retrieves Chi-square statistics when I'm testing a >> model the significance of a multivariate .gls model >> gls(x~1+2+3+x,corBrownian(phy=tree), ...). >> Is this Chi-square a two-sided test? >> >> Thank you. >> >> Best, >> Sérgio. >> >> [[alternative HTML version deleted]] >> >> __ >> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see >> https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > > and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Chi-square test
Dear all, Anova() for .car package retrieves Chi-square statistics when I'm testing a model the significance of a multivariate .gls model gls(x~1+2+3+x,corBrownian(phy=tree), ...). Is this Chi-square a two-sided test? Thank you. Best, Sérgio. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Anova() Chi-square test
Dear all, Anova() for .car package retrieves Chi-square statistics when I'm testing a model the significance of a multivariate .gls model gls(x~1+2+3+x,corBrownian(phy=tree), ...). Is this Chi-square a two-sided test? Thank you. Best, Sérgio. -- Com os melhores cumprimentos, Sérgio Ferreira Cardoso. Best regards, Sérgio Ferreira Cardoso MSc. Paleontology Museu da Lourinhã, GEAL. LATR/IST/CTN - Campus Tecnológico e Nuclear. Lisboa, Portugal [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.