4

Boa tarde.

É a minha primeira intervenção no stack, sou iniciante no R, e as minhas dúvidas são bastantes básicas.

Preciso gerar uma amostra de 1000 observações de uma função de distribuição de W.
W é uma variável aleatória discreta, que assume os valores de 1 a 6, representada pelos lados de um dado, dado este que está viciado e cuja função probabilidade é p1=0.25, p2=0.16, p3=0.18, p4=0.17, p5=0.14, p6=0.10.

Como posso escrever essa função em R?

Obrigado

2 Respostas 2

4

Para amostrar de uma variável discreta que toma um número finito de valores, pode-se usar a função base sample.
Antes de correr a função sample ou outra função que gere números pseudo-aleatórios, é sempre melhor chamar set.seed.

set.seed(7228)

W <- 1:6
p <- c(0.25, 0.16, 0.18, 0.17, 0.14, 0.10)
w <- sample(W, 1000, replace = TRUE, prob = p)

head(w, n = 20)
#[1] 4 4 1 3 6 3 5 4 1 1 2 4 4 4 4 2 6 6 5 6

Ver se as proporções do resultado são semelhantes às probabilidades dadas.

tw <- table(w)
print(tw/sum(tw), digits = 2)
#w
#   1    2    3    4    5    6 
#0.23 0.16 0.16 0.19 0.14 0.12

Não parecem ser muito diferentes. Se for preciso, sempre se pode correr um teste de Kolmogorov-Smirnov, uma vez que tanto p como as proporções da amostra vêm de uma distribuição contínua.

ks.test(p, tw/sum(tw))

Com um p-value = 0.8928 concui-se que sim, que as distribuições não são significativamente diferentes.

2

A forma de você resolver esse problema é utilizando a função de distribuição acumulada inversa dessa variável.

Você sabe as probabilidades da variável X assumir cada um dos possíveis valores (de 1 a 6), então é possível construir a função de distribuição acumulada (F(x) = P(X <= x)):

F(x) =

  • 0,25 se X <= 1;
  • 0,41 se X <= 2;
  • 0,59 se X <= 3;
  • 0,76 se X <= 4;
  • 0,90 se X <= 5;
  • 1,00 se X <= 6

Então, utilizando o conceito de função inversa F(x)^(-1), gere uma probabilidade e tenha o valor da variável.

F(x)^(-1):

  • 1 se 0 < p < 0,25;
  • 2 se 0,25 < p < 0,41;
  • 3 se 0,41 < p < 0,59;
  • 4 se 0,59 < p < 0,76;
  • 5 se 0,76 < p < 0,90;
  • 6 se 0,90 < p < 1,00

No R você gera "uma probabilidade" com a função runif(). Sendo assim, segue código que gera valores dessa sua variável aleatória:

gerar.VA <- function(n = 1){
  #P(X=1) = 0,25; P(X=2) = 0,16; P(X=3) = 0,18 ;P(X=4) = 0,17; P(X=5) = 0,14; P(X=6) = 0,10
  resultado <- NULL
  aux <- runif(n)
  for(i in 1:n){
    if(aux[i] < 0.25){
      resultado[i] <- 1
    } else{
      if(aux[i] < 0.25 + 0.16){
        resultado[i] <- 2
      } else{
        if(aux[i] < 0.25 + 0.16 + 0.18){
          resultado[i] <- 3
        } else{
          if(aux[i] < 0.25 + 0.16 + 0.18 + 0.17){
            resultado[i] <- 4
          } else{
            if(aux[i] < 0.25 + 0.16 + 0.18 + 0.17 + 0.14){
              resultado[i] <- 5
            } else{
              resultado[i] <- 6
            }
          }
        }
      }
    }
  }
  return(resultado)
}
gerar.VA(n = 1000)

Você deve fazer log-in para responder a esta pergunta.

Esta não é a resposta que você está procurando? Pesquise outras perguntas com a tag .