sample scipy distributions using the same probability array

Question

I would like to generate variates from different Beta distributions using a set of user-generated probabilities. I understand that scipy.stats.beta.rvs(a,b) generates random variates using probabilities generated internally but I would like to use my own probabilities instead. I suspect this must be trivial but have not been able to work it out myself. Thanks in advance.

Answer 1

After poking around a bit, I think the solution is to use scipy.special.betainc(a, b, x). This seems to give the same CDF as the following script shows:

import scipy.stats
import scipy.special
import matplotlib.pyplot as plt
import numpy as np

a = 3.0
b = 5.0

ntrial = 100000

p = np.linspace(0.0, 1.0, 1000)
r = np.random.rand(ntrial)
v = scipy.special.betainc(a, b, r)

# Plot the CDFs

plt.plot(np.arange(ntrial) / float(ntrial), np.sort(v), '.', label='scipy.special')
plt.plot(p, scipy.stats.beta.cdf(p, a, b), label='scipy.stats')
plt.ylabel("cumul. prob.")
plt.legend()