Primeira vez aqui? Seja bem vindo e cheque o FAQ!
x

Como conferir os resultados do gerador de números distribuídos uniformemente de seu pacote estatístico?

0 votos
29 visitas
perguntada Set 23 em Ciência da Computação por Sergio Costa (6 pontos)  
editado Set 23 por Sergio Costa

Conferindo o gerador de números distribuídos uniformemente do pacote numpy no Python. Questão 2.1 do livro "Monte Carlo Satistical Methods", por Christian Robert and George Casella (Capítulo 2, pág 62).

Compartilhe

1 Resposta

+1 voto
respondida Set 23 por Sergio Costa (6 pontos)  
editado Set 29 por Sergio Costa

Olá, amigos,

Trata-se de resolução da questão 2.1 do livro "Monte Carlo Statistical Methods", por Christian Robert e George Casella.

A questão solicita plotar em histograma 1000 variáveis aleatórias distribuídas uniformemente. Posteriormente, pede-se para verificar a correlação entre as variáveis plotando \((X_{i}, X_{i+1})\).

Neste exemplo, \(X \sim U[0,1)\) e o pacote utilizado para gerar os números aleatórios foi numpy, conforme código a seguir:

import matplotlib.pyplot as plt
import numpy as np

def dist_uniforme(amostra, a, b, n):
    #gerando números aleatórios X ~ U[0,1)
    np.random.seed(amostra)
    x = np.random.uniform(a,b,n)
    y = x[0:n-1]
    z = x[1:n]
    densidade=1/(b-a)

    #desenho histograma
    plt.subplot(211)
    plt.ylabel('f.d.p.')
    n, bins, patches = plt.hist(x, bins='auto', density=True, facecolor='b', alpha=0.9, edgecolor='black')
    plt.title('Histograma e Estudo Correlação')
    c=len(bins)
    plt.plot(bins,np.repeat(densidade,c),linewidth=2, color='r')


    #desenho correlação
    plt.subplot(212)
    plt.plot(y, z, 'go', markersize=3)
    plt.xlabel('Xi+1')
    plt.ylabel('Xi')
    plt.grid(True)
    plt.show()

    return x, y, z

if __name__=="__main__":
    dist_uniforme(1, 0, 1, 1000)

A figura resultante dá uma ideia, embora superficial, em relação a qualidade do gerador de números aleatórios.

No painel superior, o histograma é apresentado junto com a f.d.p. aparentando um ajuste apropriado e no inferior não se apresenta uma relação clara entre \((X_{i}, X_{i+1})\).

Para efeito de comparação, fixei a seed utilizada para os números pseudo-aleatórios.

A imagem será apresentada aqui.

comentou Set 29 por Thiago Trafane (21 pontos)  
Sergio, parabéns pela solução! Tenho apenas dois comentários bem pontuais:

1) Quando você escreveu \(X \sim U[0,1)\), acredito que o correto seria \(X \sim U(0,1)\), já que aqui \(U\) representa uma função cujos argumentos são \(0\) e \(1\).

2) Quando você faz o histograma, talvez fosse interessante usar as opções bins='auto' e density=True. Com a primeira opção, os intervalos para cálculo das frequências passam a ser definidos automaticamente. Com a segunda opção, as frequências são ajustadas para representarem uma função de densidade empírica. Nesse caso, você poderia incluir a função densidade da Uniforme no gráfico do histograma, o que mostraria mais claramente a qualidade do gerador de números distribuídos uniformemente.
comentou Set 29 por Sergio Costa (6 pontos)  
Muito obrigado pelas sugestões, Thiago, certamente melhoraram o código. Corrigi, então, a programação conforme indicado, exceto a notação da distribuição. Posso estar errado e mudarei sem problemas, mas tenho a impressão de que a notação para a distribuição uniforme representa o intervalo. Como o numpy, utiliza um intervalo fechado no limite inferior e aberto no superior, escrevi daquela maneira.
comentou Set 29 por Thiago Trafane (21 pontos)  
editado Set 29 por Thiago Trafane
Sergio, legal! Acho que ficou melhor assim.

Quanto à notação da distribuição, acredito que você esteja definindo os parâmetros de uma função, que é a função de densidade da Uniforme, não um intervalo. E mesmo que estivesse definindo um intervalo, o fato do numpy usar um intervalo fechado no limite inferior e aberto no superior não mudaria a função densidade da Uniforme, pois trata-se de uma variável aleatória contínua e, assim, a probabilidade de qualquer elemento do suporte é igual a 0. Desse modo, podemos,  sem perda de generalidade, definir qualquer intervalo tanto com limites fechados quanto com limites abertos, dado que, para quaisquer \(a\) e \(b\) tais que \(a \lt b\),

\(\mathbb{P}(a \lt x \lt b) = \mathbb{P}(a \leq x \lt b) = \mathbb{P}(a \lt x \leq b) = \mathbb{P}(a \leq x \leq b)\)
...