6. Amostragem

library(tidyverse)
── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.2     ✔ readr     2.1.4
✔ forcats   1.0.0     ✔ stringr   1.5.0
✔ ggplot2   3.4.3     ✔ tibble    3.2.1
✔ lubridate 1.9.2     ✔ tidyr     1.3.0
✔ purrr     1.0.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(rio)
df <- import("https://docs.google.com/spreadsheets/d/18aXD_2ISvzB8h8_kgOfSBbr9a9d9pT0QVazt-KjVLRw/edit#gid=1590128876")
df
   id                                 aluno crm altura
1   1                Aline Kliauga Ferrante  xx     NA
2   2                     André Cuneo Sagaz  xy     NA
3   3              Beatriz Haensel Teixeira  xx     NA
4   4                        Bruna Waltrich  xx    156
5   5       Carlos Eduardo Forcelini Assoni  xy     NA
6   7                Danielli Zangalli Kern  xx    169
7  11                 Gabriela Araujo Catto  xx    162
8   8                 Felipe Ricci Westphal  xy     NA
9  13         Helena dos Santos Vanderlinde  xx    180
10 10                   Gabriel Sbardelotto  xy     NA
11 14               Isabela Martins Ghizoni  xx    166
12 20        Juliana Eduarda Oliveira Gomes  xx    165
13 21        Kamilly Vitoria Siqueira Tonet  xx    164
14 29          Maria Eduarda Wendt Coutinho  xx    170
15 15        João Pedro dos Santos Angulski  xy     NA
16 30         Maria Laura Faustino Monteiro  xx    164
17 35          Rafaela Rodrigues dos Santos  xx    158
18 37                        Tassie Turcato  xx    165
19 38         Thalita Maria Gomes Rodrigues  xx    163
20 39        Wanessa Pedrinha do Nascimento  xx    177
21  6              Daniel Schmechel Affeldt  xy    169
22 22                  Leticia Herbert Post  xx     NA
23  9                 Gabriel Loche Lopasso  xy    188
24 12            Guilherme Antonio Ferreira  xy    173
25 16          Joao Pedro Lantyer Marcelino  xy    190
26 26                Luis Bortoluzzi Sobral  xy     NA
27 27                   Luiz Paulo da Silva  xy     NA
28 28              Magalhães Antonio Saxico  xy     NA
29 17                    João Vitor Germano  xy    166
30 18           João Vitor Zeferino Madeira  xy    168
31 31               Mateus Zunino Espindola  xy     NA
32 32             Matheus de Oliveira Mussi  xy     NA
33 33 Natacha Micheline de Oliveira da Rosa  xx     NA
34 34           Pierre Marcel Bruno Boisson  xy     NA
35 19         José Eduardo Pimentel e Silva  xy    169
36 23                   Lucas Lopes Ribeiro  xy    175
37 24                    Lucas Paz Claudino  xy    178
38 25               Lucas Sodre de Oliveira  xy    174
39 36              Renan Guilherme da Silva  xy    179
40 40                 Wesley Castilhos Drun  xy     NA
41 41                      Wilson Rosa Neto  xy    173

A função sample_random() retirada do pacote metan pode ser utilizada para amostrar n linhas aleatoriamente do conjunto de dados data. Utilizando prop, uma proporção dos dados é amostrada. Este último é útil ao realizar amostragens estratificadas informando o argumento by, onde cada estrato possui diferentes tamanhos de amostra.

# Função auxiliar
sample_random <- function(data,
                          n,
                          prop,
                          by = NULL,
                          weight = NULL){
  if(!missing(by)){
    data <- data |> group_by({{by}})
  }
  dplyr::slice_sample(data, n = n, prop = prop, weight_by = weight) |>
    ungroup()
}

Amostragem aleatória simples

A amostragem aleatória simples é o método mais básico de amostragem que tanto pode ser utilizado diretamente na seleção de uma amostra, quando se conhece os indivíduos da população, como ser parte de outros planos amostrais como, por exemplo, da Amostragem Estratificada.

O número possíveis de amostras com n indivíduos de uma população com N elementos, é dado por:

\[ \mathop C\nolimits_N^n = \frac{N!}{n!(N-n)!} \]

No próximo exemplo, veremos as possíveis amostras (120) com 3 indivíduos, tomadas de uma população com 10 indivíduos.

