13 Complementos sobre distribuições de probabilidade

Agora que já nos familiarizamos com o uso das distribuições de probabilidade vamos ver alguns detalhes adicionais sobre seu funcionamento.

13.1 Probabilidades e integrais

A probabilidade de um evento em uma distribuição contínua é uma área sob a curva da distribuição. Vamos reforçar esta idéia revisitando um exemplo visto na aula anterior.

Seja X uma v.a. com distribuição N(100, 100). Para calcular a probabilidade P[X < 95] usamos o comando:

  > pnorm(95, 100, 10)

  [1] 0.3085375

Vamos agora “esquecer” o comando pnorm() e ver uma outra forma de resolver usando integração numérica. Lembrando que a normal tem a função de densidade dada por

          1           1
f(x) = √------ exp{- ---2(x - μ)2}
         2π σ2       2σ
vamos definir uma função no R para a densidade normal deste problema:
  > fn <- function(x) {
  +     fx <- (1/sqrt(2 * pi * 100)) * exp((-1/200) * (x - 100)^2)
  +     return(fx)
  + }

Para obter o gráfico desta distribuição mostrado na Figura 24 usamos o fato que a maior parte da função está no intervalo entre a média +/- três desvios padrões, portanto entre 70 e 130. Podemos então fazer como nos comandos que se seguem. Para marcar no gráfico a área que corresponde a probabilidade pedida criamos um polígono com coordenadas ax e ay definindo o perímetro desta área.

  > x <- seq(70, 130, l = 200)
  > fx <- fn(x)
  > plot(x, fx, type = "l")
  > ax <- c(70, 70, x[x < 95], 95, 95)
  > ay <- c(0, fn(70), fx[x < 95], fn(95), 0)
  > polygon(ax, ay, dens = 10)


PIC


Figura 24: Funções de densidade da N(100, 100) com a área correspondente à P[X 95].


Para calcular a área pedida sem usar a função pnorm() podemos usar a função de integração numérica. Note que esta função, diferentemente da pnorm() reporta ainda o erro de aproximação numérica.

  > integrate(fn, -Inf, 95)

  0.3085375 with absolute error < 2.1e-06

Portanto para os demais ítens do problema P[90 < X < 110] e P[X > 95] fazemos:

  > integrate(fn, 90, 110)

  0.6826895 with absolute error < 7.6e-15

  > integrate(fn, 95, +Inf)

  0.6914625 with absolute error < 8.1e-05

e os resultados acima evidentemente coincidem com os obtidos anterioriormente usando pnorm().

Note ainda que na prática não precisamos definir e usar a função fn pois ela fornece o mesmo resultado que a função dnorm().

13.2 Distribuição exponencial

