I want to sample from the probability distribution function (PDF)
for $0<t<1$. The function looks like this:
The PDF can not be analytically inverted, so I tried to use Metropolis Monte Carlo sampling with the following python code:
pdf = lambda t: 3*np.exp(-t)*(0.5+0.5*np.exp(-8*t)*np.cos(40*t)) lw = 0 up = 1 sig = 0.1 ini_t = 0.01 new_t = ini_t samples =  for i in range(100000): cand_t = scipy.stats.truncnorm.rvs((lw-new_t)/sig, (up-new_t)/sig, loc=new_t, scale=sig) if np.random.uniform(0, 1) < pdf(cand_t)/pdf(new_t): new_t = cand_t samples.append(new_t) burn_in = 10000 samples = np.array(samples[burn_in:])
The result looks like the green plot of the following image:
As it is obvious, the performance of sampling is rather poor. How can I improve it?