Exercícios - Capitulo 11 - Estimadores de Riqueza
Antônio Carlos do Nascimento Mendonça
2023-01-30
#11.1 Carregue os dados - Cap11_exercicio1 - que está no pacote ecodados. Este conjunto de dados representa a abundância de 50 espécies de besouros coletados em 30 amostras. Calcule os estimadores de riqueza - Chao1 e ACE - e faça um gráfico contendo a riqueza observada e os dois estimadores de riqueza. Qual a sua interpretação sobre o esforço amostral?
## Pacotes
library(iNEXT)
library(devtools)## Carregando pacotes exigidos: usethislibrary(ecodados)
library(ggplot2)
library(vegan)## Carregando pacotes exigidos: permute##
## Attaching package: 'permute'## The following object is masked from 'package:devtools':
##
## check## Carregando pacotes exigidos: lattice## This is vegan 2.6-2library(nlme)
library(dplyr)##
## Attaching package: 'dplyr'## The following object is masked from 'package:nlme':
##
## collapse## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(piecewiseSEM)## Warning: package 'piecewiseSEM' was built under R version 4.2.2##
## This is piecewiseSEM version 2.1.0.
##
##
## Questions or bugs can be addressed to <LefcheckJ@si.edu>.## Dados
dados_Cap11_exercicio1 <- Cap11_exercicio1
head(dados_Cap11_exercicio1)## sp1 sp2 sp3 sp4 sp5 sp6 sp7 sp8 sp9 sp10 sp11 sp12 sp13 sp14 sp15
## Amostra_1 0 0 0 18 1 0 0 0 1 0 0 0 0 0 7
## Amostra_2 2 0 0 7 0 0 0 0 0 0 0 0 0 0 0
## Amostra_3 0 0 1 17 3 0 0 2 2 0 0 0 0 0 4
## Amostra_4 0 0 0 9 2 0 0 1 0 0 0 0 0 0 0
## Amostra_5 0 1 1 9 8 0 1 4 0 0 0 0 0 0 0
## Amostra_6 2 0 1 8 11 0 0 0 0 0 0 1 0 0 0
## sp16 sp17 sp18 sp19 sp20 sp21 sp22 sp23 sp24 sp25 sp26 sp27 sp28 sp29
## Amostra_1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## Amostra_2 0 0 0 0 0 0 0 0 0 4 2 0 0 0
## Amostra_3 0 1 0 0 0 0 1 0 0 3 1 0 8 0
## Amostra_4 0 0 0 0 0 0 0 0 0 1 0 0 2 0
## Amostra_5 0 0 0 0 0 0 0 0 0 2 1 1 1 0
## Amostra_6 0 0 0 0 1 0 1 0 1 0 0 0 0 0
## sp30 sp31 sp32 sp33 sp34 sp35 sp36 sp37 sp38 sp39 sp40 sp41 sp42 sp43
## Amostra_1 0 0 3 0 0 38 0 0 0 0 0 0 0 0
## Amostra_2 0 0 0 0 0 7 0 0 0 1 1 0 0 0
## Amostra_3 0 3 2 0 0 30 0 0 0 0 0 0 0 0
## Amostra_4 0 0 0 0 0 6 0 1 0 3 0 0 0 0
## Amostra_5 0 0 0 0 1 9 0 0 0 2 0 0 0 0
## Amostra_6 0 0 0 0 0 7 1 0 0 0 0 0 0 0
## sp44 sp45 sp46 sp47 sp48 sp49 sp50
## Amostra_1 0 0 0 0 0 2 0
## Amostra_2 0 0 0 0 0 0 0
## Amostra_3 0 0 0 0 0 0 0
## Amostra_4 0 0 0 0 0 3 1
## Amostra_5 0 0 0 0 0 0 0
## Amostra_6 0 0 0 0 0 0 0Cálculo do estimador de riqueza - Chao1
est_chao1 <- estaccumR(dados_Cap11_exercicio1, permutations = 100)
summary(est_chao1, display = "chao")## $chao
## N Chao 2.5% 97.5% Std.Dev
## Amostra_19 1 17.66275 7.333333 52.00000 9.256830
## Amostra_11 2 26.56590 15.095000 47.76250 9.485774
## Amostra_17 3 31.93790 20.291071 63.00000 11.150712
## Amostra_4 4 33.83035 22.722143 54.52500 7.848929
## Amostra_5 5 36.73261 24.106250 62.00000 9.471023
## Amostra_26 6 39.31002 24.892857 66.52500 9.942017
## Amostra_12 7 42.85705 27.618750 77.46250 12.859738
## Amostra_14 8 45.16965 28.868750 72.18125 12.003300
## Amostra_20 9 46.84459 31.500000 71.05000 12.046004
## Amostra_9 10 49.80979 32.596875 77.80000 15.001312
## Amostra_7 11 52.60299 33.530000 90.55000 15.933787
## Amostra_1 12 55.25254 35.500000 97.01250 15.265231
## Amostra_24 13 59.03005 36.618750 122.53750 23.423126
## Amostra_16 14 60.47538 37.106250 117.02500 17.244390
## Amostra_8 15 65.71382 39.223750 116.78750 22.315632
## Amostra_15 16 68.32893 40.285000 130.52500 23.089787
## Amostra_23 17 74.15586 43.950000 170.45000 29.987081
## Amostra_27 18 77.80712 46.000000 179.85000 31.581034
## Amostra_3 19 81.07162 51.428571 206.45000 34.837923
## Amostra_21 20 81.70574 50.425000 170.92500 32.503576
## Amostra_25 21 83.01479 50.079167 179.00000 31.145668
## Amostra_28 22 82.95100 53.095000 166.65000 29.143576
## Amostra_13 23 81.76183 54.647857 137.01250 23.934211
## Amostra_22 24 80.88017 57.237500 131.07500 22.476575
## Amostra_30 25 82.15726 58.000000 120.31250 18.079593
## Amostra_29 26 82.01276 60.166667 123.20000 16.125927
## Amostra_18 27 80.62783 60.166667 102.50000 13.620504
## Amostra_2 28 79.73790 63.500000 102.50000 11.134429
## Amostra_10 29 79.68833 67.000000 102.50000 8.775593
## Amostra_6 30 78.50000 78.500000 78.50000 0.000000
##
## attr(,"class")
## [1] "summary.poolaccum"## Preparando os dados para fazer o gráfico
resultados <- summary(est_chao1, display = c("S", "chao"))
res_chao <- cbind(resultados$chao[, 1:4], resultados$S[, 2:4])
res_chao <- as.data.frame(res_chao)
colnames(res_chao) <- c("Amostras", "Chao", "C_inferior", "C_superior",
"Riqueza", "R_inferior", "R_superior")
## Gráfico
ggplot(res_chao, aes(y = Riqueza, x = Amostras)) +
geom_point(aes(y = Chao, x = Amostras + 0.1), size = 4,
color = "darkorange", alpha = 0.7) +
geom_point(aes(y = Riqueza, x = Amostras), size = 4,
color = "cyan4", alpha = 0.7) +
geom_point(y = 27, x = 18, size = 4, color = "darkorange", alpha = 0.7) +
geom_point(y = 15, x = 18, size = 4, color = "cyan4", alpha = 0.7) +
geom_label(y = 27, x = 24, label = "Riqueza estimada - Chao 1") +
geom_label(y = 15, x = 23, label = "Riqueza observada") +
geom_line(aes(y = Chao, x = Amostras), color = "darkorange") +
geom_line(aes(y = Riqueza, x = Amostras), color = "cyan4") +
geom_linerange(aes(ymin = C_inferior, ymax = C_superior,
x = Amostras + 0.1), color = "darkorange") +
geom_linerange(aes(ymin = R_inferior, ymax = R_superior,
x = Amostras), color = "cyan4") +
scale_x_continuous(limits = c(1, 30), breaks = seq(1, 30, 1)) +
labs (x = "Número de amostras", y = "Riqueza estimada - Chao 1") +
tema_livro()## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 rows containing missing values (geom_segment).
#Interpretação dos resultados
Com base no número de espécies raras (singletons e doubletons), o estimador Chao 1 indicou a possibilidade de encontrarmos mais trinta e quatro espécies, caso o esforço amostral fosse maior e não estima tendência de estabilização da curva em uma assíntota
Cálculo do estimador de riqueza - ACE
## Análise
est_ace <- estaccumR(dados_Cap11_exercicio1, permutations = 100)
summary(est_ace, display = "ace")## $ace
## N ACE 2.5% 97.5% Std.Dev
## Amostra_6 1 23.00309 9.202797 57.84375 12.292898
## Amostra_20 2 30.20008 16.415936 57.25210 10.873999
## Amostra_18 3 34.67280 20.852239 61.75696 11.349790
## Amostra_4 4 37.24352 23.122557 60.46951 10.102926
## Amostra_10 5 38.80151 25.508697 57.06055 8.412430
## Amostra_29 6 40.36547 26.953680 60.15565 8.637590
## Amostra_15 7 42.78275 29.000625 63.86624 9.171576
## Amostra_7 8 44.24025 30.981634 63.50177 8.883794
## Amostra_23 9 44.87618 32.543983 63.25846 8.532278
## Amostra_1 10 45.54128 32.400963 60.27710 7.923353
## Amostra_12 11 47.00584 33.761033 61.54869 7.914637
## Amostra_13 12 48.40302 36.939446 60.83055 6.946311
## Amostra_5 13 49.80590 37.703859 64.37006 7.252190
## Amostra_16 14 51.44375 39.055540 67.52435 7.595130
## Amostra_27 15 52.75818 42.253980 69.72344 7.374545
## Amostra_28 16 54.52894 42.964711 69.89677 7.215556
## Amostra_24 17 55.99964 44.147478 74.64423 7.648846
## Amostra_11 18 58.26749 45.475459 78.05493 8.450264
## Amostra_8 19 60.26916 46.961040 78.66354 8.014234
## Amostra_21 20 62.01850 49.555628 78.76174 7.663162
## Amostra_2 21 63.95507 52.772392 81.92788 7.664199
## Amostra_19 22 66.09825 52.077497 82.55826 7.718167
## Amostra_25 23 69.19620 55.239339 82.81007 7.666600
## Amostra_17 24 70.87717 57.248908 82.77164 7.187296
## Amostra_9 25 72.92841 56.896208 85.97251 7.672992
## Amostra_3 26 74.61950 60.488545 91.80936 8.022045
## Amostra_14 27 76.52908 64.859025 90.08911 7.013733
## Amostra_30 28 79.14173 67.045672 90.96613 6.274410
## Amostra_22 29 81.90901 72.329125 90.05634 5.225274
## Amostra_26 30 84.82758 84.827576 84.82758 0.000000
##
## attr(,"class")
## [1] "summary.poolaccum"## Preparando os dados para fazer o gráfico
resultados_ace <- summary(est_ace, display = c("S", "ace"))
res_ace <- cbind(resultados_ace$ace[, 1:4], resultados_ace$S[, 2:4])
res_ace <- as.data.frame(res_ace)
colnames(res_ace) <- c("Amostras", "ACE", "ACE_inferior", "ACE_superior",
"Riqueza", "R_inferior", "R_superior")
## Gráfico
ggplot(res_ace, aes(y = Riqueza, x = Amostras)) +
geom_point(aes(y = ACE, x = Amostras + 0.1), size = 4,
color = "darkorange", alpha = 0.7) +
geom_point(aes(y = Riqueza, x = Amostras), size = 4,
color = "cyan4", alpha = 0.7) +
geom_point(y = 27, x = 19, size = 4, color = "darkorange", alpha = 0.7) +
geom_point(y = 20, x = 19, size = 4, color = "cyan4", alpha = 0.7) +
geom_label(y = 27, x = 25, label = "Riqueza estimada - ACE") +
geom_label(y = 20, x = 24, label = "Riqueza observada") +
geom_line(aes(y = ACE, x = Amostras), color = "darkorange") +
geom_line(aes(y = Riqueza, x = Amostras), color = "cyan4") +
geom_linerange(aes(ymin = ACE_inferior, ymax = ACE_superior,
x = Amostras + 0.1), color = "darkorange") +
geom_linerange(aes(ymin = R_inferior, ymax = R_superior,
x = Amostras), color = "cyan4") +
scale_x_continuous(limits = c(1, 30), breaks = seq(1, 30, 1)) +
labs(x = "Número de amostras", y = "Riqueza estimada - ACE") +
tema_livro()## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 rows containing missing values (geom_segment).
#Interpretação dos resultados
Com base no número de espécies raras (abundância menor que 10 indivíduos - default), o estimador ACE indica a possibilidade de encontrarmos mais trinta e quatro espécies, caso o esforço amostral fosse maior e não estimou tendência de estabilização da curva em uma assíntota.
#11.2 Utilize o mesmo conjunto de dados do exercício anterior. Calcule os estimadores de riqueza - Jackknife 1 e bootstrap. Faça um gráfico contendo a riqueza observada e os dois estimadores de riqueza. Qual a sua interpretação sobre o esforço amostral? Compare com os resultados do exercício anterior que utilizam estimadores baseados na abundância das espécies.
Análise
Cálculo do estimador de riqueza - Jackknife 1
## Análise
est_jack1 <- poolaccum(dados_Cap11_exercicio1, permutations = 100)
summary(est_jack1, display = "jack1")## $jack1
## N Jackknife 1 2.5% 97.5% Std.Dev
## [1,] 3 28.98000 21.82500 37.00000 4.153837
## [2,] 4 33.27250 25.72500 40.76250 3.737284
## [3,] 5 36.18600 29.20000 42.90500 3.934643
## [4,] 6 38.80167 32.50000 46.25417 3.973869
## [5,] 7 41.03143 33.92500 48.16429 3.990020
## [6,] 8 43.22625 36.41563 51.87500 4.228917
## [7,] 9 44.51556 36.59167 52.68333 4.572087
## [8,] 10 45.53600 37.30000 54.30500 4.736779
## [9,] 11 47.15727 37.40682 56.36364 4.698167
## [10,] 12 48.41167 39.28750 57.06458 4.764470
## [11,] 13 49.23846 39.40192 58.18077 4.769137
## [12,] 14 50.56286 41.46250 60.20179 5.032622
## [13,] 15 51.78533 44.40000 61.29000 4.759317
## [14,] 16 52.93187 43.94531 60.87500 4.763969
## [15,] 17 54.24353 45.42353 62.46618 4.796936
## [16,] 18 55.38111 45.92778 64.88889 4.792320
## [17,] 19 56.61947 46.52632 64.49737 4.656486
## [18,] 20 57.74800 47.95125 66.95250 4.850356
## [19,] 21 58.83810 48.97619 67.54762 4.977551
## [20,] 22 59.85591 50.95341 68.54659 4.718236
## [21,] 23 61.05783 50.54239 68.60978 4.614832
## [22,] 24 62.17958 53.95521 69.60833 4.103583
## [23,] 25 63.42840 54.97600 70.16000 4.069788
## [24,] 26 64.66615 57.91827 70.19231 3.401600
## [25,] 27 65.72704 59.91944 70.22222 2.884020
## [26,] 28 66.94357 62.42857 70.25000 2.457625
## [27,] 29 67.64034 64.41379 70.27586 1.889546
## [28,] 30 68.36667 68.36667 68.36667 0.000000
##
## attr(,"class")
## [1] "summary.poolaccum"## Preparando os dados para fazer o gráfico
resultados_jack1 <- summary(est_jack1, display = c("S", "jack1"))
res_jack1 <- cbind(resultados_jack1$jack1[, 1:4], resultados_jack1$S[, 2:4])
res_jack1 <- as.data.frame(res_jack1)
colnames(res_jack1) <- c("Amostras", "JACK1", "JACK1_inferior", "JACK1_superior",
"Riqueza", "R_inferior", "R_superior")
## Gráfico
ggplot(res_jack1, aes(y = Riqueza, x = Amostras)) +
geom_point(aes(y = JACK1, x = Amostras + 0.1), size = 4,
color = "darkorange", alpha = 0.7) +
geom_point(aes(y = Riqueza, x = Amostras), size = 4,
color = "cyan4", alpha = 0.7) +
geom_point(y = 30, x = 18, size = 4, color = "darkorange", alpha = 0.7) +
geom_point(y = 23, x = 18, size = 4, color = "cyan4", alpha = 0.7) +
geom_label(y = 30, x = 25, label = "Riqueza estimada - Jackknife 1") +
geom_label(y = 23, x = 23, label = "Riqueza observada") +
geom_line(aes(y = JACK1, x = Amostras), color = "darkorange") +
geom_line(aes(y = Riqueza, x = Amostras), color = "cyan4") +
geom_linerange(aes(ymin = JACK1_inferior, ymax = JACK1_superior,
x = Amostras + 0.1), color = "darkorange") +
geom_linerange(aes(ymin = R_inferior, ymax = R_superior,
x = Amostras), color = "cyan4") +
scale_x_continuous(limits = c(1, 30), breaks = seq(1, 30, 1)) +
labs(x = "Número de amostras", y = "Riqueza estimada - Jackknife 1") +
tema_livro()## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 rows containing missing values (geom_segment).
Interpretação dos resultados
Com base no número de espécies raras, o estimador Jackknife 1 estimou possibilidade de encontrarmos mais dezoito espécies, caso o esforço amostral fosse maior e não estimou tendência de estabilização da curva em uma assíntota.
##Cálculo do estimador de riqueza - Bootstrap.
## Análise
est_boot <- poolaccum(dados_Cap11_exercicio1, permutations = 100)
summary(est_boot, display = "boot")## $boot
## N Bootstrap 2.5% 97.5% Std.Dev
## [1,] 3 25.20667 20.38796 32.65463 3.197963
## [2,] 4 28.77977 23.13613 34.85908 3.414053
## [3,] 5 31.74015 26.31829 38.12750 3.083115
## [4,] 6 34.01219 28.59613 39.70312 3.048885
## [5,] 7 35.70686 29.39662 42.82751 3.338443
## [6,] 8 37.19898 30.68968 43.50737 3.247133
## [7,] 9 38.27121 32.11289 44.73793 3.229169
## [8,] 10 39.53516 32.92050 46.45817 3.358871
## [9,] 11 40.55637 34.58530 47.48214 3.338014
## [10,] 12 41.58096 35.46552 48.38686 3.394945
## [11,] 13 42.64600 36.34463 49.36374 3.295806
## [12,] 14 43.65039 37.06502 50.19546 3.280122
## [13,] 15 44.54658 38.20861 50.49432 3.147611
## [14,] 16 45.47768 39.23246 50.92845 3.014081
## [15,] 17 46.50108 40.50869 52.31463 3.044959
## [16,] 18 47.42861 41.34898 52.41223 3.100862
## [17,] 19 48.32933 42.90684 53.62822 2.948990
## [18,] 20 49.20643 43.28912 54.02016 2.822500
## [19,] 21 49.93582 44.21846 54.81146 2.665284
## [20,] 22 50.78819 44.94312 55.57678 2.750594
## [21,] 23 51.57938 46.98759 56.70861 2.668254
## [22,] 24 52.43851 47.48114 56.89705 2.511289
## [23,] 25 53.29411 48.34594 57.78879 2.237423
## [24,] 26 54.16320 49.32094 58.01045 2.232683
## [25,] 27 55.04518 52.08408 58.00550 1.869151
## [26,] 28 55.76167 52.09470 57.91311 1.649001
## [27,] 29 56.75245 54.80833 58.01775 1.124164
## [28,] 30 57.53366 57.53366 57.53366 0.000000
##
## attr(,"class")
## [1] "summary.poolaccum"## Preparando os dados para fazer o gráfico
resultados_boot <- summary(est_boot, display = c("S", "boot"))
res_boot <- cbind(resultados_boot$boot[,1:4], resultados_boot$S[,2:4])
res_boot <- as.data.frame(res_boot)
colnames(res_boot) <- c("Amostras", "BOOT", "BOOT_inferior", "BOOT_superior",
"Riqueza", "R_inferior", "R_superior")
## Gráfico
ggplot(res_boot, aes(y = Riqueza, x = Amostras)) +
geom_point(aes(y = BOOT, x = Amostras + 0.1), size = 4,
color = "darkorange", alpha = 0.7) +
geom_point(aes(y = Riqueza, x = Amostras), size = 4,
color = "cyan4", alpha = 0.7) +
geom_point(y = 30, x = 17, size = 4, color = "darkorange", alpha = 0.7) +
geom_point(y = 25, x = 17, size = 4, color = "cyan4", alpha = 0.7) +
geom_label(y = 30, x = 24, label = "Riqueza estimada - Bootstrap") +
geom_label(y = 25, x = 22, label = "Riqueza observada") +
geom_line(aes(y = BOOT, x = Amostras), color = "darkorange") +
geom_line(aes(y = Riqueza, x = Amostras), color = "cyan4") +
geom_linerange(aes(ymin = BOOT_inferior, ymax = BOOT_superior,
x = Amostras + 0.1), color = "darkorange") +
geom_linerange(aes(ymin = R_inferior, ymax = R_superior,
x = Amostras), color = "cyan4") +
scale_x_continuous(limits = c(1, 30), breaks = seq(1, 30, 1)) +
labs(x = "Número de amostras", y = "Riqueza estimada - Bootstrap") +
tema_livro()## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 rows containing missing values (geom_segment).
#Interpretação dos resultados
Com base na frequência de ocorrência das espécies, o estimador bootstrap estimou a possibilidade de encontrarmos mais sete espécies, caso o esforço amostral fosse maior e não estimou tendência de estabilização da curva em uma assíntota.
O esforço amostral parece gerar menos retorno quando calculada a frequencia de especies estimadas com Jackknife 1 e Bootstrap ao compararmos os resultados por serem menores que o exercício anterior.