Desired probability distribution

Question

I have a data "zcosmo" that follows a certain distribution. It has many fluctuations (i.e it increases slightly at some point and then flattens etc, but overall it increases from 0 to 0.5). So I fit the data using a function that follows the distribution of the data.

This function is the probability distribution lets say of the data. Here in this image, the fit gives me the blue line, which I want to use as the probabilty distribution function.

def fit(x,p1,p2):
    return (p1*x)+(p2*(x**2))

Now I theoretically know how to create data that follows this function!

• 1) compute the selection function in redshift from the distribution of cluster redshifts (the fitted function) and normalize it to unity to have a probability distribution, say p(y)

• 2) compute the indefinite integral of the probability distribution F(y)= sum (from 0 to y) p(y)dy

• 3) throw a random number uniform deviate x between 0 and 1

• 4) compute y = F-1(x)

"y" will then have the desired probability distribution p(y)

I am struggling to find a way to do it in Python. Is there a package that does the above??

Answer 1

You should take a look at the SciPy stack. It includes great libraries for probability and statistics, numerical integration and plotting.