Visualising a sine random variable in python

Question

Given a probability distribution with probability density f(t) = 0.5sin(t) for t in [0,π], plot a graph of the probability density on the interval [0,π] and in the same figure make a histogram of 10^4 random independent draws. The answer should look like this:

So far I have figured out how to plot the pdf function by doing the following:

import numpy as np
import matplotlib.pyplot as plt
import math

# Define the probability distribution
def f(t):
    return 0.5 * np.sin(t)

# Plot the probability density
t = np.linspace(0, math.pi, 10000)
pdf = f(t)
plt.plot(t, pdf, label='density')

plt.xlabel('value')
plt.ylabel('density / relative frequency')
plt.title('Visualisation of sine random variable')
plt.legend()
plt.grid(axis='y')
plt.show()

This gives me the blue line that I should have, but I dont know how to plot the histogram

I had a similar problem but with a exp(1/2) distibution on the interval [0;5.5]. There I did it like this:

import numpy as np
import matplotlib.pyplot as plt
import math

# Define exponential distribution 
rate = 0.5
def pdf(x):
    return rate*math.exp(-rate*x) #PDF = λe^(-λx), with λ the rate of the distribution

# Plot pdf
x = np.linspace(0, 5.5, 10000) 
pdf = np.array([pdf(i) for i in x])
plt.plot(x, pdf, label='density')

# Create random variables
stochast = np.random.exponential(1/rate, size=10000) 

# Make the histogram
plt.hist(stochast, bins=np.arange(0, max(stochast)+1, 0.2), density=True,
label='relative frequency 10^4 draws', edgecolor='black')

plt.xlabel('value')
plt.ylabel('density / relative frequency')
plt.title('Visualisation of a Exp(1/2) random variable')
plt.legend()
plt.xlim(0, 5.5)
plt.grid(axis='y')
plt.show()

Here I knew the rate λ of this exponential distribution which, after looking at wikipedia for a while, made it quite easy. But I don't know how to do this for the f(t) = 0.5sin(t) distribution.

Does this answer your question? [How to generate random numbers with predefined probability distribution?](https://stackoverflow.com/questions/51050658/how-to-generate-random-numbers-with-predefined-probability-distribution) — Peter O., Mar 17 '23 at 13:11
See also: https://stackoverflow.com/questions/55724501/generate-a-random-number-according-to-a-specific-distribution-with-a-function/64961053#64961053 — Peter O., Mar 17 '23 at 13:13

score 0 · Answer 1 · answered Mar 16 '23 at 20:53

You already know plotting, so question is how to sample from your PDF. And the answer is Inverse Transform Sampling

So, having

PDF(x) = sin(x) / 2

Find CDF(x) integrating PDF

CDF(x) = S PDF(x) dx from 0 to x = 1/2 - 1/2 cos(x)

Now find inverse CDF function, use U01 (uniform in 0...1 numbers) and you have your sampling routine

x = arccos(1 - 2 u01)

Code for illustation

import numpy as np
import matplotlib.pyplot as plt

rng = np.random.default_rng()

def sample_sin(N):
    u01 = rng.random(N)
    q = np.arccos(1.0 - 2.0*u01)
    return q

x = np.linspace(0.0, np.pi, 101)
cdf = 0.5 - 0.5*np.cos(x)
pdf = 0.5*np.sin(x)

plt.plot(x, cdf, "g-")
plt.plot(x, pdf, "r-")

plt.hist(sample_sin(100000), bins = 20, density = True)
plt.title("1/2 sin(x)")
plt.show()

should produce something like

Visualising a sine random variable in python

1 Answers1