Есть ли более эффективный способ перечислить вероятность для каждого из возможных результатов дискретной случайной величины в Python или R? - PullRequest
Новая
       12

Есть ли более эффективный способ перечислить вероятность для каждого из возможных результатов дискретной случайной величины в Python или R?

1 голос
/ 13 марта 2019

Я вычисляю pmf теоретически на Python.вот код. 

>>> a_coin = np.array([0,1])
>>> three_coins = np.array(np.meshgrid(a_coin,a_coin,a_coin)).T.reshape(-1,3)
>>> heads = np.sum(three_coins, axis = 1)
>>> df = pd.DataFrame({'heads': heads, 'prob': 1/8})
>>> np.array(df.groupby('heads').sum()['prob'])
array([0.125, 0.375, 0.375, 0.125])

этот фрагмент кода имитирует 1 бросок 3-х монет.возможные результаты {0,1,2,3}.Последняя строка кода вычисляет вероятность для каждого из возможных результатов соответственно.

Я должен поместить 10 'a_coin' в np.meshgrid (a_coin, ..., a_coin), если я хочу вычислить pmf дляподбрасывание 10 честных монет, что кажется скучным и неэффективным.

вопрос в том, есть ли более эффективный способ сделать это в Python или R?

Ответы [ 3 ]

1 голос
/ 13 марта 2019

Вот решение на основе fft numpy: 

import numpy as np
from scipy import fftpack

def toss(n=10, p=0.5):
    t1 = np.zeros(fftpack.next_fast_len(n+1))
    t1[:2] = 1-p, p
    f1 = fftpack.rfft(t1)
    c1 = f1[1:(len(t1) - 1) // 2 * 2 + 1].view(f'c{2*t1.itemsize}')
    c1 **= n
    f1[::(len(t1) + 1) // 2 * 2 - 1] **= n
    return fftpack.irfft(f1)[:n+1]

Например: 

>>> toss(3)
array([0.125, 0.375, 0.375, 0.125])
>>> toss(10)
array([0.00097656, 0.00976562, 0.04394531, 0.1171875 , 0.20507813,
       0.24609375, 0.20507813, 0.1171875 , 0.04394531, 0.00976562,
       0.00097656])
1 голос
/ 13 марта 2019

Вот как это сделать в R: 

> sapply(0:3, choose, n=3)/sum(sapply(0:3, choose, n=3))
[1] 0.125 0.375 0.375 0.125

Функция choose дает вам биномиальные коэффициенты.Чтобы превратить их в вероятности, просто поделите их на суммы: 

sapply(0:10, choose, n=10)
 [1]   1  10  45 120 210 252 210 120  45  10   1

sapply(0:10, choose, n=10)/ sum( sapply(0:10, choose, n=10))
 [1] 0.0009765625 0.0097656250 0.0439453125 0.1171875000 0.2050781250 0.2460937500 0.2050781250
 [8] 0.1171875000 0.0439453125 0.0097656250 0.0009765625

Не похоже, что вы действительно хотите перечислить столько, сколько рассчитать.Если вам нужно перечислить результаты 10 последовательных "честных" биномиальных розыгрышей, то вы можете использовать combn 11 раз.

0 голосов
/ 13 марта 2019

Используя стандартные библиотеки Python, вы можете получить вероятности в виде рациональных чисел (это точное решение), например, 

from fractions import Fraction
from math import factorial

n=30
[Fraction(factorial(n), factorial(n - j)) * Fraction(1, factorial(j) * 2 ** n) for j in range(0, n + 1)]

Это можно легко преобразовать в числа с плавающей запятой, например, 

list(map(float, [Fraction(factorial(n), factorial(n - j)) * Fraction(1, factorial(j) * 2 ** n) for j in range(0, n + 1)]))
...