Search code examples
pythonnumpyiterationbinomial-coefficients

Iterative Binomial Update without Loop


Can this be done without a loop?

import numpy as np
n = 10
x = np.random.random(n+1)
a, b = 0.45, 0.55

for i in range(n):
    x = a*x[:-1] + b*x[1:]

I came across this setup in another question. There it was a covered by a little obscure nomenclature. I guess it is related to Binomial options pricing model but don't quite understand the topic to be honest. I just was intrigued by the formula and this iterative update / shrinking of x and wondered if it can be done without a loop. But I can not wrap my head around it and I am not sure if this is even possible.

What makes me think that it might work is that this vatiaton

n = 10
a, b = 0.301201, 0.59692
x0 = 123
x = x0
for i in range(n):
    x = a*x + b*x
# ~42

is actually just x0*(a + b)**n

print(np.allclose(x, x0*(a + b)**n))
# True

Solution

  • You are calculating:

    sum( a ** (n - i) * b ** i * x[i] * choose(n, i) for 0 <= i <= n)
    

    [That's meant to be pseudocode, not Python.] I'm not sure of the best way to convert that into Numpy.

    choose(n, i) is n!/ (i! (n-i)!), not the numpy choose function.


    Using @mathfux's comment, one can do

    import numpy as np
    from scipy.stats import binom
    
    binomial = binom(p=p, n=n)
    pmf = binomial(np.arange(n+1))
    res = np.sum(x * pmf)
    

    So

    res = x.copy()
    for i in range(n):
        res = p*res[1:] + (p-1)*res[:-1]
    

    is just the expected value of a binomial distributed random variable x.