Dimensionamento de amostragem via simulação bootstrap
Walmes Zeviani, Ana Paula Dalla Corte e colaboradores
Definições da sessão
#-----------------------------------------------------------------------
# Pacotes.
rm(list = objects())
library(gridExtra)
library(data.table)
library(tidyverse)
#-----------------------------------------------------------------------
# Funções.
# Função para contar os valores únicos de uma variável.
count_levels <- function(x) {
x <- x[!is.na(x)]
length(unique(x))
}
# Índices de diversidade.
shannon <- function(x) {
p <- prop.table(table(x))
-sum(p * log2(p))
}
mcintosh_u <- function(x) {
f <- table(x)
sqrt(sum(f^2))
}Importação do censo
## Parsed with column specification:
## cols(
## .default = col_character(),
## CAP = col_double(),
## Alt = col_double(),
## QF = col_double(),
## Quebrada = col_logical(),
## Classe_Dia = col_double(),
## N_arv = col_double(),
## Sanidade = col_logical(),
## vol = col_double(),
## Patio = col_logical(),
## x_ok = col_double(),
## y_ok = col_double(),
## Placa = col_double(),
## index = col_double()
## )## See spec(...) for full column specifications.## Warning: 47623 parsing failures.
## row col expected actual file
## 1134 Quebrada 1/0/T/F/TRUE/FALSE x 'tb_arvores.csv'
## 1137 Quebrada 1/0/T/F/TRUE/FALSE x 'tb_arvores.csv'
## 1805 Quebrada 1/0/T/F/TRUE/FALSE x 'tb_arvores.csv'
## 2103 Quebrada 1/0/T/F/TRUE/FALSE x 'tb_arvores.csv'
## 2964 Quebrada 1/0/T/F/TRUE/FALSE x 'tb_arvores.csv'
## .... ........ .................. ...... ................
## See problems(...) for more details.## tibble [82,706 × 23] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ UPA : chr [1:82706] NA "UPA-04" "UPA-04" "UPA-04" ...
## $ UT : chr [1:82706] NA "UT-01" "UT-01" "UT-01" ...
## $ Nome_Verna: chr [1:82706] NA "Marupa" "Garapeira" "Castanheira" ...
## $ Nome_Cient: chr [1:82706] NA "Jacaranda copaia" "Apuleia leiocarpa" "Bertholletia excelsa" ...
## $ CAP : num [1:82706] NA 215 295 280 298 269 285 182 220 219 ...
## $ Alt : num [1:82706] NA 14 15 17 16 16 18 18 16 7 ...
## $ QF : num [1:82706] NA 2 1 1 1 1 1 1 1 1 ...
## $ Morta : chr [1:82706] NA NA NA NA ...
## $ Caída : chr [1:82706] NA NA NA NA ...
## $ Quebrada : logi [1:82706] NA NA NA NA NA NA ...
## $ Ocada : chr [1:82706] NA NA NA NA ...
## $ DMC : chr [1:82706] NA "Acima" "Acima" "Abaixo" ...
## $ Classe_Dia: num [1:82706] NA 65 95 85 95 85 95 55 75 65 ...
## $ Seleção : chr [1:82706] NA "Rara" "Explorável" "Proibida de Corte" ...
## $ Categoria : chr [1:82706] NA "Rara" "Corte" "Proibida de Corte" ...
## $ N_arv : num [1:82706] NA 1 2 3 4 5 6 7 8 9 ...
## $ Sanidade : logi [1:82706] NA NA NA NA NA NA ...
## $ vol : num [1:82706] NA 4.14 8.39 7.49 8.58 ...
## $ Patio : logi [1:82706] NA NA NA NA NA NA ...
## $ x_ok : num [1:82706] NA 583158 583106 583053 583008 ...
## $ y_ok : num [1:82706] NA 8969988 8969989 8969997 8970001 ...
## $ Placa : num [1:82706] NA 0 0 0 0 0 0 0 0 0 ...
## $ index : num [1:82706] 0 1 2 3 4 5 6 7 8 9 ...
## - attr(*, "problems")= tibble [47,623 × 5] (S3: tbl_df/tbl/data.frame)
## ..$ row : int [1:47623] 1134 1137 1805 2103 2964 9698 11309 12270 12687 14338 ...
## ..$ col : chr [1:47623] "Quebrada" "Quebrada" "Quebrada" "Quebrada" ...
## ..$ expected: chr [1:47623] "1/0/T/F/TRUE/FALSE" "1/0/T/F/TRUE/FALSE" "1/0/T/F/TRUE/FALSE" "1/0/T/F/TRUE/FALSE" ...
## ..$ actual : chr [1:47623] "x" "x" "x" "x" ...
## ..$ file : chr [1:47623] "'tb_arvores.csv'" "'tb_arvores.csv'" "'tb_arvores.csv'" "'tb_arvores.csv'" ...Cálculo dos valores populacionais
# Determina os valores populacionais.
censo_val <- censo %>%
summarise(Alt_medio = mean(Alt, na.rm = TRUE),
vol_medio = mean(vol, na.rm = TRUE),
vol_sum = sum(vol, na.rm = TRUE),
CAP_medio = mean(CAP, na.rm = TRUE),
n = sum(!is.na(Nome_Verna)),
n_esp = count_levels(Nome_Verna),
margalef = (n_esp - 1)/log(n),
menhinick = n_esp/sqrt(n),
shannon = shannon(Nome_Verna),
mcintosh = (n - mcintosh_u(Nome_Verna))/(n - sqrt(n)))
censo_val## # A tibble: 1 x 10
## Alt_medio vol_medio vol_sum CAP_medio n n_esp margalef menhinick shannon
## <dbl> <dbl> <dbl> <dbl> <int> <int> <dbl> <dbl> <dbl>
## 1 14.3 4.69 387494. 214. 82705 125 11.0 0.435 5.75
## # … with 1 more variable: mcintosh <dbl># ATTENTION: colocar na agregação por amostra bootstrap os valores
# máximos de número de espécies. Os índices devem ser calculados como
# máximo no número de especies determinado em uma amostragem, correto?
str(censo)## tibble [82,706 × 5] (S3: tbl_df/tbl/data.frame)
## $ index : num [1:82706] 0 1 2 3 4 5 6 7 8 9 ...
## $ Nome_Verna: chr [1:82706] NA "Marupa" "Garapeira" "Castanheira" ...
## $ Alt : num [1:82706] NA 14 15 17 16 16 18 18 16 7 ...
## $ CAP : num [1:82706] NA 215 295 280 298 269 285 182 220 219 ...
## $ vol : num [1:82706] NA 4.14 8.39 7.49 8.58 ...
## - attr(*, "problems")= tibble [47,623 × 5] (S3: tbl_df/tbl/data.frame)
## ..$ row : int [1:47623] 1134 1137 1805 2103 2964 9698 11309 12270 12687 14338 ...
## ..$ col : chr [1:47623] "Quebrada" "Quebrada" "Quebrada" "Quebrada" ...
## ..$ expected: chr [1:47623] "1/0/T/F/TRUE/FALSE" "1/0/T/F/TRUE/FALSE" "1/0/T/F/TRUE/FALSE" "1/0/T/F/TRUE/FALSE" ...
## ..$ actual : chr [1:47623] "x" "x" "x" "x" ...
## ..$ file : chr [1:47623] "'tb_arvores.csv'" "'tb_arvores.csv'" "'tb_arvores.csv'" "'tb_arvores.csv'" ...censo_long <- censo %>%
select(Alt, CAP, vol) %>%
rename("Height" = Alt, "Circunference" = CAP, "Volume" = vol) %>%
gather(key = "var", value = "val", everything())
ggplot(censo_long,
mapping = aes(x = val)) +
facet_wrap(facets = ~var, scale = "free", ncol = 1) +
geom_histogram(bins = 25)## Warning: Removed 3 rows containing non-finite values (stat_bin).
ggplot(censo_long,
mapping = aes(x = val)) +
facet_wrap(facets = ~var, scale = "free", ncol = 1) +
stat_ecdf()## Warning: Removed 3 rows containing non-finite values (stat_ecdf).
ggplot(censo_long,
mapping = aes(x = val)) +
facet_wrap(facets = ~var, scale = "free", ncol = 1) +
geom_density()## Warning: Removed 3 rows containing non-finite values (stat_density).
ggplot(censo_long,
mapping = aes(y = val)) +
facet_wrap(facets = ~var, scale = "free", ncol = 1) +
geom_boxplot()## Warning: Removed 3 rows containing non-finite values (stat_boxplot).
# Amplitudes.
xlim <- censo_long %>%
group_by(`var`) %>%
do(a = {
extendrange(.$val)
}) %>%
pull(a)
xlim## [[1]]
## [1] 54.15 908.85
##
## [[2]]
## [1] 2.7 31.3
##
## [[3]]
## [1] -3.688602 82.204219gg1 <-
grid.arrange(
ggplot(data = censo,
mapping = aes(x = CAP)) +
facet_wrap(facets = ~"Circunferece") +
geom_histogram(binwidth = 20,
mapping = aes(y = ..density..)) +
geom_density(bw = 20) +
lims(x = xlim[[1]]) +
labs(x = NULL, y = NULL),
ggplot(data = censo,
mapping = aes(x = CAP)) +
geom_boxplot() +
scale_y_continuous(breaks = NULL) +
lims(x = xlim[[1]]) +
labs(x = NULL) +
theme(plot.margin = margin(l = 1.2, r = 0.2, unit = "cm")),
ncol = 1,
heights = c(0.8, 0.2))## Warning: Removed 1 rows containing non-finite values (stat_bin).## Warning: Removed 1 rows containing non-finite values (stat_density).## Warning: Removed 2 rows containing missing values (geom_bar).## Warning: Removed 1 rows containing non-finite values (stat_boxplot).
gg2 <-
grid.arrange(
ggplot(data = censo,
mapping = aes(x = Alt)) +
facet_wrap(facets = ~"Height") +
geom_histogram(binwidth = 1,
mapping = aes(y = ..density..)) +
geom_density(bw = 1) +
lims(x = xlim[[2]]) +
labs(x = NULL, y = NULL),
ggplot(data = censo,
mapping = aes(x = Alt)) +
geom_boxplot() +
scale_y_continuous(breaks = NULL) +
lims(x = xlim[[2]]) +
labs(x = NULL) +
theme(plot.margin = margin(l = 1.05, r = 0.2, unit = "cm")),
ncol = 1,
heights = c(0.8, 0.2))## Warning: Removed 1 rows containing non-finite values (stat_bin).## Warning: Removed 1 rows containing non-finite values (stat_density).## Warning: Removed 2 rows containing missing values (geom_bar).## Warning: Removed 1 rows containing non-finite values (stat_boxplot).
gg3 <-
grid.arrange(
ggplot(data = censo,
mapping = aes(x = vol)) +
facet_wrap(facets = ~"Volume") +
geom_histogram(binwidth = 1,
mapping = aes(y = ..density..)) +
geom_density(bw = 1) +
lims(x = xlim[[3]]) +
labs(x = NULL, y = NULL),
ggplot(data = censo,
mapping = aes(x = vol)) +
geom_boxplot() +
scale_y_continuous(breaks = NULL) +
lims(x = xlim[[3]]) +
labs(x = NULL) +
theme(plot.margin = margin(l = 1.0, r = 0.2, unit = "cm")),
ncol = 1,
heights = c(0.8, 0.2))## Warning: Removed 1 rows containing non-finite values (stat_bin).## Warning: Removed 1 rows containing non-finite values (stat_density).## Warning: Removed 2 rows containing missing values (geom_bar).## Warning: Removed 1 rows containing non-finite values (stat_boxplot).
Importação dos resultados de simulação
## [1] "B_cond_100x100_1.csv" "B_cond_100x100_11.csv"
## [3] "B_cond_100x100_2.csv" "B_cond_100x100_3.csv"
## [5] "B_cond_100x100_5.csv" "B_cond_100x100_8.csv"
## [7] "B_cond_25x100_1.csv" "B_cond_25x100_11.csv"
## [9] "B_cond_25x100_2.csv" "B_cond_25x100_3.csv"
## [11] "B_cond_25x100_5.csv" "B_cond_25x100_8.csv"
## [13] "B_cond_50x100_1.csv" "B_cond_50x100_11.csv"
## [15] "B_cond_50x100_2.csv" "B_cond_50x100_3.csv"
## [17] "B_cond_50x100_5.csv" "B_cond_50x100_8.csv"
## [19] "B_cond_50x200_1.csv" "B_cond_50x200_11.csv"
## [21] "B_cond_50x200_2.csv" "B_cond_50x200_3.csv"
## [23] "B_cond_50x200_5.csv" "B_cond_50x200_8.csv"
## [25] "B_cond_50x50_1.csv" "B_cond_50x50_11.csv"
## [27] "B_cond_50x50_2.csv" "B_cond_50x50_3.csv"
## [29] "B_cond_50x50_5.csv" "B_cond_50x50_8.csv"
## [31] "B_cond_70.71x70.71_1.csv" "B_cond_70.71x70.71_11.csv"
## [33] "B_cond_70.71x70.71_2.csv" "B_cond_70.71x70.71_3.csv"
## [35] "B_cond_70.71x70.71_5.csv" "B_cond_70.71x70.71_8.csv"fls <- dir(path = "agrega_B_parc",
pattern = "^B_.*\\.csv",
full.names = TRUE,
ignore.case = FALSE)
fls <- lapply(fls, FUN = fread, sep = ";")
boot <- invoke(bind_rows, fls)
str(boot)## Classes 'data.table' and 'data.frame': 108000 obs. of 17 variables:
## $ B : int 1 2 3 4 5 6 7 8 9 10 ...
## $ n_parc : int 85 85 85 85 85 85 85 85 85 85 ...
## $ prop_vazia : num 0 0.0118 0.0118 0 0 ...
## $ n_esp_parc : num 7.54 7.55 7.4 8.56 7.87 ...
## $ mediaw_media_cap: num 213 216 213 213 213 ...
## $ mediaw_media_alt: num 14.3 14.8 14.1 14.6 14.1 ...
## $ media_sum_vol : num 43.2 43.2 41.3 48.1 45.1 ...
## $ desvio_sum_vol : num 25.2 26.9 24.1 26.2 24.5 ...
## $ li_vol : num 37.8 37.4 36.1 42.5 39.8 ...
## $ lu_vol : num 48.7 49 46.5 53.8 50.4 ...
## $ cond : chr "cond_100x100_1" "cond_100x100_1" "cond_100x100_1" "cond_100x100_1" ...
## $ n_esp : int 84 87 91 88 93 86 87 89 89 80 ...
## $ shannon : num 5.57 5.62 5.7 5.63 5.71 ...
## $ margalef : num 12.4 12.9 13.5 12.8 13.7 ...
## $ menhinick : num 2.96 3.12 3.28 2.97 3.24 ...
## $ gleason : num 6.15 6.37 6.67 6.45 6.81 ...
## $ mcintosh : num 0.855 0.857 0.865 0.858 0.861 ...
## - attr(*, ".internal.selfref")=<externalptr># Quebra string em informações individuais.
boot <- boot %>%
separate(col = "cond",
into = c("u", "larg", "comp", "inten"),
sep = "[x_]",
remove = FALSE,
convert = TRUE) %>%
select(-u) %>%
unite(col = "dimen", larg, comp, sep = "x", remove = FALSE) %>%
mutate(cond = str_replace(cond, "cond_", ""),
area = larg * comp/10000,
vol_sum_area = media_sum_vol/area,
li_vol_area = li_vol/area,
lu_vol_area = lu_vol/area)
boot %>%
group_by(cond) %>%
summarise(area = first(area),
vol = mean(media_sum_vol),
vol_area = mean(media_sum_vol/area)) %>%
arrange(area) %>%
print(n = Inf)## `summarise()` ungrouping output (override with `.groups` argument)## # A tibble: 36 x 4
## cond area vol vol_area
## <chr> <dbl> <dbl> <dbl>
## 1 25x100_1 0.25 11.3 45.4
## 2 25x100_11 0.25 11.3 45.4
## 3 25x100_2 0.25 11.3 45.4
## 4 25x100_3 0.25 11.3 45.4
## 5 25x100_5 0.25 11.3 45.3
## 6 25x100_8 0.25 11.3 45.4
## 7 50x50_1 0.25 11.3 45.4
## 8 50x50_11 0.25 11.3 45.4
## 9 50x50_2 0.25 11.3 45.4
## 10 50x50_3 0.25 11.3 45.4
## 11 50x50_5 0.25 11.3 45.4
## 12 50x50_8 0.25 11.3 45.4
## 13 70.71x70.71_1 0.500 22.7 45.4
## 14 70.71x70.71_11 0.500 22.7 45.4
## 15 70.71x70.71_2 0.500 22.7 45.4
## 16 70.71x70.71_3 0.500 22.7 45.4
## 17 70.71x70.71_5 0.500 22.7 45.3
## 18 70.71x70.71_8 0.500 22.7 45.4
## 19 50x100_1 0.5 22.7 45.3
## 20 50x100_11 0.5 22.7 45.3
## 21 50x100_2 0.5 22.7 45.4
## 22 50x100_3 0.5 22.7 45.4
## 23 50x100_5 0.5 22.7 45.4
## 24 50x100_8 0.5 22.7 45.4
## 25 100x100_1 1 45.4 45.4
## 26 100x100_11 1 45.4 45.4
## 27 100x100_2 1 45.4 45.4
## 28 100x100_3 1 45.4 45.4
## 29 100x100_5 1 45.4 45.4
## 30 100x100_8 1 45.4 45.4
## 31 50x200_1 1 45.4 45.4
## 32 50x200_11 1 45.4 45.4
## 33 50x200_2 1 45.4 45.4
## 34 50x200_3 1 45.3 45.3
## 35 50x200_5 1 45.4 45.4
## 36 50x200_8 1 45.4 45.4boot <- boot %>%
mutate(cond = factor(cond),
cond = fct_reorder(cond, n_parc),
dimen = factor(dimen),
dimen = fct_reorder(dimen, larg * comp))
# Área total do inventário.
# censo_val$area <- with(boot[1, ], n_parc * area/(inten/100))
censo_val <- boot %>%
select(n_parc, larg, comp, inten) %>%
distinct() %>%
summarise(area = mean(n_parc * larg * comp/(inten/100)/10000)) %>%
bind_cols(censo_val) %>%
mutate(vol_med_area = vol_sum/area)
str(censo_val)## 'data.frame': 1 obs. of 12 variables:
## $ area : num 8476
## $ Alt_medio : num 14.3
## $ vol_medio : num 4.69
## $ vol_sum : num 387494
## $ CAP_medio : num 214
## $ n : int 82705
## $ n_esp : int 125
## $ margalef : num 11
## $ menhinick : num 0.435
## $ shannon : num 5.75
## $ mcintosh : num 0.835
## $ vol_med_area: num 45.7# Discrimina intervalos que contém o valor populacional.
boot <- boot %>%
mutate(cobre = {
li_vol_area < censo_val$vol_med_area &
lu_vol_area > censo_val$vol_med_area
})
str(boot)## Classes 'data.table' and 'data.frame': 108000 obs. of 26 variables:
## $ B : int 1 2 3 4 5 6 7 8 9 10 ...
## $ n_parc : int 85 85 85 85 85 85 85 85 85 85 ...
## $ prop_vazia : num 0 0.0118 0.0118 0 0 ...
## $ n_esp_parc : num 7.54 7.55 7.4 8.56 7.87 ...
## $ mediaw_media_cap: num 213 216 213 213 213 ...
## $ mediaw_media_alt: num 14.3 14.8 14.1 14.6 14.1 ...
## $ media_sum_vol : num 43.2 43.2 41.3 48.1 45.1 ...
## $ desvio_sum_vol : num 25.2 26.9 24.1 26.2 24.5 ...
## $ li_vol : num 37.8 37.4 36.1 42.5 39.8 ...
## $ lu_vol : num 48.7 49 46.5 53.8 50.4 ...
## $ cond : Factor w/ 36 levels "100x100_1","50x200_1",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ dimen : Factor w/ 6 levels "25x100","50x50",..: 5 5 5 5 5 5 5 5 5 5 ...
## $ larg : num 100 100 100 100 100 100 100 100 100 100 ...
## $ comp : num 100 100 100 100 100 100 100 100 100 100 ...
## $ inten : int 1 1 1 1 1 1 1 1 1 1 ...
## $ n_esp : int 84 87 91 88 93 86 87 89 89 80 ...
## $ shannon : num 5.57 5.62 5.7 5.63 5.71 ...
## $ margalef : num 12.4 12.9 13.5 12.8 13.7 ...
## $ menhinick : num 2.96 3.12 3.28 2.97 3.24 ...
## $ gleason : num 6.15 6.37 6.67 6.45 6.81 ...
## $ mcintosh : num 0.855 0.857 0.865 0.858 0.861 ...
## $ area : num 1 1 1 1 1 1 1 1 1 1 ...
## $ vol_sum_area : num 43.2 43.2 41.3 48.1 45.1 ...
## $ li_vol_area : num 37.8 37.4 36.1 42.5 39.8 ...
## $ lu_vol_area : num 48.7 49 46.5 53.8 50.4 ...
## $ cobre : logi TRUE TRUE TRUE TRUE TRUE TRUE ...Gráficos
Densidade
# lbs_den <- labs(y = "Densidade",
# subtitle = "Intensidade amostral (%) x dimensão de parcela")
# lbs_box <- labs(x = "Dimensão da parcela",
# subtitle = "Intensidade amostral (%)")
# thm_box <- theme(axis.text.x = element_text(angle = 90,
# hjust = 1,
# vjust = 0.5))
lbs_den <- labs(y = "Density",
subtitle = "Sampling intensity (%) x plot dimensions")
lbs_box <- labs(x = "Plot dimensions",
subtitle = "Sampling intensity (%)")
thm_box <- theme(axis.text.x = element_text(angle = 90,
hjust = 1,
vjust = 0.5))
ggplot(boot) +
aes(x = prop_vazia) +
facet_grid(facets = dimen ~ inten) +
geom_density() +
# labs(x = "Proporção de pacelas vazias",
# title = "Distribuição de parcelas vazias") +
labs(x = "Proportion of empty plots (3000 bootstrap samples)",
title = "Proportion of empty plots") +
lbs_den
ggplot(boot) +
aes(x = dimen, y = prop_vazia) +
facet_grid(facets = ~inten) +
geom_boxplot() +
# labs(y = "Proporção de pacelas vazias",
# title = "Distribuição de parcelas vazias") +
labs(y = "Proportion of empty plots (3000 bootstrap samples)",
title = "Proportion of empty plots") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = mediaw_media_cap) +
facet_grid(facets = dimen ~ inten) +
geom_density() +
geom_vline(xintercept = censo_val$CAP_medio,
size = 0.25,
linetype = 2) +
# labs(x = "Circunferência à altura do peito média",
# title = "Distribuição da circunferência à altura do peito") +
labs(x = "Mean circunference at breast height (3000 bootstrap samples)",
title = "Distribution of circunference at breast height") +
lbs_den
ggplot(boot) +
aes(x = dimen, y = mediaw_media_cap) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$CAP_medio,
size = 0.25,
linetype = 2) +
# labs(y = "Circunferência à altura do peito média",
# title = "Distribuição da circunferência à altura do peito") +
labs(y = "Mean circunference at breast height (3000 bootstrap samples)",
title = "Distribution of circunference at breast height") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = mediaw_media_alt) +
facet_grid(facets = dimen ~ inten) +
geom_density() +
geom_vline(xintercept = censo_val$Alt_medio,
size = 0.25,
linetype = 2) +
# labs(x = "Altura total média",
# title = "Distribuição da altura total") +
labs(x = "Mean total height (3000 bootstrap samples)",
title = "Distribution of total height") +
lbs_den
ggplot(boot) +
aes(x = dimen, y = mediaw_media_alt) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$Alt_medio,
size = 0.25,
linetype = 2) +
# labs(y = "Altura total média",
# title = "Distribuição da altura total") +
labs(y = "Mean total height (3000 bootstrap samples)",
title = "Distribution of total height") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = vol_sum_area) +
facet_grid(facets = dimen ~ inten) +
geom_density() +
geom_vline(xintercept = censo_val$vol_med_area,
size = 0.25,
linetype = 2) +
# labs(x = "Volume total médio por hectare",
# title = "Distribuição do volume total") +
labs(x = "Mean total volume per hectare (3000 bootstrap samples)",
title = "Distribution of total volume per hectare") +
lbs_den
ggplot(boot) +
aes(x = dimen, y = vol_sum_area) +
# aes(x = dimen, y = media_sum_vol) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$vol_med_area,
size = 0.25,
linetype = 2) +
# labs(y = "Volume total médio por hectare",
# title = "Distribuição do volume total") +
labs(y = "Mean total volume per hectare (3000 bootstrap samples)",
title = "Distribution of total volume per hectare") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = n_esp) +
facet_grid(facets = dimen ~ inten) +
geom_histogram(binwidth = 1) +
geom_vline(xintercept = censo_val$n_esp,
linetype = 2,
color = "tomato") +
lbs_den +
labs(y = "Frequência",
x = "Total de espécies",
title = "Número de espécies")
ggplot(boot) +
aes(x = dimen, y = n_esp) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$n_esp,
size = 0.25,
linetype = 2) +
labs(y = "Total de espécies",
title = "Número de espécies") +
lbs_box +
thm_box
Índices ecológicos
ggplot(boot) +
aes(x = dimen, y = shannon) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$shannon,
size = 0.25,
linetype = 2) +
labs(y = "Shannon",
title = "Shannon") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = dimen, y = margalef) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$margalef,
size = 0.25,
linetype = 2) +
labs(y = "Margalef",
title = "Margalef") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = dimen, y = menhinick) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$menhinick,
size = 0.25,
linetype = 2) +
labs(y = "Menhinick",
title = "Menhinick") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = dimen, y = gleason) +
facet_grid(facets = ~inten) +
geom_boxplot() +
labs(y = "Gleason",
title = "Gleason") +
lbs_box +
thm_box
ggplot(boot) +
aes(x = dimen, y = mcintosh) +
facet_grid(facets = ~inten) +
geom_boxplot() +
geom_hline(yintercept = censo_val$mcintosh,
size = 0.25,
linetype = 2) +
labs(y = "Mcintosh",
title = "Mcintosh") +
lbs_box +
thm_box
Intervalos de confiança
ggplot(boot) +
facet_grid(facets = dimen ~ inten) +
aes(sample = vol_sum_area) +
geom_qq(pch = 1) +
# geom_qq_line(color = "orange") +
theme()
# str(boot)
boot %>%
filter(cond == "100x100_1") %>%
ggplot() +
aes(x = reorder(B, vol_sum_area),
y = vol_sum_area,
color = cobre) +
geom_point() +
geom_errorbar(aes(ymin = li_vol_area, ymax = lu_vol_area)) +
geom_line(aes(group = 1), color = "black") +
geom_hline(yintercept = censo_val$vol_med_area,
linetype = 2,
color = "orange") +
scale_x_discrete(breaks = NULL) +
scale_color_manual(values = c("orange", "gray30")) +
guides(color = FALSE) +
expand_limits(x = extendrange(c(1, 3000), f = 0.01)) +
labs(x = "Estimativa bootstrap",
y = "Volume de madeira da parcela",
subtitle = "Parcela de 100 x 100 m · 1% da área total")
Cobertura dos intervalos de confiança
# Cobertura dos intervalos de confiança.
cober <- boot %>%
group_by(cond) %>%
summarise(cober = mean(cobre),
larg = first(larg),
comp = first(comp),
inten = first(inten),
dimen = first(dimen))## `summarise()` ungrouping output (override with `.groups` argument)## # A tibble: 36 x 6
## cond cober larg comp inten dimen
## <fct> <dbl> <dbl> <dbl> <int> <fct>
## 1 100x100_1 0.938 100 100 1 100x100
## 2 50x200_1 0.941 50 200 1 50x200
## 3 100x100_2 0.940 100 100 2 100x100
## 4 50x100_1 0.940 50 100 1 50x100
## 5 50x200_2 0.946 50 200 2 50x200
## 6 70.71x70.71_1 0.944 70.7 70.7 1 70.71x70.71
## 7 100x100_3 0.95 100 100 3 100x100
## 8 50x200_3 0.940 50 200 3 50x200
## 9 25x100_1 0.946 25 100 1 25x100
## 10 50x100_2 0.943 50 100 2 50x100
## # … with 26 more rowsggplot(cober) +
aes(x = inten, y = cober, color = dimen) +
geom_point() +
geom_line(aes(group = dimen)) +
geom_hline(yintercept = 0.95, color = "black", linetype = 2) +
labs(x = "Intensidade amostral (%)",
y = "Taxa de cobertura do intervalo de confiança",
color = "Dimensão (m)") +
theme(legend.position = c(0.025, 0.025),
legend.justification = c(0, 0))