N <- 10
n <- 3
d <- combn(N, n)
t(d)
       [,1] [,2] [,3]
  [1,]    1    2    3
  [2,]    1    2    4
  [3,]    1    2    5
  [4,]    1    2    6
  [5,]    1    2    7
  [6,]    1    2    8
  [7,]    1    2    9
  [8,]    1    2   10
  [9,]    1    3    4
 [10,]    1    3    5
 [11,]    1    3    6
 [12,]    1    3    7
 [13,]    1    3    8
 [14,]    1    3    9
 [15,]    1    3   10
 [16,]    1    4    5
 [17,]    1    4    6
 [18,]    1    4    7
 [19,]    1    4    8
 [20,]    1    4    9
 [21,]    1    4   10
 [22,]    1    5    6
 [23,]    1    5    7
 [24,]    1    5    8
 [25,]    1    5    9
 [26,]    1    5   10
 [27,]    1    6    7
 [28,]    1    6    8
 [29,]    1    6    9
 [30,]    1    6   10
 [31,]    1    7    8
 [32,]    1    7    9
 [33,]    1    7   10
 [34,]    1    8    9
 [35,]    1    8   10
 [36,]    1    9   10
 [37,]    2    3    4
 [38,]    2    3    5
 [39,]    2    3    6
 [40,]    2    3    7
 [41,]    2    3    8
 [42,]    2    3    9
 [43,]    2    3   10
 [44,]    2    4    5
 [45,]    2    4    6
 [46,]    2    4    7
 [47,]    2    4    8
 [48,]    2    4    9
 [49,]    2    4   10
 [50,]    2    5    6
 [51,]    2    5    7
 [52,]    2    5    8
 [53,]    2    5    9
 [54,]    2    5   10
 [55,]    2    6    7
 [56,]    2    6    8
 [57,]    2    6    9
 [58,]    2    6   10
 [59,]    2    7    8
 [60,]    2    7    9
 [61,]    2    7   10
 [62,]    2    8    9
 [63,]    2    8   10
 [64,]    2    9   10
 [65,]    3    4    5
 [66,]    3    4    6
 [67,]    3    4    7
 [68,]    3    4    8
 [69,]    3    4    9
 [70,]    3    4   10
 [71,]    3    5    6
 [72,]    3    5    7
 [73,]    3    5    8
 [74,]    3    5    9
 [75,]    3    5   10
 [76,]    3    6    7
 [77,]    3    6    8
 [78,]    3    6    9
 [79,]    3    6   10
 [80,]    3    7    8
 [81,]    3    7    9
 [82,]    3    7   10
 [83,]    3    8    9
 [84,]    3    8   10
 [85,]    3    9   10
 [86,]    4    5    6
 [87,]    4    5    7
 [88,]    4    5    8
 [89,]    4    5    9
 [90,]    4    5   10
 [91,]    4    6    7
 [92,]    4    6    8
 [93,]    4    6    9
 [94,]    4    6   10
 [95,]    4    7    8
 [96,]    4    7    9
 [97,]    4    7   10
 [98,]    4    8    9
 [99,]    4    8   10
[100,]    4    9   10
[101,]    5    6    7
[102,]    5    6    8
[103,]    5    6    9
[104,]    5    6   10
[105,]    5    7    8
[106,]    5    7    9
[107,]    5    7   10
[108,]    5    8    9
[109,]    5    8   10
[110,]    5    9   10
[111,]    6    7    8
[112,]    6    7    9
[113,]    6    7   10
[114,]    6    8    9
[115,]    6    8   10
[116,]    6    9   10
[117,]    7    8    9
[118,]    7    8   10
[119,]    7    9   10
[120,]    8    9   10

Aplicação

Vamos considerar uma variável x, distribuida normalmente com média \(\bar X = 10\) e desvio padrão \(S = 2\), avaliada em população com N = 10.

set.seed(1)
N <- 10
df2 <- data.frame(id = 1:N,
                  x = rnorm(n = N, mean = 10, sd = 2))
df2
   id         x
1   1  8.747092
2   2 10.367287
3   3  8.328743
4   4 13.190562
5   5 10.659016
6   6  8.359063
7   7 10.974858
8   8 11.476649
9   9 11.151563
10 10  9.389223

Considerando uma amostragem com n = 3, as 120 amostras possíveis são

n <- 3
amostras <- combn(N, n) |> t()
amostras |> head()
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    1    2    4
[3,]    1    2    5
[4,]    1    2    6
[5,]    1    2    7
[6,]    1    2    8
amostras |> tail()
       [,1] [,2] [,3]
[115,]    6    8   10
[116,]    6    9   10
[117,]    7    8    9
[118,]    7    8   10
[119,]    7    9   10
[120,]    8    9   10

Médias amostrais

A seguinte função, computa a média das 120 amostras. Assim, obtém-se a distribuição das médias amostrais.

medias <- NULL
# abordagem com for-loop
for (i in 1:nrow(amostras)) {
  individ <- amostras[i,]
  valores <- df2$x[individ]
  medias <- append(medias, mean(valores))
}