A função de densidade de probabilidade da distribuição exponencial com parâmetro λ e denotada Exp(λ) é dada por:

       {  1 -x∕λ
f(x) =    λe       para x ≥ 0
          0       para x < 0

Seja uma variável X com distribuição exponencial de parâmetro λ = 500. Calcular a probabilidade P[X 400].

A solução analítica pode ser encontrada resolvendo

              ∫  ∞          ∫  ∞ 1
P [X  ≥ 400 ] =     f(x)dx =      --e-x∕λdx
                400            400 λ
que é uma integral que pode ser resolvida analiticamente. Fica como exercício encontrar o valor da integral acima.

Para ilustrar o uso do R vamos também obter a resposta usando integração numérica. Para isto vamos criar uma função com a expressão da exponencial e depois integrar no intervalo pedido e este resultado deve ser igual ao encontrado com a solução analítica.

  > fexp <- function(x, lambda = 500) {
  +     fx <- ifelse(x < 0, 0, (1/lambda) * exp(-x/lambda))
  +     return(fx)
  + }
  > integrate(fexp, 400, Inf)

  0.449329 with absolute error < 5e-06

Note ainda que poderíamos obter o mesmo resultado simplesmente usando a função pexp() com o comando pexp(400, rate=1/500, lower=F), onde o argumento corresponde a 1∕λ na equação da exponencial.

A Figura 25 mostra o gráfico desta distribuição com indicação da área correspondente à probabilidade pedida. Note que a função é positiva no intervalo (0, +) mas para fazer o gráfico consideramos apenas o intervalo (0, 2000).

  > x <- seq(0, 2000, l = 200)
  > fx <- dexp(x, rate = 1/500)
  > plot(x, fx, type = "l")
  > ax <- c(400, 400, x[x > 400], 2000, 2000)
  > ay <- c(0, dexp(c(400, x[x > 400], 2000), 1/500), 0)
  > polygon(ax, ay, dens = 10)


  > x <- seq(0, 2000, l = 200)
  > fx <- dexp(x, rate = 1/500)
  > plot(x, fx, type = "l")
  > ax <- c(400, 400, x[x > 400], 2000, 2000)
  > ay <- c(0, dexp(c(400, x[x > 400], 2000), 1/500), 0)
  > polygon(ax, ay, dens = 10)

PIC


Figura 25: Função de densidade da Exp(500) com a área correspondente à P[X 400].


13.3 Esperança e Variância

Sabemos que para a distribuição exponencial a esperança E[X] = 0xf(x)dx = λ e a variância V ar[X] = 0(x - E[X])2f(x)dx = λ2 pois podem ser obtidos analiticamente.

Novamente para ilustrar o uso do R vamos “esquecer” que conhecemos estes resultados e vamos obtê-los numericamente. Para isto vamos definir funções para a esperança e variância e fazer a integração numérica.

  > e.exp <- function(x, lambda = 500) {
  +     ex <- x * (1/lambda) * exp(-x/lambda)
  +     return(ex)
  + }
  > integrate(e.exp, 0, Inf)

  500 with absolute error < 0.00088

  > ex <- integrate(e.exp, 0, Inf)$value
  > ex

  [1] 500

  > v.exp <- function(x, lambda = 500, exp.x) {
  +     vx <- ((x - exp.x)^2) * (1/lambda) * exp(-x/lambda)
  +     return(vx)
  + }
  > integrate(v.exp, 0, Inf, exp.x = ex)

  250000 with absolute error < 6.9

13.4 Gerador de números aleatórios

A geração da amostra depende de um gerador de números aleatórios que é controlado por uma semente (seed em inglês). Cada vez que o comando rnorm() é chamado diferentes elementos da amostra são produzidos, porque a semente do gerador é automaticamente modificada pela função. Em geral o usuário não precisa se preocupar com este mecanismo. Mas caso necessário set.seed() pode ser usada para controlar o comportamento do gerador de números aleatórios. Esta função define o valor inicial da semente que é mudado a cada geração subsequente de números aleatórios. Portanto para gerar duas amostras idênticas basta usar set.seed() conforme ilustrado abaixo.

  > set.seed(214)
  > rnorm(5)

  [1] -0.46774980  0.04088223  1.00335193  2.02522505  0.30640096

  > rnorm(5)

  [1]  0.4257775  0.7488927  0.4464515 -2.2051418  1.9818137

  > set.seed(214)
  > rnorm(5)

  [1] -0.46774980  0.04088223  1.00335193  2.02522505  0.30640096

Nos comandos acima mostramos que depois da primeira amostra ser retirada a semente é mudada e por isto os elementos da segunda amostra são diferentes dos da primeira. Depois retornamos a semente ao seu estado original a a próxima amostra tem portanto os mesmos elementos da primeira.

Para saber mais sobre geração de números aleatórios no R veja |help(.Random.seed)| e |help(set.seed)|

13.5 Argumentos vetoriais e lei da reciclagem

As funções de probabilidades aceitam também vetores em seus argumentos conforme ilustrado nos exemplo abaixo.

  > qnorm(c(0.05, 0.95))

  [1] -1.644854  1.644854

  > rnorm(4, mean = c(0, 10, 100, 1000))

  [1]   0.4257775  10.7488927 100.4464515 997.7948582

  > rnorm(4, mean = c(10, 20, 30, 40), sd = c(2, 5))

  [1] 13.963627  6.872238 28.553964 35.584654

Note que no último exemplo a lei da reciclagem foi utilizada no vetor de desvios padrão, i.e. os desvios padrão utilizados foram (2, 5, 2, 5).

13.6 Aproximação pela Normal

Nos livros texto de estatística podemos ver que as distribuições binomial e Poisson podem ser aproximadas pela normal. Isto significa que podemos usar a distribuição normal para calcular probabilidades aproximadas em casos em que seria “trabalhoso” calcular as probabilidades exatas pela binomial ou Poisson. Isto é especialmente importante no caso de usarmos calculadoras e/ou tabelas para calcular probabilidades. Quando usamos um computador esta aproximação é menos importante, visto que é fácil calcular as probabilidades exatas com o auxílio do computador. De toda forma vamos ilustrar aqui este resultado.

Vejamos como fica a aproximação no caso da distribuição binomial. Seja X ~ B(n,p). Na prática, em geral a aproximação é considerada aceitável quando np 5 e n(1 - p) 5 e sendo tanto melhor quanto maior for o valor de n. A aproximação neste caso é de que X ~ B(n,p) N(np,np(1 - p)).

Seja X ~ B(10, 12) e portanto com a aproximação X N(5, 2.5). A Figura 26 mostra o gráfico da distribuição binomial e da aproximação pela normal.

  > xb <- 0:10
  > px <- dbinom(xb, 10, 0.5)
  > plot(xb, px, type = "h")
  > xn <- seq(0, 10, len = 100)
  > fx <- dnorm(xn, 5, sqrt(2.5))
  > lines(xn, fx)


PIC


Figura 26: Função de probabilidade da B(10, 12) e a aproximação pela N(5, 2.5).


Vamos também calcular as seguintes probabilidades exatas e aproximadas, lembrando que ao usar a aproximação pela normal devemos usar a correção de continuidade e/ou somando e subtraindo 0.5 ao valor pedido.

13.7 Exercícios

1.
(Bussab & Morettin, pag. 198, ex. 51)
A função de densidade de probabilidade de distribuição Weibull é dada por:
        {          λ
          λx λ-1e-x    para x ≥ 0
f (x) =   0           para x <  0

(a)
Obter E[X] para λ = 2. Obter o resultado analitica e computacionalmente.
Dica: para resolver você vai precisar da definição da função Gama:
       ∫ ∞
Γ (a) =    xa- 1e-xdx
        0
(b)
Obter E[X] para λ = 5.
(c)
Obter as probabilidades:
  • P[X > 2]
  • P[1.5 < X < 6]
  • P[X < 8]