## rm(list=ls(all=TRUE)) ##LIMPA A MEMORIA
## cat("\014") #limpa o console
## dev.off() #apaga os graficos, se houver algum
## install.packages("tidyverse")
## install.packages("openxlsx")
## install.packages("vegan")
## install.packages("gplots")
## install.packages("psych")
## install.packages("ggplot2")
library(tidyverse)
## getwd()
## setwd("C:/Seu/Diretório/De/Trabalho")
library(openxlsx)
ppbio <- read.xlsx("D:/Elvio/OneDrive/Disciplinas/_EcoNumerica/5.Matrizes/ppbio06p.xlsx",
rowNames = T,
colNames = T,
sheet = "Sheet1")
ppbio_a <- read.xlsx("D:/Elvio/OneDrive/Disciplinas/_EcoNumerica/5.Matrizes/ppbio06h.xlsx",
rowNames = T,
colNames = T,
sheet = "Sheet1")
str(ppbio)
ppbio_ma <- as.matrix(ppbio) #lê ppbio como uma matriz
str(ppbio_ma)
#ppbio
#ppbio_ma
## getwd()
## ppbio <- read.xlsx(file.choose(),
## rowNames = T, colNames = T,
## sheet = "Sheet1")
## #Lista as colunas
## colnames(ppbio_a)
## #Escolher quais colunas usar por nome
## colnames(ppbio_a)[rev(order(colSums(ppbio_a)))] #ordena por maior soma
## #Usar a função subset()
## m_part <- subset(ppbio_a[, c("a.veloc", "a.temp", "a.do", "a.transp")]) #escolhe variáveis por nome
## m_part <- subset(ppbio_a[, 18:26]) #escolhe as colunas de 18 a 26
## #m_part
m_trab <- (ppbio)
#m_trab <- (m_part)
#Dendrograma
library(vegan)
#relativização/transformação da matriz comunitária
m_trns <- asin(sqrt(decostand
(m_trab, method="total", MARGIN = 2)))
#transformação da matriz ambiental
#m_trns <- sqrt(m_trab)
vegdist <- vegdist(m_trns, method = "bray",
diag = TRUE,
upper = FALSE)
cluster_uas <- hclust(vegdist, method = "average")
plot (cluster_uas, main = "Cluster Dendrogram - Bray-Curtis",
hang = 0.1) #testar com -.01
rect.hclust(cluster_uas, k = 3, h = NULL)
#h = 0.8 fornece os grupos formados na altura h
as.matrix(vegdist)[1:6, 1:6]
#Heatmap
library("gplots")
heatdist <- as.matrix(vegdist)
col <- rev(heat.colors(999)) #rev() reverte as cores do heatmap
heatmap.2(x=(as.matrix(vegdist)), #objetos x objetos
Rowv = as.dendrogram(cluster_uas),
Colv = as.dendrogram(cluster_uas),
key = T, tracecol = NA, revC = T,
col = heat.colors, #dissimilaridade = 1 - similaridade
density.info = "none",
xlab = "UA´s", ylab = "UA´s",
mar = c(6, 6) + 0.2)
cluster_spp <- hclust((vegdist(t(m_trns), method = "bray",
diag = TRUE,
upper = FALSE)), method = "average")
plot (cluster_spp, main = "Dendrograma dos atributos")
heatmap.2(t(as.matrix(m_trns)), #objetos x atributos
Colv = as.dendrogram(cluster_uas),
Rowv = as.dendrogram(cluster_spp),
key = T, tracecol = NA, revC = T,
col = col,
density.info = "none",
xlab = "Unidades amostrais", ylab = "Espécies",
mar = c(6, 6) + 0.1) # adjust margin size
library(gt)
merge <- as.data.frame(cluster_uas$merge)
merge[nrow(merge)+1,] = c("0","0")
height <- as.data.frame(round(cluster_uas$height, 2))
height[nrow(height)+1,] = c("1.0")
fusoes <- data.frame(Cluster_uas = merge, Height = height)
colnames(fusoes) <- c("Cluster1", "Cluster2", "Height")
UA <- rownames_to_column(as.data.frame(m_trns[, 0]))
colnames(UA) <- c("UAs")
No.UA <- 1:nrow(fusoes)
fusoes <- cbind(No.UA, UA, fusoes)
fusoes$Histórico <- 1:nrow(fusoes)
#fusoes
gt(fusoes)
pca <- prcomp(m_trab)
pca
plot(pca, type = "l")
summary(pca)
## #Remover a primeira coluna
## pca_part1 <- prcomp(m_trns[,-1])
##
## #Usar apenas as 5 primeiras colunas
## colnames(m_trns) #lista as colunas
## pca_part2 <- prcomp(m_trns[,1:5])
##
## #Escolher quais colunas usar por nome
## colnames(m_trns)[rev(order(colSums(m_trns)))] #ordena por maior soma
## pca_part3 <- prcomp(~as-bimac + ge-brasi, data = m_trns) #usa apenas as colunas listadas
## #o "-" deve ser substituido, o R não o reconhece como texto
##
## #Usar a função subset()
## pca_part4 <- subset(m_trns[,1:5])
## prcomp(pca_part4, scale = TRUE)
plot(m_trns$"as-bimac", m_trns$"ge-brasi")
#plot(m_trns$"m.elev", m_trns$"m.river")
plot(scale(m_trns$"as-bimac"), scale(m_trns$"ge-brasi"))
#plot(scale(m_trns$"m.elev"), scale(m_trns$"m.river"))
plot((m_trns$"as-bimac" - mean(m_trns$"as-bimac")) / sd(m_trns$"as-bimac"))
#plot((m_trns$"m.elev" - mean(m_trns$"m.elev")) / sd(m_trns$"m.elev"))
library(psych)
pairs.panels(m_trns[,1:5],
method = "pearson", # correlation method
scale = FALSE, lm = FALSE,
hist.col = "#00AFBB", pch = 19,
density = TRUE, # show density plots
ellipses = TRUE, # show correlation ellipses
alpha = 0.5
)
pca_ce <- prcomp(m_trns, center = TRUE, scale = F)
#pca_ce
pca_cs <- prcomp(m_trns, center = TRUE, scale = TRUE)
#pca_cs
pca_sc <- prcomp(m_trns, scale = TRUE)
#pca_sc
par(mfrow = c(3,1))
plot(pca_ce, type = "l")
plot(pca_cs, type = "l")
plot(pca_sc, type = "l")
par(mfrow = c(1,1))
pca_sc #fornece os desvios padrões, Rotação dos eixos e Eigenvetores
str(pca_sc)
summary(pca_sc)
plot(pca_sc, type = "l", main = "Scree Plot (scaled data - m_trns)") #scree plot
add.col <- rownames_to_column(m_trns, var = "UAs")
#add.col
agrup <- substr(add.col[, 1], 5,6) #descendo os nomes
#agrup
m_pca_agrup <- add.col %>% mutate(Agrupamentos=c(agrup),.before=UAs)
m_pca_agrup[1:5, 1:5]
pca_sc <- prcomp(m_trns, scale = TRUE)
pca_sc #fornece os desvios padrões, Rotação dos eixos e Eigenvetores
str(pca_sc)
summary(pca_sc)
plot(pca_sc, type = "l", main = "Scree Plot (scaled data - m_trns)") #scree plot
biplot(pca_sc, choices=1:2, scale = 1,
main="Biplot da PCA. PPBio Comunidade",
xlab = "PC1 UAs",
ylab = "PC2 UAs",
cex = 0.8) #PCA plot
# Calculando o percentual de variação explicado
var_exp <- round((pca_sc$sdev^2/sum(pca_sc$sdev^2))*100, 2)
var_exp
str(pca_sc)
pca_sc$x
pca_scores <- cbind(Agrupamento = m_pca_agrup[,1], m_trns, pca_sc$x[,1:2])
head(pca_scores)
pc1_uas <- round(pca_sc$x[,1],2)
#pc1_uas
pc2_uas <- round(pca_sc$x[,2],2)
#pc2_uas
pc_uas <- data.frame(PC1 = pc1_uas, PC2 = pc2_uas)
sorted_pc_uas <- pc_uas[order(-pc_uas$PC1), ]
sorted_pc_uas <- rownames_to_column(sorted_pc_uas, var = "UAs")
sorted_pc_uas
pc1_spp <- round(pca_sc$rotation[,1],2)
#pc1_spp
pc2_spp <- round(pca_sc$rotation[,2],2)
#pc2_spp
pc_spp <- data.frame(PC1 = pc1_spp, PC2 = pc2_spp)
sorted_pc_spp <- pc_spp[order(-pc_spp$PC1), ]
sorted_pc_spp <- rownames_to_column(sorted_pc_spp, var = "Spp")
sorted_pc_spp
#nrow(sorted_pc_uas)-nrow(sorted_pc_spp)
#sorted_pc_uas[nrow(sorted_pc_uas) + 12,] <- c("<NA>", "<NA>", "<NA>")
#cbind(sorted_pc_uas, sorted_pc_spp)
library(ggplot2)
ggplot(pca_scores, aes(PC1, PC2, col = Agrupametno, fill = Agrupamento)) +
stat_ellipse(geom = "polygon", col = "black", alpha = 0.5) +
geom_point(shape = 21, col = "black")
cor(m_trns, pca_scores[, c("PC1", "PC2")])
## library(FactoMineR)
## library(factoextra)
## pca_fm <- PCA(m_trns, graph = FALSE)
## eigenval <- get_eigenvalue(pca_fm)
## eigenval
## fviz_eig(pca_fm, addlabels = TRUE)
## var <- get_pca_var(pca_fm)
## ind <- get_pca_ind(pca_fm)
## fviz_pca_var(pca_fm, col.var = "red", repel = TRUE)
## fviz_pca_ind(pca_fm, col.ind = "blue", repel = TRUE)
## grupo <- as.factor(m_pca_agrup[,1])
## fviz_pca_biplot(pca_fm, habillage = grupo, title="PCA no FactoMineR", repel = TRUE, addEllipses = FALSE,
## geom.ind = c("point","text"),
## pointshape = 21,
## pointsize = 2,
## fill.ind = m_pca_agrup$Agrupamentos)+
## theme_bw()+
## labs(title = "Biplot PCA no FactoMineR", #substitui o title
## fill = "Grupos") #substitui o habillage
## var$cos2
## library(corrplot)
## corrplot(pca_fm$var$cos2, is.corr = T)
## corrplot(pca_fm$var$contrib, is.corr = F)
## fviz_contrib(pca_fm, choice = "var", axes = 1, top = 10, title = "10 mais para a PC1")
## fviz_contrib(pca_fm, choice = "var", axes = 2, top = 10, title = "10 mais para a PC2")
## pca_sign <- dimdesc(pca_fm, axes = c(1,2,3), proba = 0.05)
## pca_signScript do R Módulo 9
Script limpo
Aqui apresento o scrip na íntegra sem os textos ou outros comentários. Você pode copiar e colar no R para executa-lo. Lembre de remover os # ou ## caso necessite executar essas linhas.