Introdução

Esse artigo é uma explicação do código descrito nesse link referente a uma palestra de valor esperado do curso de Inferência Estatística da especialização em Ciências de Dados do Coursera.

Primeiro, definiremos as bibliotecas que vamos usar, o número de simulações nosim e os tamanhos das amostras para calcular as médias em amostra, e a semente da simulação para que seja reprodutível.

library(ggplot2)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
nosim <- 10000
amostra <- c(10,20,30)
set.seed(123)

Os dataframes serão criados com duas colunas, uma com a média dos lançamentos e outra com o tamanho da amostra. A coluna com a média será determinada pelos os resultados dos lançamentos simples e dos lançamentos baseadas nos tamanhos dados em amostra. Por exemplo, se o primeiro elemento de amostra é 10, cada elemento relativa a 10 terá a média de 10 lançamentos simulados. Para fazer esses cálculos mais rapidamente, criaremos uma matriz com um número de linhas igual a nosim e o número de colunas igual ao valor correspondente na amostra, com cada entrada referente a um lançamento com reposição, determinado pela função sample. Após isso, utilizaremos a função apply em cada linha de cada matriz para calcular a média. Faremos um exemplo abaixo.

Vamos supor que 0 represente “cara” e 1 “coroa” num lançamento de moedas. Vamos criar uma matriz com 10 linhas, onde cada linha representa 10 lançamentos aleatórios de moeda, e armazenar na variável exemplo.

exemplo <- matrix(sample(0 : 1, 10 * amostra[1], replace = TRUE),
                           10)
exemplo
##       [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
##  [1,]    0    1    0    0    0    1    0    1    1     1
##  [2,]    0    1    1    1    1    0    1    0    1     1
##  [3,]    0    1    0    0    1    0    1    0    0     0
##  [4,]    1    0    0    1    0    0    0    0    1     1
##  [5,]    0    1    0    1    0    0    0    0    1     1
##  [6,]    1    0    0    0    0    1    1    1    1     0
##  [7,]    1    1    1    0    0    1    0    1    0     0
##  [8,]    1    0    1    0    1    0    0    0    0     1
##  [9,]    0    0    0    0    0    1    0    1    1     0
## [10,]    0    0    1    1    0    0    0    1    0     0

Agora, vamos aplicar a função apply nessa matriz, aplicando a função mean na linhas (número 1 na fórmula).

apply(exemplo,1,mean)
##  [1] 0.5 0.7 0.3 0.4 0.4 0.5 0.5 0.4 0.3 0.3

Utilizaremos esse procedimento para cada valor em amostra, com cada matriz tendo o tamanho de nosim. Esses valores serão juntados na variável x do dataframe. Para terminar o dataframe, criaremos uma coluna do tipo factor para dizer qual o tamanho da amostra de cada lançamento. Assim, para nosso exemplo, o dataframe terá 400000 linhas e duas colunas.

Lançamentos de moedas

Com isso, o dataframe para as simulações de lançamentos de moedas vai ser dado mostrado abaixo. Adicionaremos uma amostra de linhas desse dataframe e vamos criar o dataframe mu que será a média total das amostras.

dat <- data.frame(
        #Lançamentos simples
        x = c(sample(0 : 1, nosim, replace = TRUE),
        #Médias de lançamentos para primeiro elemento de amostra
        apply(matrix(sample(0 : 1, nosim * amostra[1], replace = TRUE),
        nosim), 1, mean),
        #Médias de lançamentos para segundo elemento de amostra
        apply(matrix(sample(0 : 1, nosim * amostra[2], replace = TRUE),
        nosim), 1, mean),
        #Médias de lançamentos para terceiro elemento de amostra
        apply(matrix(sample(0 : 1, nosim * amostra[3], replace = TRUE),
        nosim), 1, mean)
        ),
        #fatores
        size = factor(rep(c(1, amostra[1], amostra[2], amostra[3]),
                          rep(nosim, 4)))
        )