# criar um data frame com as médias
df_medias <- data.frame(amostras) |> mutate(media = medias)
head(df_medias)
  X1 X2 X3     media
1  1  2  3  9.147707
2  1  2  4 10.768314
3  1  2  5  9.924465
4  1  2  6  9.157814
5  1  2  7 10.029746
6  1  2  8 10.197009
tail(df_medias)
    X1 X2 X3     media
115  6  8 10  9.741645
116  6  9 10  9.633283
117  7  8  9 11.201023
118  7  8 10 10.613577
119  7  9 10 10.505215
120  8  9 10 10.672478

Ao computar a média das medias amostrais, obtém-se a média populacional

med_amostral <- mean(df_medias$media)
med_pop <-  mean(df2$x)

identical(med_amostral, med_pop)
[1] TRUE
ggplot(df_medias, aes(x = media)) +
  geom_histogram(bins = 8, color = "black", fill = "gray") +
  geom_vline(xintercept = med_pop, color = "red", size = 1)
Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
ℹ Please use `linewidth` instead.

Tamanho da amostra vs acurácia

No seguinte exemplo, vamos investigar o impacto do tamanho da amostra na acurácia da média. Para isso, serão amostradas aleatoriamente 1:120 médias amostrais do conjunto df_medias e calculado o desvio em relação a média populacional. O processo é repetido nboot vezes utilizando a técnica bootstrap.

nboot <- 200

samples <- list()
for(j in 1:nboot){
  tmp <- 
    map_dbl(1:nrow(df_medias), function(x){
      rows <- sample(1:nrow(df_medias), x)
      mean(df_medias[rows,]$media)
    })
  samples[[j]] <- tmp
}

# cada coluna contém as médias amostrais de um procedimento bootstrap
samples <- do.call(cbind, lapply(samples, data.frame))
colnames(samples) <- paste0("v", 1:ncol(samples))

# criar dados longo
samples_long <- 
  samples |> 
  pivot_longer(everything()) |> 
  mutate(desvio  = value -  mean(df2$x),
         x = rep(1:nrow(df_medias), each = nboot))

# média dos procedimentos bootstrap
samples_mean <- 
  samples_long |> 
  group_by(x) |> 
  summarise(mu = mean(desvio))


# criar o gráfico
ggplot(samples_long, aes(x = x,
                         y = desvio,
                         group = name)) +
  geom_line(alpha = 0.1) +
  geom_line(aes(x = x, y = mu, group = 1),
            data = samples_mean,
            color = "red") +
  scale_x_continuous(breaks = seq(0, 120, by = 10)) +
  theme_bw(base_size = 16) +
  theme(panel.grid.minor = element_blank()) +
  labs(x = "Número de médias amostrais incluídas",
       y = "Desvio em relação a média populacional",
       title = "Resultado de 200 procedimentos bootstrap",
       caption = glue::glue("Média: {round(mean(df2$x), 3)}"))

Exemplo da altura da turma

df_turma <- 
  import("https://docs.google.com/spreadsheets/d/18aXD_2ISvzB8h8_kgOfSBbr9a9d9pT0QVazt-KjVLRw/edit#gid=1590128876") |> 
  metan::remove_rows_na()
Registered S3 method overwritten by 'GGally':
  method from   
  +.gg   ggplot2
Warning: Row(s) 1, 2, 3, 5, 8, 10, 15, 22, 26, 27, 28, 31, 32, 33, 34, 40 with
NA values deleted.
linhas <- sample(1:25, 4)
df_turma[linhas,]
   id                         aluno crm altura
39 36      Renan Guilherme da Silva  xy    179
24 12    Guilherme Antonio Ferreira  xy    173
19 38 Thalita Maria Gomes Rodrigues  xx    163
9  13 Helena dos Santos Vanderlinde  xx    180
set.seed(4)

nboot <- 1000

samples <- list()
for(j in 1:nboot){
  tmp <- 
    map_dbl(1:nrow(df_turma), function(x){
      rows <- sample(1:nrow(df_turma), x)
      mean(df_turma[rows,]$altura)
    })
  samples[[j]] <- tmp
}

# cada coluna contém as médias amostrais de um procedimento bootstrap
samples <- do.call(cbind, lapply(samples, data.frame))
colnames(samples) <- paste0("v", 1:ncol(samples))

# criar dados longo
samples_long <- 
  samples |> 
  pivot_longer(everything()) |> 
  mutate(desvio  = value -  mean(df_turma$altura),
         x = rep(1:nrow(df_turma), each = nboot))

# média dos procedimentos bootstrap
samples_mean <- 
  samples_long |> 
  group_by(x) |> 
  summarise(mu = mean(desvio))


