Neste documento faremos algumas simulações e análises dos lançamentos de Moeda.
Vamos começar carregando os pacotes necessários.
#Descomente se for necessário instalar os pacotes a seguir:
#install.packages("dplyr") install.packages("Rmpfr")
suppressMessages(library(dplyr))
library(Rmpfr)
## Loading required package: gmp
##
## Attaching package: 'gmp'
## The following objects are masked from 'package:base':
##
## %*%, apply, crossprod, matrix, tcrossprod
## C code of R package 'Rmpfr': GMP using 64 bits per limb
##
## Attaching package: 'Rmpfr'
## The following objects are masked from 'package:stats':
##
## dbinom, dnorm, dpois, pnorm
## The following objects are masked from 'package:base':
##
## cbind, pmax, pmin, rbind
Vamos começar criando uma função flips. A função flips(v) simula o lançamento de v moedas. As saídas de nossas moedas serão \(\{0,1\}\).
sample.space <- c(0,1)
theta <- 0.5 # moeda honesta
flips <- function(v){ sample(sample.space,
size = v,
replace = TRUE,
prob = c(theta, 1 - theta))}
# Jogar uma moeda 100 vezes
flips(100)
## [1] 0 0 0 0 0 1 0 1 1 0 1 0 1 0 0 0 1 0 0 0 0 0 1 0 1 0 1 0 1 0 1 0 0 1 0
## [36] 1 0 1 0 1 0 1 1 1 0 0 1 0 1 1 1 1 0 1 1 0 1 0 1 1 0 1 0 0 1 0 0 0 0 0
## [71] 0 0 1 1 1 0 0 1 0 1 0 0 1 1 0 0 0 1 0 0 1 1 1 1 1 0 0 1 0 0
O conjunto dos nímeros normais em \((0,1]\) é definido como
\[ \begin{aligned} N &:=\{\omega\in(0,1]: \lim_{n\rightarrow \infty} {s_n(\omega)}/n =0\}\\ &\phantom{:}=\{\omega\in(0,1]: \lim_{n\rightarrow \infty} \frac{1}{n}\displaystyle\sum_{k=1}^n d_k(\omega) = \frac{1}{2}\}. \end{aligned} \]
O conjunto dos números anormais é definido como \(A:=(0,1]-N\).
** Conjuntos Negligencíaveis**
Um subconjunto \(B\subset (0,1]\) é dito negligenciável se para todo \(\epsilon>0\), existem subconjuntos \(B_1, B_2,\ldots\) de \(\mathcal{B}_0((0,1])\) tal que \[ \begin{aligned} &B\subset \bigcup_{k=1}^\infty B_k\quad\text{e}\quad\sum_{k=1}^\infty P[B_k]\leq \epsilon. \end{aligned} \]
** Teorema dos Números Normais de Borel**
O conjunto dos números anormais, \(A\), é negligenciável.
Vamos simular o comportamento de cem somas cumulativas no mesmo gráfico
# Vamos criar uma lista com todas as cores disponíveis
color = grDevices::colors()[grep('gr(a|e)y', grDevices::colors(), invert = T)]
# plotando a soma cumulativa
plot(cummean(flips(20000)), type = 'l',col="red", main="Normalidade de 50 números",ylab ="Médias Parciais",ylim=c(0,1) )
#adicionando outras 49 somas
for (i in 1:50) {
lines(cummean(flips(20000)), type = 'l',col= color[7*i])
}
Vamos ver o que acontece com a média após 20000 dígitos de 1000 números
milnumeros <-replicate(1000,mean(flips(20000)))
hist(milnumeros,col="gray",xlim=c(0.48,0.52), main="Histograma da média após 20000 dígitos")
Vamos ver o que acontece com a média após 200000 dígitos de 1000 números
milnumeros <-replicate(1000,mean(flips(200000)))
hist(milnumeros, col="gray",xlim=c(0.48,0.52),, main="Histograma da média após 200000 dígitos")
Não é conhecido se \(\pi\) é normal, muito menos absolutamente normal.
Para tornar as coisas mais rápidas usamos o arquivo pi-10million.txt que contém 10 milhões de dígito de \(\pi\)
#
myPi = readLines("pi-10million.txt")[1]
## Warning in readLines("pi-10million.txt"): linha final incompleta encontrada
## em 'pi-10million.txt'
# I think this is how I managed to get Pi into a vector, it wasn't easy.
piVector = unlist(strsplit(myPi,""))
piVector = unlist(lapply(piVector,as.numeric))
Fazemos um histograma da frequência dos dígitos
hist(piVector[1:1000000],main="Distribuição dos dígitos de pi",col="blue", breaks=c(-0.1,0.9,1.9,2.9,3.9,4.9,5.9,6.9,7.9,8.9,9.9))