O programa R inclui funcionalidade para operações com distribuições de probabilidades. Para cada distribuição há 4 operações básicas indicadas pelas letras:
Para usar as funções deve-se combinar uma das letras acima com uma abreviatura do nome da distribuição. Por exemplo, para calcular probabilidades usamos: pnorm() para distribuição normal, pexp() para exponencial, pbinom() para binomial, ppois() para Poisson e assim por diante. Vamos ver com mais detalhes algumas distribuições de probabilidades.
A funcionalidade para distribuição normal é implementada por argumentos que combinam as letras
acima com o termo norm. Vamos ver alguns exemplos com a distribuição normal padrão. Por default
as funções assumem a distribuição normal padrão \(N(\mu =0, \sigma ^2=1)\).
> dnorm(-1)
[1] 0.2419707
> pnorm(-1)
[1] 0.1586553
> qnorm(0.975)
[1] 1.959964
> rnorm(10)
[1] 0.20655549 -0.64716079 -0.25471700 0.03371804 0.11033425 0.98218218 [7] -0.60410167 0.22780774 0.66930840 1.07551727
O primeiro valor corresponde ao valor da densidade da normal \[ f(x) = \frac {1}{\sqrt {2 \pi \sigma ^2}}\exp \{-\frac {1}{2\sigma ^2}(x - \mu )^2\}\] com parâmetros \((\mu =0, \sigma ^2=1)\) no ponto \(-1\). Portanto, o mesmo valor seria obtido substituindo \(x\) por \(-1\) na expressão da normal padrão:
> (1/sqrt(2*pi)) * exp((-1/2)*(-1)^2)
[1] 0.2419707
A função pnorm(-1) calcula a probabilidade acumulada \(P(X \leq -1)\). O comando qnorm(0.975) calcula o valor de um quantil \(a\) tal que \(P(X \leq a)=0.975\). Finalmente, o comando rnorm(10) gera uma amostra de 10 elementos da normal padrão. Note que os valores que você obtém rodando este comando podem ser diferentes dos mostrados acima.
As funções acima possuem argumentos adicionais, para os quais valores padrão (default) foram assumidos, e podem ser modificados. Usamos args() para ver os argumentos de uma função e help() para visualizar a documentação detalhada:
> args(rnorm)
function (n, mean = 0, sd = 1) NULL
As funções relacionadas à distribuição normal possuem os argumentos mean e sd para definir média e desvio padrão da distribuição que podem ser modificados como nos exemplos a seguir. Note nestes exemplos que os argumentos podem ser passados de diferentes formas.
> qnorm(0.975, mean = 100, sd = 8)
[1] 115.6797
> qnorm(0.975, m = 100, s = 8)
[1] 115.6797
> qnorm(0.975, 100, 8)
[1] 115.6797
Para informações mais detalhadas pode-se usar help(). O comando
> help(rnorm)
irá exibir em uma janela a documentação da função que pode também ser chamada com
?rnorm. Note que ao final da documentação são apresentados exemplos que podem ser
rodados pelo usuário e que auxiliam na compreensão da funcionalidade. Note também
que as 4 funções relacionadas à distribuição normal são documentadas conjuntamente,
portanto help(rnorm), help(qnorm), help(dnorm) e help(pnorm) irão exibir a mesma
documentação.
Problemas de cálculo de probabilidades usuais, para os quais utilizávamos tabelas estatísticas podem ser facilmente obtidos como no exemplo a seguir.
Seja \(X\) uma v.a. com distribuição \(N(100, 100)\). Calcular as probabilidades:
\( P[X < 95] \)
\( P[90 < X < 110] \)
\( P[X > 95] \)
Calcule estas probabilidades de forma usual, usando a tabela da normal. Depois compare com os resultados fornecidos pelo R. Os comandos do R para obter as probabilidades pedidas são:
> pnorm(95, 100, 10)
[1] 0.3085375
> pnorm(110, 100, 10) - pnorm(90, 100, 10)
[1] 0.6826895
> 1 - pnorm(95, 100, 10)
[1] 0.6914625
> pnorm(95, 100, 10, lower=FALSE)
[1] 0.6914625
Note que a última probabilidade foi calculada de duas formas diferentes, a segunda usando o argumento lower que implementa um algoritmo de cálculo de probabilidades mais estável numericamente.
A seguir vamos ver comandos para fazer gráficos de distribuições de probabilidade. Vamos fazer gráficos de funções de densidade e de probabilidade acumulada. Estude cuidadosamente os comandos a seguir e verifique os gráficos por eles produzidos. A Figura 20 mostra gráficos da densidade (esquerda) e probabilidade acumulada (direita) da normal padrão. Para fazer o gráfico consideramos valores de \(X\) entre -3 e 3 que correspondem a +/- três desvios padrões da média, faixa que concentra 99,73% da massa de probabilidade da distribuição normal.
> plot(dnorm, -3, 3) > plot(pnorm, -3, 3)
A Figura 21 mostra gráficos da densidade (esquerda) e probabilidade acumulada (direita) da \(N(100, 64)\). Para fazer estes gráficos tomamos uma sequência de valores de \(x\) entre 70 e 130 e para cada um deles calculamos o valor das funções \(f(x)\) e \(F(x)\). Depois unimos os pontos \((x,f(x))\) em um gráfico e \((x,F(x))\) no outro.
> x <- seq(70, 130, len=100) > fx <- dnorm(x, 100, 8) > plot(x, fx, type='l') > Fx <- pnorm(x, 100, 8) > plot(x, Fx, type='l')
Note que, alternativamente, os mesmos gráficos poderiam ser produzidos com os comandos a seguir.
> plot(function(x) dnorm(x, 100, 8), 70, 130) > plot(function(x) pnorm(x, 100, 8), 70, 130)
ou ainda
> curve(dnorm(x, 100, 8), 70, 130) > curve(pnorm(x, 100, 8), 70, 130)
Comandos usuais do R podem ser usados para modificar a aparência dos gráficos. Por exemplo, podemos incluir títulos e mudar texto dos eixos conforme mostrado no gráfico da esquerda da Figura 22 e nos dois primeiros comandos abaixo. Os demais comandos mostram como colocar diferentes densidades em um mesmo gráfico como ilustrado à direita da mesma Figura.
> plot(dnorm, -3, 3, xlab='valores de X', ylab='densidade de probabilidade') > title('Distribuição Normal\nX ~ N(100, 64)') > > curve(dnorm(x, 100, 8), 60, 140, ylab='f(x)') > curve(dnorm(x, 90, 8), 60, 140, add=T, col=2) > curve(dnorm(x, 100, 15), 60, 140, add=T, col=3) > legend("topright", c("N(100,64)","N(90,64)","N(100,225)"), fill=1:3)
Cálculos para a distribuição binomial são implementados combinando as letras básicas vistas
acima com o termo binom. Vamos primeiro investigar argumentos e documentação com args() e
dbinom().
> args(dbinom)
function (x, size, prob, log = FALSE) NULL
> help(dbinom)
Seja \(X\) uma v.a. com distribuição Binomial com \(n=10\) e \(p=0.35\). Vamos ver os comandos do R para:
fazer o gráfico da função de densidade
idem para a função de probabilidade
calcular \(P[X = 7]\)
calcular \(P[X < 8] = P[X \leq 7]\)
calcular \(P[X \geq 8] = P[X > 7]\)
calcular \(P[3 < X \leq 6] = P[4 \leq X < 7]\)
Note que sendo uma distribuição discreta de probabilidades os gráficos são diferentes dos obtidos para distribuição normal e os cálculos de probabilidades devem considerar as probabilidades nos pontos. Os gráficos das funções de densidade e probabilidade são mostrados na Figura 23.
> x <- 0:10 > fx <- dbinom(x, 10, 0.35) > plot(x, fx, type='h') > Fx <- pbinom(x, 10, 0.35) > plot(x, Fx, type='s')
As probabilidades pedidas são obtidas com os comandos a seguir.
> dbinom(7, 10, 0.35)
[1] 0.02120302
> pbinom(7, 10, 0.35)
[1] 0.9951787
> sum(dbinom(0:7, 10, 0.35))
[1] 0.9951787
> 1-pbinom(7, 10, 0.35)
[1] 0.004821265
> pbinom(7, 10, 0.35, lower=F)
[1] 0.004821265
> pbinom(6, 10, 0.35) - pbinom(3, 10, 0.35)
[1] 0.4601487
> sum(dbinom(4:6, 10, 0.35))
[1] 0.4601487
Para a distribuição uniforme contínua usa-se as funções *unif() onde * deve ser \(p\), \(q\), \(d\) ou \(r\) como mencionado anteriormente. Nos comandos a seguir inspecionamos os argumentos, sorteamos 5 valores da \(U(0,1)\) e calculamos a probabilidade acumulada até 0,75.
> args(runif)
function (n, min = 0, max = 1) NULL
> runif(5)
[1] 0.7163657 0.5931028 0.9449806 0.1846381 0.5385011
> punif(0.75)
[1] 0.75
Portanto, o default é uma distribuição uniforme no intervalo \([0,1]\) e os argumentos opcionais são min e max. Por exemplo, para simular 5 valores de \(X \sim U(5, 20)\) usamos:
> runif(5, min=5, max=20)
[1] 8.737943 10.087838 17.835511 13.626632 16.895792
Não há entre as funções básicas do R uma função específica para a distribuição uniforme discreta com opções de prefixos \(r,d,p\) e \(q\). Provavelmente isso é devido a sua simplicidade, embora algumas outras funções possam ser usadas. Por exemplo para sortear números pode-se usar sample(). No exemplo a seguir são sorteados 15 valores de uma uniforme discreta com valores (inteiros) entre 1 e 10 (\(X \sim U_d(1,10)\)).
> sample(1:10, 15, rep=T)
[1] 3 9 7 6 1 6 6 5 8 3 1 2 5 6 9
A função sample() não é restrita à distribuição uniforme discreta, podendo ser usada para sorteios, com ou sem reposição (conforme indicado pelo argumento replace, com default sem reposição) e ainda com a possibilidade de associar diferentes probabilidades a cada elemento (argumento prob, default probabilidades iguais para os elementos).
> args(sample)
function (x, size, replace = FALSE, prob = NULL) NULL
Vejamos alguns exemplos:
sorteio de 3 números entre os inteiros de 0 a 20
> sample(0:20, 3)
[1] 16 14 4
sorteio de 5 números entre os elementos de um certo vetor
> x <- c(23, 34, 12, 22, 17, 28, 18, 19, 20, 13, 18) > sample(x, 5)
[1] 12 34 23 20 28
sorteio de 10 números entre os possíveis resultados do lançamento de um dado, com reposição
> sample(1:6, 10, rep=T)
[1] 5 4 6 5 1 2 1 4 4 5
idem ao anterior, porém agora com a probabilidade de cada face proporcional ao valor da face.
> sample(1:6, 10, prob=1:6, rep=T)
[1] 6 4 6 2 4 6 2 3 5 6
O último exemplo ilustra ainda que os valores passados para o argumento prob não precisam ser probabilidades, são apenas entendidos como pesos. A própria função trata isto internamente fazendo a ponderação adequada.
Nos exercícios abaixo iremos também usar o R como uma calculadora estatística para resolver alguns exemplos/exercícios de probabilidade tipicamente apresentados em um curso de estatística básica.
Os exercícios abaixo com indicação de página foram retirados de:
Magalhães, M.N. & Lima, A.C.P. (2001) Noções de Probabilidade e Estatística. 3 ed. São
Paulo, IME-USP. 392p.
(Ex 1, pag 67) Uma moeda viciada tem probabilidade de cara igual a 0.4. Para quatro lançamentos independentes dessa moeda, estude o comportamento da variável número de caras e faça um gráfico de sua função de distribuição.
(Ex 5, pag 77) Sendo \(X\) uma variável seguindo o modelo Binomial com parâmetro \(n = 15\) e \(p = 0.4\), pergunta-se:
(Ex 8, pag 193) Para \(X \sim N(90, 100)\), obtenha:
Faça os seguintes gráficos:
A probabilidade de indivíduos nascerem com certa característica é de 0,3. Para o nascimento de 5 indivíduos e considerando os nascimentos como eventos independentes, estude o comportamento da variável número de indivíduos com a característica e faça um gráfico de sua função de distribuição.
Sendo \(X\) uma variável seguindo o modelo Normal com média \(\mu = 130\) e variância \(\sigma ^2 = 64\), pergunta-se: (a) \(P(X \geq 120)\) (b) \(P(135 < X \leq 145)\) (c) \(P(X < 120 \;{\rm ou }\; X \geq 150)\)
(Ex 3.6, pag 65) Num estudo sobre a incidência de câncer foi registrado, para cada paciente com este diagnóstico o número de casos de câncer em parentes próximos (pais, irmãos, tios, filhos e sobrinhos). Os dados de 26 pacientes são os seguintes:
| Paciente | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 | 13 |
| Incidência | 2 | 5 | 0 | 2 | 1 | 5 | 3 | 3 | 3 | 2 | 0 | 1 | 1 |
| Paciente | 14 | 15 | 16 | 17 | 18 | 19 | 20 | 21 | 22 | 23 | 24 | 25 | 26 |
| Incidência | 4 | 5 | 2 | 2 | 3 | 2 | 1 | 5 | 4 | 0 | 0 | 3 | 3 |
Estudos anteriores assumem que a incidência de câncer em parentes próximos pode ser modelada pela seguinte função discreta de probabilidades:
| Incidência | 0 | 1 | 2 | 3 | 4 | 5 |
| \(p_i\) | 0.1 | 0.1 | 0.3 | 0.3 | 0.1 | 0.1 |
A distribuição da soma de duas variáveis aleatórias uniformes não é uniforme. Verifique isto gerando dois vetores \(x\) e \(y\) com distribuição uniforme \([0,1]\) com 3000 valores cada e fazendo \(z=x+y\). Obtenha o histograma para \(x\), \(y\) e \(z\). Descreva os comandos que utilizou.
(extraído de Magalhães e Lima, 2001) A resistência (em toneladas) de vigas de concreto produzidas por uma empresa, comporta-se como abaixo:
| Resistência | 2 | 3 | 4 | 5 | 6 |
| pi | 0,1 | 0,1 | 0,4 | 0,2 | 0,2 |
Simule a resistência de 5000 vigas a partir de valores gerados de uma uniforme [0,1]. (Dica: Use o comando ifelse() do R). Verifique o histograma.