dat[sample(1:40000,10),]
##               x size
## 16229 0.6000000   10
## 16742 0.5000000   10
## 33643 0.4666667   30
## 1488  1.0000000    1
## 5228  1.0000000    1
## 14061 0.6000000   10
## 8860  1.0000000    1
## 14447 0.6000000   10
## 35467 0.4333333   30
## 21383 0.3500000   20
#Dataframe com a média por grupo
mu <- dat %>% group_by(size) %>% summarise(grp.mean=mean(x),.groups="drop")
mu
## # A tibble: 4 x 2
##   size  grp.mean
##   <fct>    <dbl>
## 1 1        0.499
## 2 10       0.498
## 3 20       0.500
## 4 30       0.499

O valor da média grp.mean corresponde ao valor esperado calculado:

\[E[X] = 0.5 \times 0 + 0.5 \times 1 = 0.5 \]

Agora faremos um gráfico com a biblioteca ggplot2, contendo um histograma das médias para cada fator, e uma linha tracejada em cada histograma com a média total de cada fator.

g <- ggplot(dat, aes(x = x, fill = size)) + 
        geom_histogram( binwidth = 1/12, colour = "black")+
        geom_density(alpha=0.2)+ facet_grid(. ~ size)+
        geom_vline(data=mu, aes(xintercept=grp.mean),
                   linetype="dashed")+
        labs( x="Média das amostras",y = "Número de lançamentos",fill="Tamanho \n amostra")

g 
Distribuição de médias de amostras de 
 lançamentos de moedas

Distribuição de médias de amostras de lançamentos de moedas

Lançamentos de dados

Faremos o mesmo do caso anterior, só trocando a amostra, que antes era de entre 0 e 1 no caso das moedas e agora faremos entre 1 e 6 representando cada face do dado.

dat <- data.frame(
        #Lançamentos simples
        x = c(sample(1 : 6, nosim, replace = TRUE),
        #Médias de lançamentos para primeiro elemento de amostra
        apply(matrix(sample(1 : 6, nosim * amostra[1], replace = TRUE),
        nosim), 1, mean),
        #Médias de lançamentos para segundo elemento de amostra
        apply(matrix(sample(1 : 6, nosim * amostra[2], replace = TRUE),
        nosim), 1, mean),
        #Médias de lançamentos para terceiro elemento de amostra
        apply(matrix(sample(1 : 6, nosim * amostra[3], replace = TRUE),
        nosim), 1, mean)
        ),
        #fatores
        size = factor(rep(c(1, amostra[1], amostra[2], amostra[3]),
                          rep(nosim, 4)))
        )
dat[sample(1:40000,10),]
##              x size
## 36487 3.433333   30
## 12403 3.400000   10
## 17617 3.700000   10
## 17386 4.600000   10
## 355   3.000000    1
## 17040 1.700000   10
## 13165 3.600000   10
## 10373 3.200000   10
## 498   5.000000    1
## 1438  2.000000    1
#Dataframe com a média por grupo
mu <- dat %>% group_by(size) %>% summarise(grp.mean=mean(x),.groups="drop")
mu
## # A tibble: 4 x 2
##   size  grp.mean
##   <fct>    <dbl>
## 1 1         3.48
## 2 10        3.49
## 3 20        3.50
## 4 30        3.50

O valor da média grp.mean corresponde ao valor esperado calculado:

\[ E[X] = 1 \times \frac{1}{6} + 2 \times \frac{1}{6} + 3 \times \frac{1}{6} + 4 \times \frac{1}{6} + 5 \times \frac{1}{6} + 6 \times \frac{1}{6} = 3.5 \]

E o gráfico é dado por

g <- ggplot(dat, aes(x = x, fill = size)) + 
        geom_histogram(alpha=0.4, binwidth = .25, colour = "black")+
        geom_density(alpha=0.2)+ facet_grid(. ~ size)+
        geom_vline(data=mu, aes(xintercept=grp.mean),
                   linetype="dashed")+
        labs(x="Média das amostras",y = "Número de lançamentos",fill="Tamanho \n amostra")

g 
Distribuição de médias de amostras de 
 lançamentos de dados

Distribuição de médias de amostras de lançamentos de dados