Using MATLAB, how can I find the integral of a bounded CDF?

Question

Using MATLAB, I am trying to find the integral of a bounded range of a CDF. Please refer to the following code:

u = 1;
s = 1;
X = random('Normal',u,s,1,10000);
pd = makedist('Normal','mu',u,'sigma',s);
xAxis = min(X):.0001:max(X);
c_pd = cdf(pd,xAxis);
r = icdf(pd,[.3,.6]);
plot(xAxis,c_pd)

Basically, I am trying to integrate c_pd between the corresponding X values for .3 and .6 (found by using icdf). However, c_pd is a vector and not the actual cdf function. Does anyone have ideas on how I can find the integral of this regardless of the distribution type (i.e. Normal, Rician, etc.)? Please advise. Thank you.

Answer 1

First, it looks like you're off on the wrong step using random. Second, you're correct in that the documentation for Matlab's ProbDist classes is not particularly good and lacks examples.

Let's use the these parameters that match those that you used (however I'm not sure if you bounds, a and b, are meant to be to be probabilities, P(X), -do you want do integrate over a range of probabilities? -in that case you actually want to use the inverse CDF):

mu = 1; sig = 1;
a = 0.3; b = 0.6;

You have several options depending on what you need. It'll likely be easiest to use the integral function to perform numerical integration (quadrature). First you can implement the CDF yourself:

normalCDF = @(t,mu,sig)(1+erf((t-mu)./(sqrt(2)*sig)))/2;
q = integral(@(t)normalCDF(t,mu,sig),a,b)

Or use one of the older style cdf function, normcdf in this case:

q = integral(@(t)normcdf(t,mu,sig),a,b)

Or use the general cdf function (type help cdf to see a list of all of the supported distributions for the older style CDFs):

q = integral(@(t)cdf('norm',t,mu,sig),a,b)

Or use one of the new ProbDistUnivParam class methods:

normalPD = ProbDistUnivParam('normal',[mu sig]);
q = integral(@(t)normalPD.cdf(t),a,b)

See here for a list of distributions supported with this new class. Note that the .cdf(t) is not to be confused with the cdf function used just above. This one is a method of the ProbDistUnivParam class. Type help ProbDistUnivParam and help ProbDistUnivParam/cdf.

If you want to attempt to solve for symbolic solutions, then you'll likely need to implement CDF functions yourself. Most high level Matlab functions support floating point calculations only unless they are part of the Symbolic Toolbox. Here's how you might go about solving for for these symbolically using int:

syms t MU SIG A B real
normalCDFsym = (1+erf((t-MU)./(sqrt(2)*SIG)))/2;
qsym = simplify(int(normalCDFsym,t,A,B)); % Solve integral symbolically
pretty(qsym)                              % Print out result
q = subs(qsym,{MU,SIG,A,B},{mu,sig,a,b})  % Plug in numeric values

Note that for more complicated distributions, you may not always be able to get a solution. Also, here I've left MU, SIG, A, and B all as symbolic. In some cases you may not be able to get a solution with all symbolic parameter so you can try letting some of them be explicit values if you know what those values are, e.g.:

syms t MU A real
normalCDFsym = (1+erf((t-MU)./(sqrt(2)*sym(1))))/2;
qsym = simplify(int(normalCDFsym,t,A,sym(0.6))); % Solve integral symbolically
pretty(qsym)                                     % Print out result
q = subs(qsym,{MU,A},{mu,a})                     % Plug in numeric values