Estatística Básica - Bussab e Morettin (2017)
Comandos R para análises estatísticas
Capítulo 6: Variáveis Aleatórias Discretas
6.1 Introdução
Neste capítulo são descritos os comandos para utilização de variáveis aleatórias R. No R, é possível simular essas variáveis, obter sua distribuição acumulada e seus quantis. Por exemplo, para uma variável aleatória com distribuição binomial, podemos utilizar os seguintes comandos:
Usage
dbinom(x, size, prob, log = FALSE)
pbinom(q, size, prob, lower.tail = TRUE, log.p = FALSE)
qbinom(p, size, prob, lower.tail = TRUE, log.p = FALSE)
rbinom(n, size, prob)
Principais argumentos:
x : valor observador para o qual se deseja calcular o valor da distribuição
q : vetor de quantis
p : vetor de probabilidades
n : Número de observações a serem simuladas
size : parâmetro da distribuição: número de ensaios
prob : parâmetro da distribuição: probabilidade de sucesso de cada ensaio
lower.tail : Variável TRUE/FALSE. Se TRUE, P[X ≤ x] para TRUE e P[X > x].
No R, os comandos para trabalhar com variáveis aleatórias seguem sempre o mesmo padrão, i.e. sua nomenclatura e parâmetros são semelhantes. Com relação à nomenclatura, note que no exemplo anterior, a referência à distribuição é sempre precedida de uma letra, d para densidade, p para função distribuição, q para quantis e r para simulação. Dessa forma, se quisermos simular uma distribuição de poisson, o comando será rpois(...), para acessar os quantis de uma distribuição geométrica, qgeom(...) e assim por diante. Os sufixos para as distribuições discretas mais conhecidas são apresentados na Tabela 6.A.
Tabela 6.A: principais distribuições discretas e seus sufixos no R.
| Distribuição | Sufixo | Parâmetros |
|---|---|---|
| Uniforme | -unif |
min, max |
| Binomial | -binom |
size, prob |
| Geométrica | -geom |
prob |
| Poison | -pois |
lambda |
| Hipergeométrica | -hyper |
m, n, k |
Figura 6.1
Para construir a Figura 6.1, utilizaremos o pacote diagram e o comando plotmat.
Este comando constrói o diagrama de acordo com uma matriz A que representará os nós do gráfico. Temos 13 nós na Figura 6.1, logo nossa matriz terá 13 colunas e 13 linhas, com as colunas representando a origem das setas e as linhas, o destino. Assim, se identificarmos de 1 a 13 os nós que temos na figura da forma a seguir
então, para imprimir a seta do primeiro ramo com 0.8, por exemplo, colocaremos na posição (1,4) da matriz A o valor 0.8.
A<-matrix(0,13,13)
A[4,1]<-0.8
A## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13]
## [1,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [2,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [3,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [4,] 0.8 0 0 0 0 0 0 0 0 0 0 0 0
## [5,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [6,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [7,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [8,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [9,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [10,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [11,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [12,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
## [13,] 0.0 0 0 0 0 0 0 0 0 0 0 0 0
Nossa matriz A para construir o diagrama da Figura 6.1 é apresentada na próxima sequência de comandos. Dois outros parâmetros fazem-se necessários: uma matriz com duas colunas e 13 linhas, cada uma delas com a posição relativa (entre 0 e 1) das coordenadas x e y de cada um dos nós, em ordem(objeto pos), e um vetor com os nomes de cada um dos nós(objeto knots)
A<-matrix(0,13,13) # Matriz de setas: colunas são as origens e linhas, os destinos
A[2,1]<-0.1
A[3,1]<-0.1
A[4,1]<-0.8
A[5,2]<-0.1
A[6,2]<-0.2
A[7,2]<-0.7
A[8,3]<-0.1
A[9,3]<-0.2
A[10,3]<-0.7
A[11,4]<-0.1
A[12,4]<-0.2
A[13,4]<-0.7
col <- A
col[] <- "blue"
pos=rbind(c(0.05,0.5),
c(0.3,0.2),
c(0.3,0.5),
c(0.3,0.8),
c(0.7,0.1),
c(0.7,0.2),
c(0.7,0.3),
c(0.7,0.4),
c(0.7,0.5),
c(0.7,0.6),
c(0.7,0.7),
c(0.7,0.8),
c(0.7,0.9)) # Matriz de coordenadas dos nós
knots=c(" ","L","C","B","C 0.01","L 0.02","B 0.07","C 0.01","L 0.02","B 0.07","C 0.08","L 0.16","B 0.56") # labels dos nós.
pp <- plotmat(A, pos = pos, name = knots,
box.cex = c(2,2,2,2,1,1,1,1,1,1,1,1,1), box.type = "none", box.prop = 1.0, # Parâmetros do nó: tamanho da letra, tipo de borda e proporção
cex.txt = 1, dtext=-0.4, #parâmetros do texto das setas: tamanho e posição
arr.lcol = col, arr.col = col, arr.pos = 0.8, arr.type = "simple", arr.len = 0.2,
curve = 0, segment.from = 0.8, segment.to = 0.1,lwd = 2, # Parametros complementares de cor, posição da ponta da seta, tipo de seta, espessura da linha, curvatura, comprimento da linha
main = "Diagrama de árvore: Cilindro X Esfera " # titulo
)
Figura 6.1: Diagrama em árvore para o Exemplo 6.1.
Figura 6.4
Da mesma forma que fizemos na Figura 6.1, precisamos montar os parâmetros A, pos e knots
A<-matrix(0,7,7) # Matriz de setas: colunas são as origens e linhas, os destinos
A[2,1]<-"2/5"
A[3,1]<-"3/5"
A[4,2]<-"2/4"
A[5,2]<-"2/4"
A[6,3]<-"3/4"
A[7,3]<-"1/4"
col <- A
col[] <- "blue"
pos=rbind(c(0.05,0.5),
c(0.3,0.33),
c(0.3,0.67),
c(0.7,0.2),
c(0.7,0.4),
c(0.7,0.6),
c(0.7,0.8)) # Matriz de coordenadas dos nós
knots=c(" ","V","B","V","B","V","B") # labels dos nós.
pp <- plotmat(A, pos = pos, name = knots,
box.cex = c(2,2,2,1,1,1,1), box.type = "none", box.prop = 1.0, # Parâmetros do nó: tamanho da letra, tipo de borda e proporção
cex.txt = 1, dtext=0, #parâmetros do texto das setas: tamanho e posição
arr.lcol = col, arr.col = col, arr.pos = 0.7, arr.type = "simple", arr.len = 0.3,
curve = 0, segment.from = 0.8, segment.to = 0.1,lwd = 2, # Parametros complementares de cor, posição da ponta da seta, tipo de seta, espessura da linha, curvatura, comprimento da linha
main = "Diagrama de árvore" # titulo
)
Figura 6.4: Diagrama em árvore para o Exemplo 6.3.
Figura 6.5
Da mesma forma que fizemos na Figura 6.1, precisamos montar os parâmetros A, pos e knots
A<-matrix(0,7,7) # Matriz de setas: colunas são as origens e linhas, os destinos
A[2,1]<-"1/2"
A[3,1]<-"1/2"
A[4,2]<-"1/2"
A[5,2]<-"1/2"
A[6,3]<-"1/2"
A[7,3]<-"1/2"
col <- A
col[] <- "blue"
pos=rbind(c(0.05,0.5),
c(0.3,0.33),
c(0.3,0.67),
c(0.7,0.2),
c(0.7,0.4),
c(0.7,0.6),
c(0.7,0.8)) # Matriz de coordenadas dos nós
knots=c(" ","R","C","R","C","R","C") # labels dos nós.
pp <- plotmat(A, pos = pos, name = knots,
box.cex = c(2,2,2,1,1,1,1), box.type = "none", box.prop = 1.0, # Parâmetros do nó: tamanho da letra, tipo de borda e proporção
cex.txt = 1, dtext=0, #parâmetros do texto das setas: tamanho e posição
arr.lcol = col, arr.col = col, arr.pos = 0.7, arr.type = "simple", arr.len = 0.3,
curve = 0, segment.from = 0.8, segment.to = 0.1,lwd = 2, # Parametros complementares de cor, posição da ponta da seta, tipo de seta, espessura da linha, curvatura, comprimento da linha
main = "Diagrama de árvore" # titulo
)
Figura 6.5: Diagrama em árvore para o Exemplo 6.4.
Figura 6.7
O gráfico da Figura 6.7 é um gráfico de dispersão e sua construção foi discutida no Capítulo 3. Relembrando, se temos a Tabela 6.3:
| x | p(x) |
|---|---|
| 15 | 0.56 |
| 10 | 0.23 |
| 5 | 0.02 |
| -5 | 0.19 |
| Total | 1.00 |
Imprimiremos o gráfico da distribuição da v.a. X = lucro por montagem
lucro<-data.frame(x=c(15,10,5,-5), px=c(0.56 , 0.23, 0.02, 0.19))
plot(lucro$x,lucro$px, pch=16, col="darkblue", xlab="Lucro", ylab="p(x)", bty="n", ylim=c(0,0.7)) 
# Alternativamente, pode-se modificar a posição do eixo Y utilizando o comando:
# axis(2,at=seq(0,0.6,0.1),labels=seq(0,0.6,0.1),pos=0)
# Se desejar isso, é necessário adicionar a opção `yaxt="n"` no comando que gera o gráfico.Figura 6.7: Gráfico de p(x): distribuição da v.a. X = lucro por montagem.
Figura 6.8
Para imprimir o gráfico de uma função escada, utilize o comando stepfun associado com plot:
### Sintaxe
stepfun(x, y,...)
x: Vetor contendo as posições em que ocorrem os saltos (para f.d.a., x deve estar ordenado)
y: Vetor contendo a altura do degrau que começa na posição (para f.d.a., y representa as probabilidades acumuladas)
o.lucro<-lucro[order(lucro$x),] # Ordenando o dataset
o.lucro<-cbind(o.lucro,fda_x=cumsum(o.lucro$px)) # Calculando a fda de X
plot(stepfun(o.lucro$x,c(0,o.lucro$fda_x)),
main="f.d.a para a v.a. X = Lucro por montagem",
xlab="Lucro", ylab="p(x)", bty="n",col="darkblue",pch=20)
lines(x=c(5,5),y=c(0,o.lucro$fda_x[1]),lty=2,col="gray")
lines(x=c(10,10),y=c(0,o.lucro$fda_x[2]),lty=2,col="gray")
lines(x=c(15,15),y=c(0,o.lucro$fda_x[3]),lty=2,col="gray")
Figura 6.8: f.d.a. para a v.a. X = lucro por montagem.
Tabela 6.12
Para a Tabela 6.12, utilizaremos o comando dbinom(...):
x = c(0,1,2,3)
px<-dbinom(x=x,size=3,p=1/2)
table612<-data.frame(x,px)| x | px |
|---|---|
| 0 | 0.125 |
| 1 | 0.375 |
| 2 | 0.375 |
| 3 | 0.125 |
Figura 6.12
plot(table612$x,table612$px, pch=16, col="darkblue", xlab="Lucro", ylab="p(x)", bty="n") 
Figura 6.12: Gráfico da f.p. p(x) para n=3 e p=1/2.
Exemplo 6.15
A probabilidade para uma v.a. hipergeométrica pode ser calculada da seguinte forma:
# x: numero de peças defeituosas desejadas
# m: Número total de peças defeituosas
# n: Número total de peças não defeituosas
# k: Número de bolas retiradas
dhyper(x=0 , m = 10, n = 90, k = 5)## [1] 0.58375Exemplo 6.17
Para Ilustrar a qualidade da aproximação da distribuição binomial pela distribuição de Poisson, considere as duas situações apresentadas na Tabela 6.13a a seguir. Nela, as duas primeiras colunas, pbinomial.5 e ppoisson.5 foram calculadas com p=0.5 ao passo que as duas últimas pbinomial.005 e ppoisson.005 foram calculadas com p=0.005.
x<-0:9
p.5=0.5
p.005=0.005
n=10
pbinomial.5<-dbinom(x,n,p.5)
ppoisson.5<-dpois(x,lambda=n*p.5)
pbinomial.005<-dbinom(x,n,p.005)
ppoisson.005<-dpois(x,lambda=n*p.005)
table613a<-data.frame(x,pbinomial.5,ppoisson.5,pbinomial.005,ppoisson.005)
table613a<-rbind(table613a,c(10,1-sum(pbinomial.5),1-sum(ppoisson.5),1-sum(pbinomial.005),1-sum(ppoisson.005)))| x | pbinomial.5 | ppoisson.5 | pbinomial.005 | ppoisson.005 |
|---|---|---|---|---|
| 0 | 0.00098 | 0.00674 | 0.95111 | 0.95123 |
| 1 | 0.00977 | 0.03369 | 0.04779 | 0.04756 |
| 2 | 0.04395 | 0.08422 | 0.00108 | 0.00119 |
| 3 | 0.11719 | 0.14037 | 0.00001 | 0.00002 |
| 4 | 0.20508 | 0.17547 | 0.00000 | 0.00000 |
| 5 | 0.24609 | 0.17547 | 0.00000 | 0.00000 |
| 6 | 0.20508 | 0.14622 | 0.00000 | 0.00000 |
| 7 | 0.11719 | 0.10444 | 0.00000 | 0.00000 |
| 8 | 0.04395 | 0.06528 | 0.00000 | 0.00000 |
| 9 | 0.00977 | 0.03627 | 0.00000 | 0.00000 |
| 10 | 0.00098 | 0.03183 | 0.00000 | 0.00000 |
Figura 6.13a: Comparação entre aproximações da distribuição binomial pela Poisson para p=0.5 e p=0.005
Exemplo 6.19
x<-0:12
px<-dbinom(x,size=14,p=0.3)
fdax<-cumsum(px)
quadro61<-data.frame(x,px,fdax)| x | px | fdax |
|---|---|---|
| 0 | 0.00678 | 0.00678 |
| 1 | 0.04069 | 0.04748 |
| 2 | 0.11336 | 0.16084 |
| 3 | 0.19433 | 0.35517 |
| 4 | 0.22903 | 0.58420 |
| 5 | 0.19631 | 0.78052 |
| 6 | 0.12620 | 0.90672 |
| 7 | 0.06181 | 0.96853 |
| 8 | 0.02318 | 0.99171 |
| 9 | 0.00662 | 0.99833 |
| 10 | 0.00142 | 0.99975 |
| 11 | 0.00022 | 0.99997 |
| 12 | 0.00002 | 1.00000 |
x<-0:17
px<-dpois(x,lambda=5.2)
fdax<-cumsum(px)
quadro62<-data.frame(x,px,fdax)| x | px | fdax |
|---|---|---|
| 0 | 0.00552 | 0.00552 |
| 1 | 0.02869 | 0.03420 |
| 2 | 0.07458 | 0.10879 |
| 3 | 0.12928 | 0.23807 |
| 4 | 0.16806 | 0.40613 |
| 5 | 0.17479 | 0.58091 |
| 6 | 0.15148 | 0.73239 |
| 7 | 0.11253 | 0.84492 |
| 8 | 0.07314 | 0.91806 |
| 9 | 0.04226 | 0.96033 |
| 10 | 0.02198 | 0.98230 |
| 11 | 0.01039 | 0.99269 |
| 12 | 0.00450 | 0.99719 |
| 13 | 0.00180 | 0.99899 |
| 14 | 0.00067 | 0.99966 |
| 15 | 0.00023 | 0.99989 |
| 16 | 0.00008 | 0.99997 |
| 17 | 0.00002 | 0.99999 |