Re: [R-br] Dúvida - elipse de confiança
Olá, Mauro Não sei se resolve no seu caso, mas dê uma olhada na função ordiellipse() do pacote vegan. Ela funciona bem para as ordenações feitas no próprio pacote, mas não sei se é compatível com outras ordenações, como no seu caso. Mas uma possibilidade é fazer a pca pelo vegan também, usando a função rda() com o argumento scale=T. Abraços Marcos -- Marcos Vinícius Carneiro Vital Universidade Federal de Alagoas Instituto de Ciências Biológicas e da Saúde Setor de Biodiversidade e Ecologia
Mauro, boa tarde! Trabalhei um pouco em cima do seu código e cheguei a plotar as elipses, mas precisaria revisar o conteúdo de multivariada pra saber se o procedimento está correto. De qualquer modo, posto o código para que você veja se é possível aproveitar algo. ### <code r> dURL <- (' https://www.dropbox.com/s/405ggda0h6scnie/marajoara320amostras.txt') dURL2 <- (' https://dl.dropboxusercontent.com/s/405ggda0h6scnie/marajoara320amostras.txt... ') dName <- basename(dURL); dName download.file(dURL2, dName, mode='wb') file.show(dName) maraOri <- read.table('marajoara320amostras.txt', head=FALSE); head(maraOri) maraNorm <- data.frame(scale(maraOri)); head(maraNorm) ### Cluster maraNorm.h <- hclust(dist(maraNorm)^2, method="ward") maraGroups <- cutree(maraNorm.h, 3); maraGroups maraInitial <- aggregate(maraNorm, by=list(maraGroups), mean)[,-1]; maraInitial maraNorm.km1 <- kmeans(maraNorm,maraInitial) plot(maraNorm.h, main="marajoara320amostras") ### PCA maraNorm.pca <- princomp(maraNorm) maraNorm.px <- predict(maraNorm.pca) library(MASS) eqscplot(maraNorm.px[,1:2],type="n",xlab="componente 1",ylab="componente 2") text(maraNorm.px[,1:2],labels=as.character(maraNorm.km1$cluster),col=as.numeric(maraNorm.km1$cluster),cex=.5) title("marajoara320amostras") #identify(maraNorm.px[,1:2]) # identificar amostras no gráfico de PCA library(ellipse) tmp <- split(data.frame(maraNorm.px[,1:2]), maraGroups) lapply(tmp, cor) lapply(tmp, function(x) lines(ellipse(cor(x), level=0.95, centre=c(mean(x[,1]),mean(x[,2]))), lty=2)) ### </code> Éder Comunello <c <comunello.eder@gmail.com>omunello.eder@gmail.com> Dourados, MS - [22 16.5'S, 54 49'W]
Muito obrigado! Parece que resolveu o problema, vamos tentar verificar o embasamento teórico por detrás. Mauro, boa tarde! Trabalhei um pouco em cima do seu código e cheguei a plotar as elipses, mas precisaria revisar o conteúdo de multivariada pra saber se o procedimento está correto. De qualquer modo, posto o código para que você veja se é possível aproveitar algo. ### <code r> dURL <- ('https://www.dropbox.com/s/405ggda0h6scnie/marajoara320amostras.txt') dURL2 <- ('https://dl.dropboxusercontent.com/s/405ggda0h6scnie/marajoara320amostras.txt...') dName <- basename(dURL); dName download.file(dURL2, dName, mode='wb') file.show(dName) maraOri <- read.table('marajoara320amostras.txt', head=FALSE); head(maraOri) maraNorm <- data.frame(scale(maraOri)); head(maraNorm) ### Cluster maraNorm.h <- hclust(dist(maraNorm)^2, method="ward") maraGroups <- cutree(maraNorm.h, 3); maraGroups maraInitial <- aggregate(maraNorm, by=list(maraGroups), mean)[,-1]; maraInitial maraNorm.km1 <- kmeans(maraNorm,maraInitial) plot(maraNorm.h, main="marajoara320amostras") ### PCA maraNorm.pca <- princomp(maraNorm) maraNorm.px <- predict(maraNorm.pca) library(MASS) eqscplot(maraNorm.px[,1:2],type="n",xlab="componente 1",ylab="componente 2") text(maraNorm.px[,1:2],labels=as.character(maraNorm.km1$cluster),col=as.numeric(maraNorm.km1$cluster),cex=.5) title("marajoara320amostras") #identify(maraNorm.px[,1:2]) # identificar amostras no gráfico de PCA library(ellipse) tmp <- split(data.frame(maraNorm.px[,1:2]), maraGroups) lapply(tmp, cor) lapply(tmp, function(x) lines(ellipse(cor(x), level=0.95, centre=c(mean(x[,1]),mean(x[,2]))), lty=2)) ### </code> Éder Comunello <comunello.eder@gmail.com> Dourados, MS - [22 16.5'S, 54 49'W]
participantes (3)
-
Marcos Vital -
Mauro Sznelwar -
Éder Comunello