How to generate random samples for a custom probability density function

I am using 'fmincon' and 'patternsearch' of the OptimizationToolbox, R2016a. My ObjectiveFunction involves estimating the mean and variance of many arbitrary pdfs that are 1) not necessarily integrating to 1, and 2) generated in realtime by other custom functions I wrote so that the pdfs are numerically approximated rather than having any explicit analytical expression. In other words, each pdf is merely expressed as a vector 'likelihood' as a function of x, and the shape/expression of the pdf is unknown and changes as other parameters in the ObjectiveFunction vary across iterations. So I guess I should generate random samples for each pdf and estimate the stats based on these samples. However, I need the sampling to be done as fast as possible because I might need >10000 samples each time.

Answer 1

[Alan Weiss]

I think that the first thing you need to do is scale your likelihood so that it becomes a real pdf that integrates to 1. You can call integral for that purpose. Then I suggest that, instead of taking random samples, you numerically integrate moments of your newly-made pdf to find the mean and variance.

Alan Weiss

MATLAB mathematical toolbox documentation

[Yinan Cao]

Each iteration generates 2 pdfs for 2 independent choices, respectively, e.g., A or B, the integrals of pdf_A and pdf_B sum to 1. integral(pdf_A) = 0.6 and integral(pdf_B) = 0.4, for instance. Someone told me I don't have to generate samples. Instead, I can simply approximate the mean and variance as follows: (assuming x is the vector of the random variable and L consists of the likelihood values of x)

MEAN_A = sum(x.*(L_A./sum(L_A))); VAR_A = sum(x.^2.*(L_A/sum(L_A))) - MEAN_A^2;

MEAN_B = sum(x.*(L_B./sum(L_B))); VAR_B = sum(x.^2.*(L_B/sum(L_B))) - MEAN_B^2;

Just to clarify: If you plot(x, L_A) you get the probability density function for choosing A; And if you plot(x, L_B) you get the probability density function for choosing B.