# criar o gráfico
ggplot(samples_long, aes(x = x,
                         y = desvio,
                         group = name)) +
  geom_line(alpha = 0.1) +
  geom_line(aes(x = x, y = mu, group = 1),
            data = samples_mean,
            color = "red") +
  scale_x_continuous(breaks = seq(1, nrow(df_turma), by = 1)) +
  theme_bw(base_size = 16) +
  theme(panel.grid.minor = element_blank()) +
  labs(x = "Número de alunos incluídos na amostra",
       y = "Desvio em relação a média populacional",
       title = "Resultado de 200 procedimentos bootstrap",
       caption = glue::glue("Média: {round(mean(df_turma$altura), 3)}; N = {nrow(df_turma)}"))

Abordagem paralela

Quando o número de amostras cresce bastante, a abordagem for-loop não é computacionalmente eficiente. Assim, uma abordagem utilizando sapply() é mais eficiente. Quando paraleliza-se a função, a eficiência aumenta mais ainda.

# criando uma função para obter a média de um id
get_mean <- function(df, var, amostras, id){
  individ <- amostras[id,]
  mean(df[[var]][individ])
}


N <- 30
n <- 5
df2 <- data.frame(id = 1:N,
                  x = rnorm(n = N, mean = 10, sd = 2))
amostras2 <- combn(N, n) |> t()


system.time(
  medias2 <-
    map_dbl(1:nrow(amostras2), function(i){
      get_mean(df2, "x", amostras2, id = i)
    })
)


library(parallel)
clust <- makeCluster(5)
clusterExport(clust,
              varlist = c("df2", "amostras2", "get_mean"))
system.time(
  medias3 <-
    parLapply(clust, 1:nrow(amostras2), function(i){
      get_mean(df2, "x", amostras2, id = i)
    })
)
stopCluster(clust)

Amostragem Aleatória Estratificada

Número igual dentro de cada estrato

sample_random(df, n = 3, by = crm)
# A tibble: 6 × 4
     id aluno                          crm   altura
  <int> <chr>                          <chr>  <int>
1    39 Wanessa Pedrinha do Nascimento xx       177
2    35 Rafaela Rodrigues dos Santos   xx       158
3    29 Maria Eduarda Wendt Coutinho   xx       170
4    34 Pierre Marcel Bruno Boisson    xy        NA
5    40 Wesley Castilhos Drun          xy        NA
6    12 Guilherme Antonio Ferreira     xy       173

Proporção da população em cada estrato

sample_random(df,
              prop = 0.3,
              by = crm)
# A tibble: 12 × 4
      id aluno                           crm   altura
   <int> <chr>                           <chr>  <int>
 1     3 Beatriz Haensel Teixeira        xx        NA
 2    39 Wanessa Pedrinha do Nascimento  xx       177
 3     1 Aline Kliauga Ferrante          xx        NA
 4    22 Leticia Herbert Post            xx        NA
 5     4 Bruna Waltrich                  xx       156
 6     5 Carlos Eduardo Forcelini Assoni xy        NA
 7    31 Mateus Zunino Espindola         xy        NA
 8    18 João Vitor Zeferino Madeira     xy       168
 9    17 João Vitor Germano              xy       166
10    24 Lucas Paz Claudino              xy       178
11    26 Luis Bortoluzzi Sobral          xy        NA
12    40 Wesley Castilhos Drun           xy        NA

Amostragem aleatória sistemática

sample_systematic <- function(data,
                              n,
                              r = NULL,
                              by =  NULL){
  aux <- function(data, n, r = NULL){
    k <- floor(nrow(data) / n)
    message("k = ", k)
    if(is.null(r)){
      r <- sample(1:k, 1)
    }
    if(r == 1){
      rows <- sample(1:nrow(data), n)
    } else{
      rows <- seq(r, r + k*(n-1), k)
    }
    slice(data, rows) |>
      mutate(.id = rows, .before = 1)
  }
  if(!missing(by)){
    data <- data |> group_by({{by}})
  }
  if(is.grouped_df(data)){
    groups <- group_vars(data)
    data |>
      ungroup() |>
      nest(data = -c(!!!syms(groups))) |>
      mutate(sample = map(data, ~.x |>  aux(n = n, r = r))) |>
      select(-data) |>
      unnest(sample)
  } else{
    aux(data, n = n, r = r)
  }
}
set.seed(1)
sample_systematic(df, n = 4)
k = 10
  .id id                         aluno crm altura
1   9 13 Helena dos Santos Vanderlinde  xx    180
2  19 38 Thalita Maria Gomes Rodrigues  xx    163
3  29 17            João Vitor Germano  xy    166
4  39 36      Renan Guilherme da Silva  xy    179
Free Website Hit Counter
Free website hit counter