Este é um documento anexo que pertence ao artigo intitulado “Analyzing beta-diversity and the biogeographical delimitation of forest and woody savanna formations in South America”, de autoria de Ângela Lúcia Bagnatori Sartori, Geraldo Alves Damasceno Junior, Emerson Pereira da Silva, Giovana Amaral Umar, Isabela Oliveira Ferreira, Mateus César Araújo Pestana, Maxwell da Rosa Oliveira, Thomaz Ricardo Favreto Sinani e Flávio Macedo Alves
Estas análises são baseadas em 112 referencias artigos e base de dados dos autores realizados em remanescentes de Chaco e ecorregiões adjacentes. Nós compilamos dados de ocorrência de 1.628 táxons da flora em 307 localidades, que podem ser acessados no Material Suplementar “matrix_original.xlsx” e as análises replicadas seguindo este roteiro.
library(readxl)
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.1 ✔ readr 2.2.0
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ ggplot2 4.0.3 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2
## ✔ purrr 1.2.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(vegan)
## Carregando pacotes exigidos: permute
library(betapart)
library(venn)
Anotações: se houver erro, os pacotes devem ser previamente instalados.
# 2.1. Importar a planilha original
input_matrix <- read_excel("matrix_original.xlsx")
# 2.2. Processar e pivotar os dados
pivot_matrix <- input_matrix %>%
# Selecionar apenas as colunas que contêm localidade e taxa
select(locality, taxa) %>%
# Remover registros duplicados da mesma espécie na mesma localidade
distinct(locality, taxa) %>%
# Criar uma coluna de preenchimento indicando a presença (1)
mutate(presence = 1) %>%
# Pivotamos o dataframe: cada 'taxa' vira uma coluna
pivot_wider(names_from = taxa, values_from = presence, values_fill = 0)
# 2.3. Ajustar os nomes das linhas (rownames) e converter para matriz numérica
# Transformar em data.frame clássico para habilitar rownames
matriz_df <- as.data.frame(pivot_matrix)
# Definir a coluna 'locality' como o identificador de cada linha
rownames(matriz_df) <- matriz_df$locality
# Remover a coluna de texto 'locality' que ficou residual
matriz_df <- matriz_df %>% select(-locality)
# Converter formalmente para a classe 'matrix' do R
final_matrix <- as.matrix(matriz_df)
# 2.4. Checagem rápida do resultado
# Exibir as 5 primeiras linhas e colunas para conferência
as_tibble(final_matrix[1:5, 1:5], rownames = "Localidade")
## # A tibble: 5 × 6
## Localidade Acrocomia_aculeata_(Jacq.)_Lodd._ex_M…¹ Agonandra_brasiliens…²
## <chr> <dbl> <dbl>
## 1 CERDAO_CER_01 1 1
## 2 CER.S.ST_CER_01 0 1
## 3 CERDAO_CER_02 1 1
## 4 CERDAO_CER_03 0 1
## 5 CERDAO_CER_04 0 0
## # ℹ abbreviated names: ¹`Acrocomia_aculeata_(Jacq.)_Lodd._ex_Mart.`,
## # ²`Agonandra_brasiliensis_Benth._&_Hook.f.`
## # ℹ 3 more variables: Aspidosperma_subincanum_Mart. <dbl>,
## # Aspidosperma_tomentosum_Mart. <dbl>, Astronium_fraxinifolium_Schott <dbl>
# Verificar o número total de localidades (linhas) e espécies (colunas)
print(dim(final_matrix))
## [1] 307 1628
Anotações: a planilha “matrix_original.xlsx” está disponível para download no material suplementar.
# Calcular o número total de ocorrências de cada espécie (soma das colunas)
frequency_sp <- colSums(final_matrix)
# Ordenar as espécies da mais comum para a mais rara
frequency_order <- sort(frequency_sp, decreasing = TRUE)
# Converter para percentagem (proporção face às N localidades)
percent_ocur <- (frequency_order / nrow(final_matrix)) * 100
# Criar o gráfico de Frequência-Ordem
plot(percent_ocur,
type = "l", # Desenha uma linha
lwd = 2, col = "blue",
main = "Curva de Frequencia-Ordem das Angiospermas",
xlab = "Ordenacao das Especies (da mais comum a mais rara)",
ylab = "Frequencia de Ocorrencia (%)")
# Adicionar uma linha de corte visual para espécies raras (ex: presentes em menos de 5% das áreas)
abline(h = 5, col = "red", lty = 2)
# Extrair dados interessantes para a sua tese/artigo:
cat("A espécie mais comum ocorre em", round(percent_ocur[1], 1), "% das áreas.\n")
## A espécie mais comum ocorre em 35.8 % das áreas.
cat("Número de espécies hiper-raras (<5% de ocorrência):", sum(percent_ocur < 5), "\n")
## Número de espécies hiper-raras (<5% de ocorrência): 1357
Anotações: o gráfico inicia alto à esquerda e cai rapidamente para formar uma cauda em “J” invertido; o alto cruzamento da linha vermelha (espécies raras) no eixo X determina se o conjunto de espécies especialistas locais ou de transição.
# Preparar a matriz (o R desenha de baixo para cima, então invertemos para ficar intuitivo)
matrix_image <- as.matrix(final_matrix)
matrix_image <- t(matrix_image)[, nrow(matrix_image):1]
# Gerar o Heatmap
image(x = 1:ncol(final_matrix),
y = 1:nrow(final_matrix),
z = matrix_image,
col = c("white", "black"), # 0 (ausência) é branco, 1 (presença) é preto
main = "Mapa de Calor da Matriz de Presença/Ausência",
xlab = "Espécies",
ylab = "Localidades")
Anotações: Os pontos pretos no heatmap representam presença (1) e a parte branca representa ausência (0); é possível detectar agrupamentos prévios das espécies nas localidades.
# Calcular a curva de acumulação exata
# 'random' vai misturar as suas N áreas dezenas de vezes para criar uma curva suave
regional_rare <- specaccum(final_matrix, method = "random", permutations = 100)
# Plotar a curva
plot(regional_rare,
ci.type = "poly", col = "darkgreen", lwd = 2,
ci.lty = 0, ci.col = "lightgreen",
main = "Curva de Acumulacao Regional de Especies",
xlab = "Esforco Amostral (Número de Localidades Adicionadas)",
ylab = "Riqueza Acumulada de Especies")
# Adicionar uma grelha para facilitar a leitura
grid()
Anotações: Foi construída uma curva de acumulação de espécies baseada em amostras (dados de incidência em matriz binária) utilizando permutações aleatórias das localidades.
# Calcular a partição usando a sua matriz binária limpa
part_beta <- beta.multi(final_matrix, index.family = "jaccard")
# Visualizar os resultados gerais da região
cat("Dissimilaridade Jaccard Total (Regiao):", round(part_beta$beta.JAC, 3), "\n")
## Dissimilaridade Jaccard Total (Regiao): 0.997
cat("Fracao por Turnover (Substituicao - Biomas diferentes):", round(part_beta$beta.JTU, 3), "\n")
## Fracao por Turnover (Substituicao - Biomas diferentes): 0.995
cat("Fracao por Aninhamento (Perda de especies):", round(part_beta$beta.JNE, 3), "\n")
## Fracao por Aninhamento (Perda de especies): 0.002
Anotações: Quando o Turnover (JTU) é muito maior que o Aninhamento (JNE), há evidências fortes de que a área de estudo abrange diferentes domínios florísticos.
# Extrair os nomes das localidades da matriz de presença/ausência
# Transformar o nome das linhas em coluna
df_locality <- tibble(locality = rownames(final_matrix))
# Fazer o Cruzamento (Left Join) para encontrar o 'dominio_bioma' correspondente
df_biome_align <- df_locality %>%
left_join(input_matrix %>% select(locality, domain_biome) %>%
distinct(), by = "locality")
df_biome_align <- df_locality %>%
left_join(input_matrix %>% select(locality, domain_biome) %>% distinct(),
by = "locality")
# Criar o vetor final de biomas
biome_vec <- df_biome_align$domain_biome
# Descobrir quais são os biomas únicos presentes nos dados
biome_uniq <- unique(na.omit(biome_vec))
# Criar a lista de espécies para o diagrama de forma automática
# Usar lapply para varrer cada bioma e extrair as espécies presentes nele
list_venn_dinamic <- lapply(biome_uniq, function(biome_actual) {
# Filtrar a matriz apenas para as linhas deste bioma
row_biome <- final_matrix[biome_vec == biome_actual, , drop = FALSE]
# Quais espécies têm soma > 0 nestas linhas?
present_sp <- colnames(row_biome)[colSums(row_biome) > 0]
return(present_sp)
})
# Dar os nomes corretos a cada conjunto na lista
names(list_venn_dinamic) <- biome_uniq
# Colar a riqueza total de cada bioma diretamente no nome da categoria
biome_name <- names(list_venn_dinamic)
biome_rich <- sapply(list_venn_dinamic, length)
# Atualizar os nomes na lista com o valor entre parênteses
names(list_venn_dinamic) <- paste0(biome_name, "\n(", biome_rich, " spp.)")
# Gerar o gráfico
venn(list_venn_dinamic,
# --- CORES E SOBREPOSIÇÃO ---
zcolor = "style",
opacity = 0.4,
# --- FORÇAR EXIBIÇÃO DOS NÚMEROS ---
counts = TRUE,
ilcs = 0.5,
sncs = 0.8,
# --- ESTÉTICA GERAL ---
box = FALSE,
borders = FALSE)
Anotações: Os biomas indicados no diagrama de Venn foram obtidos das informações presentes do material consultado com a finalidade de obter uma visão geral desta análise exploratória de dados.
# Calcular a matriz de distância florística (Dissimilaridade de Jaccard)
dissimilarity_matrix <- vegdist(final_matrix, method = "jaccard")
# Criar a Árvore de Agrupamento (Clustering) com método 'ward.D2', excelente para criar grupos densos e bem definidos
tree_biome <- hclust(dissimilarity_matrix, method = "ward.D2")
# Desenhar a árvore para decidir o número ótimo de grupos (k)
plot(tree_biome, labels = FALSE, hang = -1,
main = "Dendrograma Florístico: Procurando os Biomas Reais",
xlab = "Localidades", ylab = "Distância de Ward")
# Baseado no dendrograma anterior, o k pode ser ajustado (aqui foi mantido k=2)
k_excell <- 2
# Cortar a árvore para gerar 2 grupos
my_data_group <- cutree(tree_biome, k = k_excell)
# Extrair os grupos da árvore
my_data_group <- factor(my_data_group,
levels = c(1, 2),
labels = c("Floristic Domain A", "Floristic Domain B"))
# Correr o NMDS
nmds_result <- metaMDS(dissimilarity_matrix, k = 2, trymax = 100)
## Run 0 stress 0.1481864
## Run 1 stress 0.1505784
## Run 2 stress 0.1545715
## Run 3 stress 0.1559888
## Run 4 stress 0.1516041
## Run 5 stress 0.1522861
## Run 6 stress 0.1537738
## Run 7 stress 0.1504466
## Run 8 stress 0.1514693
## Run 9 stress 0.1517982
## Run 10 stress 0.152757
## Run 11 stress 0.1513451
## Run 12 stress 0.1500807
## Run 13 stress 0.1539637
## Run 14 stress 0.152988
## Run 15 stress 0.1555398
## Run 16 stress 0.1519908
## Run 17 stress 0.1547997
## Run 18 stress 0.1543904
## Run 19 stress 0.1514179
## Run 20 stress 0.1515637
## Run 21 stress 0.1539893
## Run 22 stress 0.1534569
## Run 23 stress 0.1487386
## Run 24 stress 0.1541035
## Run 25 stress 0.1533714
## Run 26 stress 0.1531064
## Run 27 stress 0.1456277
## ... New best solution
## ... Procrustes: rmse 0.03291639 max resid 0.237303
## Run 28 stress 0.1498658
## Run 29 stress 0.1502705
## Run 30 stress 0.1509448
## Run 31 stress 0.1547399
## Run 32 stress 0.1451041
## ... New best solution
## ... Procrustes: rmse 0.02026104 max resid 0.1355561
## Run 33 stress 0.1517035
## Run 34 stress 0.1525734
## Run 35 stress 0.1500419
## Run 36 stress 0.1537767
## Run 37 stress 0.1478305
## Run 38 stress 0.1498785
## Run 39 stress 0.1509294
## Run 40 stress 0.1488916
## Run 41 stress 0.1563524
## Run 42 stress 0.1566497
## Run 43 stress 0.149964
## Run 44 stress 0.1506817
## Run 45 stress 0.1518778
## Run 46 stress 0.1521549
## Run 47 stress 0.152964
## Run 48 stress 0.1475119
## Run 49 stress 0.1529887
## Run 50 stress 0.1487619
## Run 51 stress 0.1455183
## ... Procrustes: rmse 0.02280899 max resid 0.1549416
## Run 52 stress 0.1523059
## Run 53 stress 0.1464137
## Run 54 stress 0.1512414
## Run 55 stress 0.1479771
## Run 56 stress 0.1541085
## Run 57 stress 0.1547149
## Run 58 stress 0.1551865
## Run 59 stress 0.1493155
## Run 60 stress 0.150717
## Run 61 stress 0.152267
## Run 62 stress 0.1506578
## Run 63 stress 0.1529446
## Run 64 stress 0.1534701
## Run 65 stress 0.1517628
## Run 66 stress 0.1557881
## Run 67 stress 0.1539189
## Run 68 stress 0.1486776
## Run 69 stress 0.1526171
## Run 70 stress 0.1504691
## Run 71 stress 0.1556569
## Run 72 stress 0.152415
## Run 73 stress 0.1559738
## Run 74 stress 0.1529853
## Run 75 stress 0.1513983
## Run 76 stress 0.1525546
## Run 77 stress 0.1522805
## Run 78 stress 0.1516075
## Run 79 stress 0.1527916
## Run 80 stress 0.152411
## Run 81 stress 0.1550271
## Run 82 stress 0.1515647
## Run 83 stress 0.1492319
## Run 84 stress 0.1546589
## Run 85 stress 0.1492685
## Run 86 stress 0.155874
## Run 87 stress 0.1479237
## Run 88 stress 0.1533925
## Run 89 stress 0.1452041
## ... Procrustes: rmse 0.02102635 max resid 0.1488202
## Run 90 stress 0.1546082
## Run 91 stress 0.1552715
## Run 92 stress 0.1501499
## Run 93 stress 0.152561
## Run 94 stress 0.1542943
## Run 95 stress 0.1525859
## Run 96 stress 0.1556779
## Run 97 stress 0.1500961
## Run 98 stress 0.1538947
## Run 99 stress 0.1505794
## Run 100 stress 0.1514099
## *** Best solution was not repeated -- monoMDS stopping criteria:
## 66: no. of iterations >= maxit
## 34: stress ratio > sratmax
# Juntar as Coordenadas com os Novos Grupos Gerados
df_nmds_blind <- data.frame(
Locality = rownames(final_matrix),
NMDS1 = nmds_result$points[, 1],
NMDS2 = nmds_result$points[, 2],
Floristic_Group = my_data_group
)
# Plotar o Gráfico Final e Cego
ggplot(df_nmds_blind, aes(x = NMDS1, y = NMDS2, color = Floristic_Group, fill = Floristic_Group)) +
geom_point(size = 3, alpha = 0.7) +
# As elipses agora são baseadas nos dados puros, não na sua planilha
stat_ellipse(geom = "polygon", alpha = 0.15, linetype = "solid", size = 0.8) +
theme_bw() +
labs(title = "Ordenacao Floristica Data-Driven (NMDS)",
subtitle = paste("Separação primaria em 2 biomas (Ward) | Stress =", round(nmds_result$stress, 3)),
x = "Eixo NMDS 1",
y = "Eixo NMDS 2",
color = "Comunidade (k)",
fill = "Comunidade (k)") +
# Usar cores distintas para os grupos identificados
scale_color_manual(values = c("#d95f02", "#1b9e77")) +
scale_fill_manual(values = c("#d95f02", "#1b9e77"))+
theme(legend.position = "bottom",
plot.title = element_text(face = "bold"))
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
# O Boxplot focando no Eixo 1 do NMDS
ggplot(df_nmds_blind, aes(x = Floristic_Group, y = NMDS1, fill = Floristic_Group)) +
# Criar as caixas. O outlier.shape destaca localidades atípicas dentro do próprio bioma
geom_boxplot(alpha = 0.8, outlier.shape = 16, outlier.size = 2, outlier.color = "red") +
theme_bw() +
labs(title = "Separacao Floristica ao Longo do Gradiente Principal",
subtitle = "Distribuicao das localidades no Eixo 1 do NMDS",
x = "Dominio Fitogeografico (k=2)",
y = "Coordenadas do Eixo 1 (NMDS1)") +
scale_fill_manual(values = c("#d95f02", "#1b9e77")) +
theme(legend.position = "none",
plot.title = element_text(face = "bold"),
axis.text.x = element_text(size = 11, face = "bold"))
Anotações: o teste NMDS seguiu uma análise não supervisionada (Data-Driven) para evitar viés tendencioso com a informação dos biomas da matriz de dados; o dendrograma apresentou 2 grandes ramos, então k=2; uma localidade foi detectada como outlier no gráfico de ordenação e mais de 3 no bloxpot; apesar de apresentar outliers, o Stress foi <0.20, que pode ser usado.
# Correr a PCoA com base na matriz de distância de Jaccard e objetos do tópico de NMDS
pcoa_result <- cmdscale(dissimilarity_matrix, k = 2, eig = TRUE)
# Calcular a variância explicada pelos dois eixos principais
percent_axis1 <- round(pcoa_result$eig[1] / sum(pcoa_result$eig) * 100, 1)
percent_axis2 <- round(pcoa_result$eig[2] / sum(pcoa_result$eig) * 100, 1)
# Montar o Data Frame com as Coordenadas da PCoA e os Grupos do Cluster
df_pcoa_blind <- data.frame(
Locality = rownames(final_matrix),
PCoA1 = pcoa_result$points[, 1],
PCoA2 = pcoa_result$points[, 2],
Floristic_Group = my_data_group
)
# Plotar o Gráfico Final da PCoA Data-Driven
ggplot(df_pcoa_blind, aes(x = PCoA1, y = PCoA2, color = Floristic_Group, fill = Floristic_Group)) +
geom_point(size = 3, alpha = 0.7) +
stat_ellipse(geom = "polygon", alpha = 0.15, linetype = "solid", size = 0.8) +
theme_bw() +
labs(title = "Ordenacao Floristica Data-Driven (PCoA)",
subtitle = "Grupos definidos cegamente pela flora (Ward Clustering)",
x = paste0("Eixo PCoA 1 (", percent_axis1, "%)"),
y = paste0("Eixo PCoA 2 (", percent_axis2, "%)"),
color = "Comunidade (k)",
fill = "Comunidade (k)") +
theme(legend.position = "bottom",
plot.title = element_text(face = "bold"))
# O Boxplot focando no Eixo 1 da PCoA
ggplot(df_pcoa_blind, aes(x = Floristic_Group, y = PCoA1, fill = Floristic_Group)) +
# Criar as caixas. O outlier.shape destaca localidades atípicas dentro do próprio bioma
geom_boxplot(alpha = 0.8, outlier.shape = 16, outlier.size = 2, outlier.color = "red") +
theme_bw() +
labs(title = "Separacao Floristica ao Longo do Gradiente Principal",
subtitle = "Distribuicao das localidades no Eixo 1 da PCoA",
x = "Dominio Fitogeografico (k=2)",
y = "Coordenadas do Eixo 1 (PCoA1)") +
scale_fill_manual(values = c("#d95f02", "#1b9e77")) +
theme(legend.position = "none",
plot.title = element_text(face = "bold"),
axis.text.x = element_text(size = 11, face = "bold"))
Anotações: diferentemente do NMDS, o PCoA não apresentou outliers e foca na distância métrica linear; poderia ser uma alternativa ao NMDS.
# Criar um dataframe simples contendo os grupos gerados pelo corte da árvore (k=2)
df_permanova_blind <- data.frame(
Floristic_Group = my_data_group)
# Executar a PERMANOVA usando a matriz de Jaccard
permanova_test <- adonis2(dissimilarity_matrix ~ Floristic_Group,
data = df_permanova_blind,
permutations = 999)
print("--- RESULTADO DA PERMANOVA ---")
## [1] "--- RESULTADO DA PERMANOVA ---"
print(permanova_test)
## Permutation test for adonis under reduced model
## Permutation: free
## Number of permutations: 999
##
## adonis2(formula = dissimilarity_matrix ~ Floristic_Group, data = df_permanova_blind, permutations = 999)
## Df SumOfSqs R2 F Pr(>F)
## Model 1 12.861 0.09353 31.47 0.001 ***
## Residual 305 124.642 0.90647
## Total 306 137.502 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Calcular a dispersão e extrair coordenadas
dispers_grupos <- betadisper(dissimilarity_matrix, df_permanova_blind$Floristic_Group)
# Criar dataframe com as coordenadas dos pontos (Localidades)
df_perm_scatter <- data.frame(
PCoA1 = dispers_grupos$vectors[, 1],
PCoA2 = dispers_grupos$vectors[, 2],
Floristic_Group = df_permanova_blind$Floristic_Group
)
# Criar dataframe com as coordenadas dos "Corações" (Centróides)
df_centroids <- data.frame(
Floristic_Group = rownames(dispers_grupos$centroids),
Centroid1 = dispers_grupos$centroids[, 1],
Centroid2 = dispers_grupos$centroids[, 2]
)
# Juntar as informações para poder desenhar as linhas conectando pontos aos centróides
df_perm_scatter <- df_perm_scatter %>% left_join(df_centroids, by = "Floristic_Group")
# Plotar a Teia (Spider Plot)
ggplot(df_perm_scatter, aes(color = Floristic_Group)) +
# Desenha a "teia" (linha ligando o ponto ao centro)
geom_segment(aes(x = Centroid1, y = Centroid2, xend = PCoA1, yend = PCoA2), alpha = 0.3) +
# Desenha os pontos normais das localidades
geom_point(aes(x = PCoA1, y = PCoA2), size = 2, alpha = 0.7) +
# Desenha os Centróides em destaque (Losango preto)
geom_point(aes(x = Centroid1, y = Centroid2), size = 6, shape = 18, color = "black") +
theme_bw() +
labs(title = "Dispersao da PERMANOVA (Teia de Centróides)",
subtitle = "Projecao multivariada das distâncias florísticas ao centro de cada grupo",
x = "Eixo Principal 1",
y = "Eixo Principal 2") +
scale_color_manual(values = c("#d95f02", "#1b9e77")) +
theme(legend.position = "bottom", plot.title = element_text(face = "bold"))
# Inserir os valores calculados no dataframe original
df_permanova_blind$Centroid_Distance <- dispers_grupos$distances
# Plotar o Boxplot
ggplot(df_permanova_blind, aes(x = Floristic_Group, y = Centroid_Distance, fill = Floristic_Group)) +
geom_boxplot(alpha = 0.8, outlier.shape = 16, outlier.size = 2, outlier.color = "red") +
theme_bw() +
labs(title = "Variância Interna dos Grupos (PERMANOVA)",
subtitle = "Mediana das distâncias de cada localidade para o núcleo do seu bioma",
x = "Domínio Fitogeográfico (k=2)",
y = "Distância Florística (Jaccard)") +
scale_fill_manual(values = c("#d95f02", "#1b9e77")) +
theme(legend.position = "none",
plot.title = element_text(face = "bold"),
axis.text.x = element_text(size = 11, face = "bold"))
Anotações: para verificar se a separação de grupos via NMDS ou PCoA ficaram estatisticamente significativas, usou-se a PERMANOVA; o p-value (Pr(>F)) foi XXX por ser >< X; O R2 (R-squared) diz a porcentagem exata de toda a variação da vegetação de cada grupo formado.
# Fazer o Cruzamento (Left Join) para buscar as coordenadas na matriz bruta
df_coords_align <- df_locality %>%
left_join(input_matrix %>%
select(locality, longitud, latitud) %>%
distinct(),
by = "locality")
# Limpeza segura das coordenadas (Tratar vírgulas e converter para números)
# O mutate do dplyr faz isso de forma muito limpa
df_coords_align <- df_coords_align %>%
mutate(
longitud = as.numeric(gsub(",", ".", longitud)),
latitud = as.numeric(gsub(",", ".", latitud))
)
# Preparar para a distância matemática
# A função vegdist só aceita números, então precisamos remover a coluna de texto 'locality'
coords_numericas <- df_coords_align %>%
select(longitud, latitud)
# Calcular a Distância Geográfica e Executar o Teste de Mantel
geo_distance <- vegdist(coords_numericas, method = "euclidean")
mantel_test <- mantel(dissimilarity_matrix, geo_distance,
method = "spearman",
permutations = 999)
# Imprimir o resultado
print("--- RESULTADO DO TESTE DE MANTEL ---")
## [1] "--- RESULTADO DO TESTE DE MANTEL ---"
print(mantel_test)
##
## Mantel statistic based on Spearman's rank correlation rho
##
## Call:
## mantel(xdis = dissimilarity_matrix, ydis = geo_distance, method = "spearman", permutations = 999)
##
## Mantel statistic r: 0.4571
## Significance: 0.001
##
## Upper quantiles of permutations (null model):
## 90% 95% 97.5% 99%
## 0.0249 0.0314 0.0398 0.0480
## Permutation: free
## Number of permutations: 999
# 'Achatar' as matrizes para formato vetorial
vec_floristic <- as.vector(dissimilarity_matrix)
vec_geo <- as.vector(geo_distance)
# Montar dataframe cruzado
df_mantel_visual <- data.frame(
Dist_Geo = vec_geo,
Dist_Flor = vec_floristic
)
# Plotar a regressão contínua
ggplot(df_mantel_visual, aes(x = Dist_Geo, y = Dist_Flor)) +
geom_point(alpha = 0.08, color = "#2c3e50", size = 0.4) +
geom_smooth(method = "lm", formula = y ~ x, color = "red", linetype = "dashed", size = 1, se = FALSE) +
theme_bw() +
labs(title = "Isolamento por Distancia (Dispersao de Mantel)",
subtitle = "Autocorrelacao espacial: Dissimilaridade floristica vs. Separacao geografica",
x = "Distancia Geografica Euclidiana",
y = "Dissimilaridade Floristica (Jaccard)") +
theme(plot.title = element_text(face = "bold"))
# Fatiar a distância geográfica contínua em 4 classes discretas
# A função 'cut' gera intervalos automáticos baseados nos seus dados reais
df_mantel_visual$Class_Geo <- cut(df_mantel_visual$Dist_Geo,
breaks = 4,
labels = c("Curta Distancia", "Media-Curta", "Media-Longa", "Longa Distancia"))
# Plotar o Boxplot Espacial
ggplot(df_mantel_visual, aes(x = Class_Geo, y = Dist_Flor, fill = Class_Geo)) +
geom_boxplot(alpha = 0.8, outlier.shape = 16, outlier.size = 1, outlier.color = "grey50") +
theme_bw() +
labs(title = "Decaimento da Similaridade em Classes (Mantel Boxplot)",
subtitle = "Evolucao da dissimilaridade da flora conforme o afastamento geografico",
x = "Intervalos de Distancia Geografica",
y = "Dissimilaridade Floristica (Jaccard)") +
# Usar uma paleta de cores gradiente (ex: azul claro para distâncias curtas, escuro para longas)
scale_fill_brewer(palette = "Blues") +
theme(legend.position = "none",
plot.title = element_text(face = "bold"),
axis.text.x = element_text(size = 11, face = "bold"))
Anotações: Se o valor for muito baixo (próximo de zero ou negativo) e/ou o p-value não for significativo, então a distância geográfica no mapa não tem poder para explicar a flora, isso prova que as duas elipses distintas que viu na PCoA representam blocos fitogeográficos verdadeiros, e não apenas um gradiente espacial contínuo.
# Calcular a riqueza observada de cada espécie para cada localidade
riqueza_sp <- specnumber(final_matrix)
riqueza_sp
## CERDAO_CER_01 CER.S.ST_CER_01 CERDAO_CER_02
## 47 28 94
## CERDAO_CER_03 CERDAO_CER_04 CERDAO_CER_05
## 61 77 57
## CERDAO_CER_06 CERDAO_CER_07 CERDAO_CER_08
## 28 74 63
## CERDAO_CER_09 CERDAO_CER_10 CERDAO_CER_11
## 112 92 77
## CERDAO_CER_23 CERDAO_CER_24 CER.S.ST_PANT_01
## 63 50 29
## CERDAO_CER_12 CERDAO_CER_13 CILIAR_CER_01
## 149 132 16
## CILIAR_CER_02 CILIAR_CER_03 CILIAR_CER_04
## 9 22 23
## CILIAR_CER_05 CAMP.MUR_CER_01 CAMP.MUR_CER_02
## 13 92 39
## CERDAO_CER_14 CER.RUPES_CER_01 CILIAR_CER_06
## 25 65 61
## CER.S.ST_CER_02 CILIAR_CER_07 CILIAR_CER_08
## 57 111 61
## CILIAR_CER_09 CERDAO_CER_15 CERDAO_CER_25
## 111 25 45
## CERDAO_CER_26 CERDAO_CER_27 FES_PANT_17
## 25 39 52
## CILIAR_CER_10 CILIAR_CER_11 CERDAO_CER_28
## 111 111 77
## CERDAO_CER_29 CERDAO_CER_30 CERDAO_CER_31
## 14 43 15
## CILIAR_CER_12 CERDAO_CER_16 CERDAO_CER_17
## 61 33 143
## CERDAO_CER_18 TRANS_CER_01 CHA_UM_CHA_01
## 74 86 62
## CHA_UM_CHA_02 CHA_UM_CHA_03 CHA_SEC_01
## 57 87 93
## CHA_SEC_02 CHA_SEC_03 CHA_SEC_04
## 64 48 48
## CHA_SEC_05 CHA_SEC_06 CHA_SEC_07
## 68 63 77
## CHA_ARID_01 CHA_SERRA_01 CHA_BOL_01
## 76 93 20
## CHA_PAMPA_01 CHA_PANT_01 CHA_PANT_02
## 52 64 37
## CHA_PANT_03 CHA_PANT_04 CHA_UM_CHA_04
## 93 43 44
## CHA_UM_CHA_05 CHA_UM_CHA_06 CHA_PANT_05
## 136 123 46
## CHA_PANT_06 CHA_UM_CHA_07 CHA_BOL_02
## 61 33 27
## CHA_SEM_SAL_01 CHA_SEM_ARID_01 CHA_SERRA_02
## 31 44 60
## CHA_PANT_07 CHA_PANT_08 FED_PANT_01
## 102 114 48
## FED_CER_01 FED_PANT_02 FED_PANT_03
## 56 18 39
## CILIAR_CER_13 CILIAR_MA_01 CILIAR_CER_14
## 169 127 34
## CILIAR_CER_15 CILIAR_CER_16 CILIAR_CER_17
## 5 1 13
## FES_CER_01 FES_AMAZ_01 FES_CER_02
## 42 63 56
## FES_CER_03 FES_MA_01 FES_MA_02
## 89 59 42
## FES_MA_03 FES_MA_04 FES_PANT_01
## 63 76 50
## FES_CHI_01 CILIAR_MA_02 FES_CER_04
## 30 51 62
## FES_CER_05 FES_CER_06 FES_CER_07
## 30 28 50
## FES_PANT_02 FES_PANT_03 FES_PANT_04
## 69 61 68
## FES_PANT_05 FES_PANT_06 FES_CER_08
## 46 56 30
## FES_CER_09 FES_CER_10 FES_CER_11
## 44 35 82
## FES_MA_05 FES_MA_06 FES_MA_07
## 51 41 57
## FES_MA_08 FED_PANT_04 PARAT_PANT_03
## 46 20 75
## CERDAO_PANT_01 CERDAO_PANT_02 CERDAO_PANT_03
## 19 47 26
## CILIAR_PANT_01 TRANS_CER_02 CAPAO_PANT_13
## 40 23 10
## CERDAO_PANT_23 CAPAO_PANT_07 CAPAO_PANT_06
## 88 76 81
## CILIAR_PANT_02 CILIAR_PANT_03 CILIAR_PANT_04
## 60 40 12
## CILIAR_PANT_05 CILIAR_PANT_06 CAMP.CERR_PANT_05
## 37 20 10
## FES_PANT_18 FED_PANT_05 CAPAO_PANT_01
## 20 33 42
## CAPAO_PANT_08 CERDAO_PANT_04 CERDAO_PANT_05
## 166 36 32
## CERDAO_PANT_06 CERDAO_PANT_07 CERDAO_PANT_08
## 92 41 45
## CERDAO_PANT_09 CERDAO_PANT_10 CERDAO_PANT_11
## 30 35 24
## CAPAO_PANT_02 CERDAO_PANT_20 CAPAO_PANT_03
## 14 84 28
## CAMP.CERR_PANT_01 CERDAO_PANT_12 CILIAR_PANT_07
## 23 25 30
## CERDAO_PANT_19 CERDAO_PANT_21 FED_PANT_07
## 81 29 57
## FED_PANT_08 FED_PANT_09 FED_PANT_10
## 58 44 50
## CAMP.CERR_PANT_02 CAPAO_PANT_14 PARAT_PANT_02
## 20 20 13
## CAMP.MUR_PANT_01 CAPAO_PANT_12 CAPAO_PANT_09
## 28 16 34
## CAMP.MUR_PANT_02 CERDAO_PANT_13 CAPAO_PANT_04
## 20 14 25
## CAPAO_PANT_05 FES_PANT_07 CAMP.CERR_PANT_03
## 25 74 12
## CERDAO_PANT_14 CAPAO_PANT_11 CILIAR_PANT_08
## 67 76 12
## FES_PANT_08 CILIAR_PANT_09 CILIAR_PANT_10
## 40 36 21
## CILIAR_PANT_11 CERDAO_PANT_22 MAT_ACURI_PANT_01
## 66 60 44
## CAPAO_PANT_10 CAMB_PANT_01 CAPAO_PANT_15
## 113 38 127
## CERDAO_PANT_15 CERDAO_CER_19 CILIAR_PANT_12
## 78 75 97
## CERDAO_PANT_16 FES_PANT_09 FES_PANT_10
## 55 33 54
## CILIAR_PANT_13 CILIAR_PANT_14 FESM_PANT_01
## 45 23 49
## CILIAR_PANT_15 FES_PANT_11 CAMP.CERR_PANT_04
## 31 24 34
## FES_PANT_12 CAMP.CERR_CER_01 CAMP.MUR_CER_03
## 91 34 81
## CER.S.ST_CER_03 CERDAO_CER_20 BABA_CER_01
## 116 55 41
## CAPAO_CER_01 CILIAR_CER_18 CILIAR_MA_03
## 32 82 52
## CERDAO_CER_21 FED_PANT_06 CERDAO_CER_22
## 75 41 40
## CHA_PANT_09 CHA_SEM_ARID_02 CHA_SEM_ARID_03
## 61 32 32
## CHA_SEM_ARID_04 CHA_SEM_ARID_05 CHA_SEM_ARID_06
## 32 32 32
## CHA_SEM_ARID_07 CHA_SEM_ARID_08 CHA_PANT_10
## 20 13 23
## CERDAO_PANT_17 CERDAO_PANT_18 FES_PANT_13
## 35 44 38
## FES_PANT_14 CER.S.ST_PANT_02 FES_PANT_15
## 40 25 27
## FES_PANT_16 CAMB_PANT_02 CAMP.MUR_PANT_03
## 35 13 9
## CHA_AREN_01 CHA_AREN_02 CHA_ENC_AREN_02
## 30 28 23
## CHA_INTERDUN_01 CHA_INTERDUN_02 CHA_INTERDUN_03
## 25 35 36
## CHA_LEQ_ALU_01 CHA_LEQ_ALU_02 CHA_LEQ_ALU_03
## 25 35 29
## CHA_LEQ_ALU_04 CHA_LEQ_ALU_05 CHA_LEQ_ALU_06
## 29 37 33
## CHA_LEQ_ALU_07 CHA_LEQ_ALU_08 CHA_LEQ_ALU_09
## 30 35 32
## CHA_LEQ_ALU_10 CHA_LEQ_ALU_11 CHA_LEQ_ALU_12
## 36 35 36
## CHA_LEQ_ALU_13 CHA_LEQ_ALU_14 CHA_LEQ_ALU_15
## 28 21 29
## CHA_LEQ_ALU_16 CHA_LEQ_ALU_17 CHA_LEQ_ALU_18
## 29 27 37
## CHA_LEQ_ALU_19 CHA_CHIQ_01 CHA_CHIQ_02
## 28 39 44
## CHA_CHIQ_03 CHA_CHIQ_04 CHA_CHIQ_05
## 37 41 36
## CHA_CHIQ_06 CHA_CHIQ_07 CHA_CHIQ_08
## 24 41 42
## CHA_CHIQ_09 CHA_CHIQ_10 CHA_CHIQ_11
## 32 19 39
## CHA_CHIQ_12 CHA_CHIQ_13 CHA_EASTERN_01
## 29 36 31
## CHA_EASTERN_02 CHA_EASTERN_03 CHA_EASTERN_04
## 27 32 33
## CHA_LEQ_ALU_20 CHA_LEQ_ALU_21 CHA_LEQ_ALU_22
## 23 21 15
## CHA_LEQ_ALU_23 CHA_LEQ_ALU_24 CHA_LEQ_ALU_25
## 26 28 19
## CHA_LEQ_ALU_26 CHA_LEQ_ALU_27 CHA_LEQ_ALU_28
## 24 19 23
## CHA_SEM_AR_LEQ_ALU_01 CHA_EASTERN_05 CHA_EASTERN_06
## 17 32 17
## CHA_EASTERN_07 CHA_CHIQ_14 CHA_CHIQ_15
## 23 17 22
## CHA_CHIQ_16 CHA_CHIQ_17 CHA_CHIQ_18
## 6 18 15
## CHA_CHIQ_19 CHA_HYGRO_01 CHA_HYGRO_02
## 13 7 8
## CHA_HYGRO_03 CHA_FLOO_PALM_01 CHA_SEC_08
## 7 4 43
## CHA_UM_CHA_08 CHA_SEC_09 CHA_SEC_10
## 84 82 43
## CHA_ESPI_01 CHA_UM_CHA_09 CHA_SEC_11
## 91 179 24
## CHA_UM_CHA_10 CHA_SEC_12 CHA_SEC_13
## 122 122 32
## CHA_SEC_14 CHA_SEC_15 CHA_UM_CHA_11
## 132 31 82
## CHA_SEC_16 CERDAO_CER_32 CERDAO_CER_33
## 63 63 10
## CERDAO_CER_34 CERDAO_CER_35 CERDAO_CER_36
## 34 10 70
## CER.S.ST_CER_04
## 27
View(riqueza_sp)
#Calcular a diversidade beta da planilha de presença e ausência (PA)
#O resultado esperado será:
#Diversidade beta total = índice de Sørensen (beta.sor)
#Componente de substituição = índice de Simpson (beta.sim)
#Componente de aninhamento = diferença na riqueza (beta.sne)
resultado_PA <- beta.pair(final_matrix, index.family = "sorensen")
#Montar um data frame com os resultados
dfresultado_PA <- data.frame(round(as.numeric(resultado_PA$beta.sor), 2),
round(as.numeric(resultado_PA$beta.sim), 2),
round(as.numeric(resultado_PA$beta.sne), 2))
colnames(dfresultado_PA) <- c("Sorensen", "Simpson", "Nestedness")
head(dfresultado_PA)
## Sorensen Simpson Nestedness
## 1 0.28 0.04 0.24
## 2 0.62 0.43 0.19
## 3 0.69 0.64 0.05
## 4 0.77 0.70 0.07
## 5 0.46 0.40 0.06
## 6 0.52 0.36 0.16
#Gerarar a matriz de dissimilaridade de Jaccard
dissimilarity_matrix <- vegdist(final_matrix, method = "jaccard")
#Converter a matrix anterior para uma matriz quadrada
square_matrix <- as.matrix(dissimilarity_matrix)
#Calcular a distância média de cada localidade para as demais restantes
distm_matrix <- rowMeans(square_matrix)
#Testar e visualizar se há outliers na matriz de distância
par(mfrow = c(1,2))
hist(distm_matrix, main = "Histograma das Distancias Medias",
xlab = "Distancia Media de Jaccard")
boxplot(distm_matrix, main = "Boxplot (Pontos isolados sao Outliers)",
ylab = "Distancia Media")
str(final_matrix)
## num [1:307, 1:1628] 1 0 1 0 0 1 1 0 1 0 ...
## - attr(*, "dimnames")=List of 2
## ..$ : chr [1:307] "CERDAO_CER_01" "CER.S.ST_CER_01" "CERDAO_CER_02" "CERDAO_CER_03" ...
## ..$ : chr [1:1628] "Acrocomia_aculeata_(Jacq.)_Lodd._ex_Mart." "Agonandra_brasiliensis_Benth._&_Hook.f." "Aspidosperma_subincanum_Mart." "Aspidosperma_tomentosum_Mart." ...