RESGATE VIA R MARKDOWN
Hernani Ramos
28 de abril de 2018
###########################
library("tinytex", lib.loc="~/R/win-library/3.4") #PARA IMPRIMIR RELATORIO EM PDF## Warning: package 'tinytex' was built under R version 3.4.4###########################
########################################################
###Começamos estabelecendo uma pasta, lendo os dados e conferindo:
setwd("c:/R") #É o endereço no Dell da Adriana
#Se preferir, troque o comando por: setwd(choose.dir())
#E escolha a pasta onde estão todos os arquivos que acompanham este script
dados<-read.table("Exemplo_Dell_R_Adriana.txt", h=T, row.names=1)
attach(dados)
#Conferindo rapidamente os dados:
summary(dados)## Estado Latitude Longitude Temperatura
## Alagoas :4 Min. :-12.100 Min. :-38.40 Min. :18.20
## Bahia :4 1st Qu.:-10.523 1st Qu.:-37.78 1st Qu.:22.07
## Paraíba :2 Median : -9.500 Median :-37.41 Median :24.50
## Pernambuco:6 Mean : -9.493 Mean :-37.27 Mean :25.13
## Sergipe :4 3rd Qu.: -7.965 3rd Qu.:-36.87 3rd Qu.:28.18
## Max. : -7.630 Max. :-35.82 Max. :32.20
## Cobertura Luminosidade Frutos Profundidade
## Min. :27.0 Min. :53.00 Min. : 5.00 Min. :32.00
## 1st Qu.:39.0 1st Qu.:61.25 1st Qu.:16.50 1st Qu.:40.00
## Median :50.5 Median :69.75 Median :25.50 Median :48.50
## Mean :52.8 Mean :68.90 Mean :25.55 Mean :53.75
## 3rd Qu.:70.0 3rd Qu.:77.75 3rd Qu.:35.25 3rd Qu.:65.00
## Max. :78.0 Max. :84.00 Max. :47.00 Max. :87.00
## sp.1 sp.2 sp.3 sp.4 sp.5
## Min. :0 Min. :0.00 Min. : 1.0 Min. :0.00 Min. :0.00
## 1st Qu.:4 1st Qu.:0.00 1st Qu.: 4.0 1st Qu.:1.75 1st Qu.:0.00
## Median :4 Median :1.00 Median : 5.0 Median :3.00 Median :1.00
## Mean :4 Mean :1.15 Mean : 5.6 Mean :3.45 Mean :1.50
## 3rd Qu.:5 3rd Qu.:1.25 3rd Qu.: 7.0 3rd Qu.:5.25 3rd Qu.:2.25
## Max. :6 Max. :5.00 Max. :10.0 Max. :8.00 Max. :5.00
## sp.6 sp.7 sp.8 sp.9
## Min. : 0.00 Min. :0.0 Min. :0.0 Min. :0.00
## 1st Qu.: 4.00 1st Qu.:0.0 1st Qu.:2.0 1st Qu.:0.00
## Median : 8.50 Median :2.0 Median :3.0 Median :2.00
## Mean : 19.80 Mean :1.9 Mean :3.2 Mean :1.95
## 3rd Qu.: 23.25 3rd Qu.:3.0 3rd Qu.:5.0 3rd Qu.:3.00
## Max. :107.00 Max. :6.0 Max. :7.0 Max. :6.00
## sp.10 sp.11 sp.12 sp.13
## Min. :0.00 Min. : 0.0 Min. :0.0 Min. :0.00
## 1st Qu.:0.75 1st Qu.: 0.0 1st Qu.:0.0 1st Qu.:1.75
## Median :1.00 Median : 2.0 Median :1.5 Median :3.50
## Mean :1.55 Mean : 2.8 Mean :2.1 Mean :3.20
## 3rd Qu.:2.25 3rd Qu.: 5.0 3rd Qu.:3.0 3rd Qu.:5.00
## Max. :6.00 Max. :11.0 Max. :9.0 Max. :7.00
## sp.14 sp.15 sp.16 sp.17
## Min. : 0.00 Min. :0.00 Min. :0.00 Min. :0.0
## 1st Qu.: 4.00 1st Qu.:0.00 1st Qu.:0.00 1st Qu.:0.0
## Median : 6.00 Median :1.00 Median :0.00 Median :0.0
## Mean : 6.25 Mean :1.40 Mean :0.75 Mean :0.6
## 3rd Qu.: 8.25 3rd Qu.:2.25 3rd Qu.:1.00 3rd Qu.:0.0
## Max. :13.00 Max. :5.00 Max. :4.00 Max. :6.0
## sp.18 sp.19 sp.20 sp.21 sp.22
## Min. :0.00 Min. :0.00 Min. :0.0 Min. :0.00 Min. :1.00
## 1st Qu.:0.00 1st Qu.:0.00 1st Qu.:0.0 1st Qu.:2.00 1st Qu.:1.75
## Median :0.00 Median :0.00 Median :0.0 Median :2.50 Median :2.00
## Mean :0.85 Mean :0.75 Mean :0.6 Mean :2.35 Mean :1.90
## 3rd Qu.:2.00 3rd Qu.:1.25 3rd Qu.:0.0 3rd Qu.:3.00 3rd Qu.:2.00
## Max. :3.00 Max. :4.00 Max. :4.0 Max. :4.00 Max. :3.00
## sp.23 sp.24 sp.25 sp.26 sp.27
## Min. :0.0 Min. :0.0 Min. :0.00 Min. :0.00 Min. :0.0
## 1st Qu.:1.0 1st Qu.:0.0 1st Qu.:0.00 1st Qu.:0.00 1st Qu.:0.0
## Median :1.0 Median :0.5 Median :1.00 Median :0.50 Median :0.0
## Mean :0.9 Mean :0.6 Mean :0.85 Mean :0.55 Mean :0.1
## 3rd Qu.:1.0 3rd Qu.:1.0 3rd Qu.:1.00 3rd Qu.:1.00 3rd Qu.:0.0
## Max. :2.0 Max. :2.0 Max. :2.00 Max. :2.00 Max. :1.0
## sp.28 sp.29 sp.30 sp.31 sp.32
## Min. :0.0 Min. :0.0 Min. :0.00 Min. :0.00 Min. :0.0
## 1st Qu.:0.0 1st Qu.:0.0 1st Qu.:0.00 1st Qu.:0.00 1st Qu.:0.0
## Median :0.0 Median :0.0 Median :0.00 Median :1.00 Median :0.0
## Mean :0.2 Mean :0.2 Mean :0.15 Mean :0.65 Mean :0.2
## 3rd Qu.:0.0 3rd Qu.:0.0 3rd Qu.:0.00 3rd Qu.:1.00 3rd Qu.:0.0
## Max. :2.0 Max. :3.0 Max. :2.00 Max. :1.00 Max. :1.0
## sp.33 sp.34 sp.35 sp.36 sp.37
## Min. :0.00 Min. :0.0 Min. :0.00 Min. :0.00 Min. :0.00
## 1st Qu.:0.00 1st Qu.:0.0 1st Qu.:0.00 1st Qu.:0.00 1st Qu.:0.00
## Median :0.00 Median :0.0 Median :0.00 Median :0.00 Median :0.00
## Mean :0.15 Mean :0.1 Mean :0.25 Mean :0.15 Mean :0.05
## 3rd Qu.:0.00 3rd Qu.:0.0 3rd Qu.:0.00 3rd Qu.:0.00 3rd Qu.:0.00
## Max. :1.00 Max. :1.0 Max. :3.00 Max. :2.00 Max. :1.00
## sp.38 sp.39 sp.40 sp.41
## Min. :0.0 Min. :0.0 Min. :0.0 Min. :0.00
## 1st Qu.:0.0 1st Qu.:0.0 1st Qu.:0.0 1st Qu.:0.00
## Median :0.0 Median :0.0 Median :0.0 Median :0.00
## Mean :0.8 Mean :0.3 Mean :0.1 Mean :0.65
## 3rd Qu.:0.0 3rd Qu.:0.0 3rd Qu.:0.0 3rd Qu.:0.00
## Max. :7.0 Max. :3.0 Max. :1.0 Max. :7.00#Conferindo dados de apenas algumas colunas.
summary(dados[,1:8])## Estado Latitude Longitude Temperatura
## Alagoas :4 Min. :-12.100 Min. :-38.40 Min. :18.20
## Bahia :4 1st Qu.:-10.523 1st Qu.:-37.78 1st Qu.:22.07
## Paraíba :2 Median : -9.500 Median :-37.41 Median :24.50
## Pernambuco:6 Mean : -9.493 Mean :-37.27 Mean :25.13
## Sergipe :4 3rd Qu.: -7.965 3rd Qu.:-36.87 3rd Qu.:28.18
## Max. : -7.630 Max. :-35.82 Max. :32.20
## Cobertura Luminosidade Frutos Profundidade
## Min. :27.0 Min. :53.00 Min. : 5.00 Min. :32.00
## 1st Qu.:39.0 1st Qu.:61.25 1st Qu.:16.50 1st Qu.:40.00
## Median :50.5 Median :69.75 Median :25.50 Median :48.50
## Mean :52.8 Mean :68.90 Mean :25.55 Mean :53.75
## 3rd Qu.:70.0 3rd Qu.:77.75 3rd Qu.:35.25 3rd Qu.:65.00
## Max. :78.0 Max. :84.00 Max. :47.00 Max. :87.00#Confereindo dados de outras colunas:
summary(dados[,41:49])## sp.33 sp.34 sp.35 sp.36 sp.37
## Min. :0.00 Min. :0.0 Min. :0.00 Min. :0.00 Min. :0.00
## 1st Qu.:0.00 1st Qu.:0.0 1st Qu.:0.00 1st Qu.:0.00 1st Qu.:0.00
## Median :0.00 Median :0.0 Median :0.00 Median :0.00 Median :0.00
## Mean :0.15 Mean :0.1 Mean :0.25 Mean :0.15 Mean :0.05
## 3rd Qu.:0.00 3rd Qu.:0.0 3rd Qu.:0.00 3rd Qu.:0.00 3rd Qu.:0.00
## Max. :1.00 Max. :1.0 Max. :3.00 Max. :2.00 Max. :1.00
## sp.38 sp.39 sp.40 sp.41
## Min. :0.0 Min. :0.0 Min. :0.0 Min. :0.00
## 1st Qu.:0.0 1st Qu.:0.0 1st Qu.:0.0 1st Qu.:0.00
## Median :0.0 Median :0.0 Median :0.0 Median :0.00
## Mean :0.8 Mean :0.3 Mean :0.1 Mean :0.65
## 3rd Qu.:0.0 3rd Qu.:0.0 3rd Qu.:0.0 3rd Qu.:0.00
## Max. :7.0 Max. :3.0 Max. :1.0 Max. :7.00##########################################################
###Iniciando com alguns calculos simples:
#Média e resumo de uma variável:
mean(Temperatura)## [1] 25.13summary(Temperatura)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 18.20 22.07 24.50 25.13 28.18 32.20#Avaliar dados por meio de histograma:
hist(Temperatura, las=1, main=NULL, ylab="Frequ?ncia", xlab="Temperaturas",
col="gray", xlim=c(15, 35))
# Incluir linha da média no gráfico
abline(v=mean(Temperatura), lty=2, col="red")
#Avaliar dados por estado com boxplot
tapply(Temperatura, Estado, mean)## Alagoas Bahia Paraíba Pernambuco Sergipe
## 22.47500 25.62500 29.45000 28.26667 20.42500boxplot(Temperatura~Estado, las=1, ylab="Temperatura", xlab="Estados")
#####################################################################
###Avaliar outros gráficos e a relação entre variáveis
#Um gráfico simples de plot:
plot(Cobertura, Temperatura, pch=16, las=1)
#Analisar a correlação entre dois dados:
cor(Cobertura, Temperatura)## [1] -0.8350364#Diferenciar com cores os dados do gráfico
plot(Cobertura, Temperatura, pch=16, las=1,
col=c("red","green3","blue", "gray", "yellow")[Estado])
#Inclusão de legenda no gráfico
legend("topright", legend=levels(Estado), pch=16,
col=c("red","green3","blue", "gray", "yellow"))
#Como fazer uma análise rápida de correlação entre dois dados
modelo1<-lm(Temperatura~Cobertura)
summary(modelo1)##
## Call:
## lm(formula = Temperatura ~ Cobertura)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.3524 -1.4926 -0.2025 1.1415 4.0627
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 36.0099 1.7696 20.349 7.15e-14 ***
## Cobertura -0.2061 0.0320 -6.439 4.64e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.352 on 18 degrees of freedom
## Multiple R-squared: 0.6973, Adjusted R-squared: 0.6805
## F-statistic: 41.46 on 1 and 18 DF, p-value: 4.644e-06abline(modelo1, col="red", lty=2)
################################################################
###Analisar muitas variáveis de uma vez:
cor(dados[,4:8])## Temperatura Cobertura Luminosidade Frutos Profundidade
## Temperatura 1.0000000 -0.8350364 0.8570506 0.2897747 -0.5719227
## Cobertura -0.8350364 1.0000000 -0.9608754 -0.2642353 0.8381537
## Luminosidade 0.8570506 -0.9608754 1.0000000 0.3411597 -0.7304425
## Frutos 0.2897747 -0.2642353 0.3411597 1.0000000 0.1211804
## Profundidade -0.5719227 0.8381537 -0.7304425 0.1211804 1.0000000#Plotar uma matriz de graficos de correlações:
pairs(dados[,4:8])
#inserir mais informações na matriz de gráficos:
source("panel.hist.R")
source("panel.cor.R")
pairs(dados[,4:8], diag.panel=panel.hist, upper.panel=panel.cor)
################################################################
###Analizando abundâncias:
#######
##Abundância de apenas uma espécie:
summary(sp.1)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0 4 4 4 5 6tabela.sp1<-table(sp.1)
tabela.sp1## sp.1
## 0 1 3 4 5 6
## 2 1 1 8 4 4barplot(tabela.sp1, las=1, ylab="Frequência", xlab="Abundância")
#######
##Avaliar a abundância de todas as espécies de uma localidade:
amostra1<-dados[1, 9:49]
amostra1## sp.1 sp.2 sp.3 sp.4 sp.5 sp.6 sp.7 sp.8 sp.9 sp.10 sp.11 sp.12 sp.13
## AL1 4 0 7 2 0 0 3 2 1 1 0 1 3
## sp.14 sp.15 sp.16 sp.17 sp.18 sp.19 sp.20 sp.21 sp.22 sp.23 sp.24
## AL1 3 0 3 0 0 0 0 0 1 1 1
## sp.25 sp.26 sp.27 sp.28 sp.29 sp.30 sp.31 sp.32 sp.33 sp.34 sp.35
## AL1 0 2 0 1 0 0 1 0 0 0 0
## sp.36 sp.37 sp.38 sp.39 sp.40 sp.41
## AL1 0 0 0 0 0 0tabela.amostra1<-as.matrix(amostra1)
barplot(tabela.amostra1, las=2, ylab="Frequência", xlab="Espécie")
#Organizar a abundancia em ordem de importância
amostra1.ordem<-sort(amostra1, decreasing=T)
amostra1.ordem## sp.3 sp.1 sp.7 sp.13 sp.14 sp.16 sp.4 sp.8 sp.26 sp.9 sp.10 sp.12
## AL1 7 4 3 3 3 3 2 2 2 1 1 1
## sp.22 sp.23 sp.24 sp.28 sp.31 sp.2 sp.5 sp.6 sp.11 sp.15 sp.17 sp.18
## AL1 1 1 1 1 1 0 0 0 0 0 0 0
## sp.19 sp.20 sp.21 sp.25 sp.27 sp.29 sp.30 sp.32 sp.33 sp.34 sp.35
## AL1 0 0 0 0 0 0 0 0 0 0 0
## sp.36 sp.37 sp.38 sp.39 sp.40 sp.41
## AL1 0 0 0 0 0 0tabela.amostra1.ordem<-as.matrix(amostra1.ordem)
barplot(tabela.amostra1.ordem, las=2, ylab="Frequência",xlab="Espécie")
#Remover as espécies ausentes:
barplot(tabela.amostra1.ordem[,1:17], las=2, ylab="Frequência",
xlab="Espécie")
#######
##Ilustrar todas as abundâncias graficamente
abundâncias<-colSums(dados[,9:49])
abundâncias## sp.1 sp.2 sp.3 sp.4 sp.5 sp.6 sp.7 sp.8 sp.9 sp.10 sp.11 sp.12
## 80 23 112 69 30 396 38 64 39 31 56 42
## sp.13 sp.14 sp.15 sp.16 sp.17 sp.18 sp.19 sp.20 sp.21 sp.22 sp.23 sp.24
## 64 125 28 15 12 17 15 12 47 38 18 12
## sp.25 sp.26 sp.27 sp.28 sp.29 sp.30 sp.31 sp.32 sp.33 sp.34 sp.35 sp.36
## 17 11 2 4 4 3 13 4 3 2 5 3
## sp.37 sp.38 sp.39 sp.40 sp.41
## 1 16 6 2 13abundâncias.ordem<-sort(abundâncias, decreasing=T)
barplot(abundâncias.ordem, las=2, ylab="Abundâncias", xlab="Espécies")
#Olhando para estes dados a fundo
dados[,c(1:3, 14)]## Estado Latitude Longitude sp.6
## AL1 Alagoas -9.90 -37.07 0
## SE1 Sergipe -11.12 -37.41 4
## SE2 Sergipe -10.07 -37.46 3
## AL2 Alagoas -9.53 -37.88 9
## AL3 Alagoas -9.12 -37.77 8
## SE3 Sergipe -9.95 -37.56 0
## BA1 Bahia -11.68 -38.40 7
## BA2 Bahia -10.53 -37.88 5
## SE4 Sergipe -10.52 -36.82 1
## BA3 Bahia -12.10 -38.02 5
## BA4 Bahia -11.90 -37.80 27
## PE1 Pernambuco -8.44 -36.07 4
## AL4 Alagoas -9.47 -36.21 107
## PB1 Paraíba -7.87 -36.46 53
## PE2 Pernambuco -7.73 -37.38 22
## PE3 Pernambuco -7.92 -37.77 17
## PE4 Pernambuco -7.63 -37.26 43
## PE5 Pernambuco -8.68 -37.40 14
## PB2 Paraíba -7.73 -36.88 57
## PE6 Pernambuco -7.98 -35.82 10#Vamos supor que os dados estejam corretos. Como resolver? Log!
log.abundâncias<-log10(abundâncias.ordem+1)
barplot(log.abundâncias, las=2, ylab="Logaritmo das abundâncias",
xlab="Espécies", ylim=c(0,3))
#Vamos melhorar colocando os valores reais nas barras
gráfico<-barplot(log.abundâncias, las=2, ylab="Logaritmo das abundâncias",
xlab="Espécies", ylim=c(0,3))
text(x=gráfico, y=log.abundâncias+0.1, abundâncias.ordem)
##################################################################
###E a biodiversidade?
#ATENÇÃO ATUALIZAR DIRETORIO DOS PACOTES
#SUGESTÃO SUBIR PARA INICIO DO SCRIPT OS COMANDOS DE ABRIR PACOTES
library("vegan", lib.loc="~/R/win-library/3.4") ## Warning: package 'vegan' was built under R version 3.4.4## Loading required package: permute## Warning: package 'permute' was built under R version 3.4.4## Loading required package: lattice## This is vegan 2.5-1
#Se não tiver o pacote, vê no menu "Pacotes" e clique em "Instalar pacotes"
#Escolha os servidor de onde vai baixar e encontre o pacote desejado na lista
#Riqueza de espécies por localidade de coleta:
riqueza<-specnumber(dados[,9:49])
riqueza## AL1 SE1 SE2 AL2 AL3 SE3 BA1 BA2 SE4 BA3 BA4 PE1 AL4 PB1 PE2 PE3 PE4 PE5
## 17 25 21 19 22 17 21 17 17 18 20 15 23 20 19 22 19 25
## PB2 PE6
## 18 19#E a riqueza total por estado:
specnumber(dados[,9:49], group=Estado)## Alagoas Bahia Paraíba Pernambuco Sergipe
## 33 29 23 28 34#E a diversidade? Vamos de Shannon:
diversity(dados[,9:49], index="shannon")## AL1 SE1 SE2 AL2 AL3 SE3 BA1 BA2
## 2.624194 2.976366 2.928375 2.697784 2.758447 2.724256 2.849312 2.647016
## SE4 BA3 BA4 PE1 AL4 PB1 PE2 PE3
## 2.628297 2.711446 2.417293 2.579156 1.894462 2.365457 2.583259 2.822355
## PE4 PE5 PB2 PE6
## 2.197376 2.940925 2.144414 2.681331#Plotar curva de rarefação.
rarefação<-specaccum(dados[,9:49])## Warning in cor(x > 0): o desvio padrão é zerorarefação## Species Accumulation Curve
## Accumulation method: exact
## Call: specaccum(comm = dados[, 9:49])
##
##
## Sites 1.00000 2.00000 3.00000 4.00000 5.00000 6.00000 7.00000
## Richness 19.70000 26.12105 29.54825 31.83633 33.55076 34.91736 36.04172
## sd 2.64764 2.51567 2.43165 2.31283 2.18803 2.05671 1.92348
##
## Sites 8.00000 9.00000 10.00000 11.00000 12.00000 13.00000 14.00000
## Richness 36.98042 37.76836 38.42997 38.98398 39.44561 39.82755 40.14061
## sd 1.79009 1.65715 1.52485 1.39273 1.26059 1.12726 0.99249
##
## Sites 15.00000 16.00000 17.00000 18.00000 19.00000 20
## Richness 40.39396 40.59546 40.75175 40.86842 40.95000 41
## sd 0.85486 0.71348 0.56641 0.40259 0.21794 0plot(rarefação, las=1, xlab="Esforço amostral",
ylab="Riqueza acumulada", ylim=c(10, 45))
#Estimando a riqueza
specpool(dados[,9:49])## Species chao chao.se jack1 jack1.se jack2 boot boot.se
## All 41 41.07917 0.3244698 41.95 0.95 37.73421 42.30226 1.263051
## n
## All 20#######
#Pensando na composição de espécies.
#Vamos criar um cluster vertical com o coeficiente de Jaccard:
matriz.jac<-vegdist(dados[,9:49], method="jac", binary=T)
agrupamento.jac<-hclust(matriz.jac, method="average")
plot(agrupamento.jac, hang=-1)
#Mudando o cluster para a horizontal:
plot(as.dendrogram(agrupamento.jac), horiz=T)
#######
#Representação mediante NMDS
#Iniciando com l?:
nmds.jac<-metaMDS(matriz.jac)## Run 0 stress 0.1603297
## Run 1 stress 0.1804793
## Run 2 stress 0.1843736
## Run 3 stress 0.1607154
## ... Procrustes: rmse 0.01440179 max resid 0.05266457
## Run 4 stress 0.1804793
## Run 5 stress 0.1607154
## ... Procrustes: rmse 0.0144174 max resid 0.0526862
## Run 6 stress 0.1716372
## Run 7 stress 0.1668308
## Run 8 stress 0.1839113
## Run 9 stress 0.1628222
## Run 10 stress 0.1716381
## Run 11 stress 0.161823
## Run 12 stress 0.1854323
## Run 13 stress 0.1745404
## Run 14 stress 0.1773553
## Run 15 stress 0.1673048
## Run 16 stress 0.167373
## Run 17 stress 0.1677702
## Run 18 stress 0.1777048
## Run 19 stress 0.1861626
## Run 20 stress 0.1694352
## *** No convergence -- monoMDS stopping criteria:
## 20: stress ratio > sratmaxplot(nmds.jac, type="t")## species scores not available
#Mudar cores por estado e incluir legenda:
plot(nmds.jac, type="n")## species scores not availablepoints(nmds.jac, pch=16,
col=c("red","green3","blue", "gray", "yellow")[Estado])
legend("bottomleft", legend=levels(Estado), pch=16,
col=c("red","green3","blue", "gray", "yellow"))
#E que tal pontos com polígonos por estado?
plot(nmds.jac, type="n")## species scores not availablepoints(nmds.jac, pch=16)
ordihull(nmds.jac, group=Estado, show="Alagoas", col="red")
ordihull(nmds.jac, group=Estado, show="Bahia", col="green3")
ordihull(nmds.jac, group=Estado, show="Para?ba", col="blue")
ordihull(nmds.jac, group=Estado, show="Pernambuco", col="gray")
ordihull(nmds.jac, group=Estado, show="Sergipe", col="yellow")
legend("bottomleft", legend=levels(Estado), lty=1,
col=c("red","green3","blue", "gray", "yellow"))
#######
#Analisar a relação com as variáveis ambientais:
#Primeiro, uma PCA:
resultado.pca<-rda(dados[,4:8], scale=T)
summary(resultado.pca)##
## Call:
## rda(X = dados[, 4:8], scale = T)
##
## Partitioning of correlations:
## Inertia Proportion
## Total 5 1
## Unconstrained 5 1
##
## Eigenvalues, and their contribution to the correlations
##
## Importance of components:
## PC1 PC2 PC3 PC4 PC5
## Eigenvalue 3.4760 1.1035 0.31331 0.08726 0.019986
## Proportion Explained 0.6952 0.2207 0.06266 0.01745 0.003997
## Cumulative Proportion 0.6952 0.9159 0.97855 0.99600 1.000000
##
## Scaling 2 for species and site scores
## * Species are scaled proportional to eigenvalues
## * Sites are unscaled: weighted dispersion equal on all dimensions
## * General scaling constant of scores: 3.121986
##
##
## Species scores
##
## PC1 PC2 PC3 PC4 PC5
## Temperatura 1.248 -0.14171 0.58134 0.18582 -0.01340
## Cobertura -1.378 -0.07164 0.11098 0.08085 -0.16019
## Luminosidade 1.360 -0.09205 0.02064 -0.28382 -0.09844
## Frutos 0.425 -1.29489 -0.28723 0.09676 -0.01204
## Profundidade -1.135 -0.66378 0.42137 -0.19780 0.05730
##
##
## Site scores (weighted sums of species scores)
##
## PC1 PC2 PC3 PC4 PC5
## AL1 -1.02957 0.255624 0.57253 -0.169273 1.577743
## SE1 -0.44258 0.794137 -1.17293 -0.009534 0.340298
## SE2 -0.96679 -0.001699 -0.55269 1.137916 -0.817939
## AL2 -0.84237 -0.373709 -0.33273 0.659017 -0.330739
## AL3 0.23145 0.619376 -0.66615 0.391840 -0.755281
## SE3 -1.02363 0.816577 0.15611 -2.003702 -0.678930
## BA1 -0.93442 -0.764020 1.34952 0.007348 -0.008326
## BA2 0.17452 0.790036 1.39590 0.429347 -0.374980
## SE4 0.28112 1.350440 -0.59549 0.026945 1.316592
## BA3 0.16423 0.684458 0.07558 1.156419 -0.787000
## BA4 -0.63077 -0.303762 0.25033 -0.207189 -0.365035
## PE1 1.06121 0.291911 0.41558 0.079573 0.617395
## AL4 -0.17658 -1.252884 -0.69737 0.614629 0.722140
## PB1 0.90465 -0.523667 0.75635 0.282218 -0.293837
## PE2 0.79302 -0.030483 0.80155 -0.276849 -0.174969
## PE3 0.92496 -0.521060 -0.71764 -0.659446 -0.115347
## PE4 0.08705 -0.854247 0.04089 0.100571 1.118842
## PE5 0.15994 -1.010178 -0.64287 -1.072086 -0.504148
## PB2 0.48247 -0.380949 -0.40184 -0.271135 -0.398108
## PE6 0.78210 0.414099 -0.03463 -0.216611 -0.088369biplot(resultado.pca)
#Incluir polígonos nesta PCA
ordihull(resultado.pca, group=Estado, show="Alagoas", col="red")
ordihull(resultado.pca, group=Estado, show="Bahia", col="green3")
ordihull(resultado.pca, group=Estado, show="Para?ba", col="blue")
ordihull(resultado.pca, group=Estado, show="Pernambuco", col="gray")
ordihull(resultado.pca, group=Estado, show="Sergipe", col="yellow")
#######
#Relacionando tudo:
nmds.env<-envfit(nmds.jac, dados[,4:8])
nmds.env##
## ***VECTORS
##
## NMDS1 NMDS2 r2 Pr(>r)
## Temperatura 0.77792 0.62836 0.5985 0.001 ***
## Cobertura -0.79494 -0.60668 0.6590 0.001 ***
## Luminosidade 0.86708 0.49817 0.6421 0.001 ***
## Frutos 0.81910 -0.57365 0.7046 0.001 ***
## Profundidade -0.56166 -0.82737 0.3194 0.049 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Permutation: free
## Number of permutations: 999plot(nmds.jac, type="n")## species scores not availablepoints(nmds.jac, pch=16,
col=c("red","green3","blue", "gray", "yellow")[Estado])
legend("bottomleft", legend=levels(Estado), pch=16,
col=c("red","green3","blue", "gray", "yellow"))
plot(nmds.env)
############################################################
#Para entender melhor basta analisar os dados em mapas
library("dismo", lib.loc="~/R/win-library/3.4")## Warning: package 'dismo' was built under R version 3.4.4## Loading required package: raster## Warning: package 'raster' was built under R version 3.4.4## Loading required package: sp## Warning: package 'sp' was built under R version 3.4.4library("maptools", lib.loc="~/R/win-library/3.4")## Warning: package 'maptools' was built under R version 3.4.4## Checking rgeos availability: FALSE
## Note: when rgeos is not available, polygon geometry computations in maptools depend on gpclib,
## which has a restricted licence. It is disabled by default;
## to enable gpclib, type gpclibPermit()#Lidando com mapas no R:
data(wrld_simpl)
plot(wrld_simpl)
#Plotando os dados no mapa:
plot(wrld_simpl, xlim=c(-90, -30), ylim=c(-60, 20),
axes=T, col="gray")
points(dados[,c(3,2)], pch=16, col="red", cex=0.5)
##########################################
# PAREI AQUI PRECISO FAZER UPDATE DE RGDAL
##########################################
# TENTAR INSTALAR VERSÕES POSTERIORES A CORRENTE
#################################################
library("rgdal", lib.loc="~/R/win-library/3.4")## Warning: package 'rgdal' was built under R version 3.4.4## rgdal: version: 1.2-18, (SVN revision 718)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 2.2.3, released 2017/11/20
## Path to GDAL shared files: C:/Users/adriana.colosio/Documents/R/win-library/3.4/rgdal/gdal
## GDAL binary built with GEOS: TRUE
## Loaded PROJ.4 runtime: Rel. 4.9.3, 15 August 2016, [PJ_VERSION: 493]
## Path to PROJ.4 shared files: C:/Users/adriana.colosio/Documents/R/win-library/3.4/rgdal/proj
## Linking to sp version: 1.2-7########################################################################
#pacotes a serem instalados um a um no ubuntu
########################################################################
#library("lattice", lib.loc="/usr/lib/R/library")
#library("sp", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4")
#library("codetools", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4")
#library("iterators", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4")
#library("RUnit", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4") # sugerido
#library("R.methodsS3", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4")
#library("Rcpp", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4")
#library("raster", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4")
#library("methods", lib.loc="/usr/lib/R/library")
#library("foreach", lib.loc="~/R/x86_64-pc-linux-gnu-library/3.4")
#library("grDevices", lib.loc="/usr/lib/R/library")
#library("graphics", lib.loc="/usr/lib/R/library")
#library("stats", lib.loc="/usr/lib/R/library")
#library("utils", lib.loc="/usr/lib/R/library")
########################################################################
#Aplicando zoom e incluir os estados:
brasil<-shapefile("uf_car.shp")## Warning in .local(x, ...): .prj file is missingplot(brasil, xlim=c(-40, -34), ylim=c(-13, -5),
axes=T, col="gray")
points(dados[,c(3,2)], pch=16, col="red")
#Plotar por proporcionalidade de riqueza de espécies
plot(brasil, xlim=c(-40, -34), ylim=c(-13, -5),
axes=T, col="gray")
points(dados[,c(3,2)], pch=16, col="black", cex=0.1)
points(dados[,c(3,2)], pch=1, col="red", cex=riqueza/5)
#Avaliar o efeito do Rio São Francisco separando os dados:
plot(brasil, xlim=c(-40, -34), ylim=c(-13, -5),
axes=T, col="gray")
points(dados[,c(3,2)], pch=16, col="red")
rio<-shapefile("SaoFrancisco.shp")## Warning in .local(x, ...): .prj file is missingplot(rio, add=T, col="blue", lwd=2)
#######
#Retornando ao NMDS, para separar os dados em "sul e norte" do rio São Francisco
#Primeiro passo, criar os dois grupos:
Posição<-ifelse(Latitude<(-10), "Sul", "Norte")
Posição<-as.factor(Posição)
#Plotando novamente para o NMDS:
plot(nmds.jac, type="n")## species scores not availablepoints(nmds.jac, pch=16,
col=c("red","green3")[Posição])
legend("bottomleft", legend=levels(Posição), pch=16,
col=c("red","green3"))
#######
#Vamos um pouco além com o SIG, e trazer alguns gridfiles:
files<-list.files(path=paste(system.file(package="dismo"),
"/ex", sep=""), pattern="grd", full.names=T)
precipitação<-stack(files[2]) #Aqui usamos apenas a precipitação, como exemplo
plot(precipitação, xlim=c(-40, -34), ylim=c(-13, -5))
plot(brasil, xlim=c(-40, -34), ylim=c(-13, -5), add=T)
points(dados[,c(3,2)], pch=16, col="red", cex=1)
#Qual a precipitação nos pontos de coleta?
precipitação.pontos<-extract(precipitação, dados[,c(3,2)])
novos.dados<-cbind(dados, precipitação.pontos)
novos.dados[,c(1,50)]## Estado bio12
## AL1 Alagoas 604
## SE1 Sergipe 1497
## SE2 Sergipe 824
## AL2 Alagoas 569
## AL3 Alagoas 615
## SE3 Sergipe 569
## BA1 Bahia 963
## BA2 Bahia 924
## SE4 Sergipe 1370
## BA3 Bahia 1587
## BA4 Bahia 1385
## PE1 Pernambuco 580
## AL4 Alagoas 1435
## PB1 Paraíba 425
## PE2 Pernambuco 636
## PE3 Pernambuco 748
## PE4 Pernambuco 636
## PE5 Pernambuco 642
## PB2 Paraíba 507
## PE6 Pernambuco 793hist(precipitação.pontos)
Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated results.