error in symbolic convolution using integral and isAlways

I want to symbolically integrate the following function using int

This is actually similar to convoution

so this is what I did

syms z alpha lam_f u_f sigma ep 
nCk = @(n,kVec,z)arrayfun(@(k)lam_f*exp(-lam_f*alpha).*(nchoosek(n,k)*(1*ep).^(k).*((1-ep)*normcdf(z,u_f,sigma)).^(n-(k)).*(normcdf(alpha-u_f-z,sigma))^k),kVec);
%expo 

%tf = isAlways(lam_f > 0 & u_f > 0  & sigma > 0)
F(z,u_f,lam_f,ep)=int(sum(nCk(n,0:n,z)),alpha,0, z,'IgnoreAnalyticConstraints',true);

However, I get the following error

Unable to prove 'z < u_f & sigma == 0' literally. Use 'isAlways' to test the statement mathematically.

I tried also using isAlways

[Torsten]

normcdf(alpha-u_f-z,sigma)

Why does this call to normcdf only has two arguments ? Especially: Why is the second argument the standard deviation ?

What's the problem from probability theory behind your question ?

[Muna Tageldin]

Yes, I made a mistake in formulating the question, it should be

normcdf(alpha-z,u_f,sigma)

Also, the aim behind this is to find the CDF function

[Torsten]

What is the underlying random variable for which your expression is the CDF ?

I'm almost certain that this complicated function does not have an analytic antiderivative.

[Muna Tageldin]

This is the random variable I'm trying to integrate where is the rest of summation part as in above

[Muna Tageldin]

Do you know if I can numerically solve this CDF because I need to find the mean/average as well

Answer 1

[Walter Roberson]

normpdf() is not designed to accept symbolic inputs.

You will need to re-implement normpdf without the line

p(sigma==0 & x<mu) = 0;

[Torsten]

It should be possible to redefine "normcdf" with the help of the error function (which accepts symbolic inputs).

[Walter Roberson]

The code hints

y = exp(-0.5 * ((x - mu)./sigma).^2) ./ (sqrt(2*pi) .* sigma)

In practice for symbolic work you would want

y = exp(sym(-0.5) * ((x - mu)./sigma).^2) ./ (sqrt(2*sym(pi)) .* sigma)

[Torsten]

OP works with the CDF, not the PDF.

[Walter Roberson]

I just noticed ;-(

z = (x-mu) ./ sigma;
p = 0.5 * erfc(-z ./ sqrt(2));

appears to be the basic code.

[Muna Tageldin]

Thank you, I tried using erf function and it works with symbolic integration

[Torsten]

@Muna Tageldin

Can you show us how you succeeded ? I don't get it (e.g. for k = 2):

syms x 
syms alpha lambda positive
syms mu sigma positive
int((0.5*erfc(-(alpha-x-mu)/sigma/sqrt(sym('2'))))^2*lambda*exp(-lambda*alpha),alpha,0,x)

ans = 

[Muna Tageldin]

n=3;
syms z positive 
syms alpha positive 
syms lam_f positive 
syms u_f positive 
syms sigma positive
syms ep positive
syms z positive

nCk = @(n,kVec,z)arrayfun(@(k)lam_f*exp(-lam_f*alpha).*(nchoosek(n,k)*(1*ep).^(k).*((1-ep)*(0.5+0.5*erf((z-u_f)/(sqrt(2)*sigma)))).^(n-(k)).*(0.5+0.5*erf((z-alpha-u_f)/(sqrt(2)*sigma)))^k),kVec);

F(z,u_f,lam_f,sigma,ep)=int(nCk(n,0:n,z),alpha,0, z)

F(z, u_f, lam_f, sigma, ep) = 

[Torsten]

If you look at what comes out for F(z,u_f,lam_f,sigma,ep), you'll see that the erf integrals are still left.

No chance for an analytical antiderivative, I guess.

Answer 2

[Steven Lord]

Are you trying to tell MATLAB that % lam_f > 0 & u_f > 0 & sigma > 0 or are you trying to ask MATLAB if those three conditions are satisfied? Using isAlways asks that question. To tell MATLAB use assume and assumeAlso or tell MATLAB when it creates those symbolic variable that they're all positive, passing the positive option to syms or sym: syms x positive.

[Muna Tageldin]

I'm trying to tell MATLAB these conditions are true. So using assume and assume Also is true. However, I get the same error when I enforce these assumptions (lam_f > 0 & u_f > 0 & sigma >). I'm not sure what other methods I can use?

[Steven Lord]

However, I get the same error when I enforce these assumptions

Can you show your updated code and the full and exact text of the error message? I'd be a bit surprised if you told MATLAB that sigma was always positive and it couldn't tell if that condition was satisfied or not (it wouldn't be.)

Please show us all the text displayed in red in the Command Window (and if there are any warning messages displayed in orange, please show us those too.) The exact text may be useful and/or necessary to determine what's going on and how to avoid the warning and/or error.

[Muna Tageldin]

n=2;
syms z positive 
syms alpha positive 
syms lam_f positive 
syms u_f positive 
syms sigma positive
syms ep positive
syms z positive
%assume(z < u_f )
assume (ep <1 & ep>0)
assumeAlso(sigma > 0)
%assume(x/2,"integer")
nCk = @(n,kVec,z)arrayfun(@(k)lam_f*exp(-lam_f*alpha).*(nchoosek(n,k)*(1*ep).^(k).*((1-ep)*normcdf(z,u_f,sigma)).^(n-(k)).*(normcdf(alpha-u_f-z,sigma))^k),kVec); 

Error using symengine

Unable to prove 'z < u_f & sigma == 0' literally. Use 'isAlways' to test the statement mathematically.

[Muna Tageldin]

The error I got :

Error using symengine

Unable to prove 'z < u_f & sigma == 0' literally. Use 'isAlways' to test the statement mathematically.

Error in sym/subsindex (line 849)

X = find(mupadmex('symobj::logical',A.s,9)) - 1;

Error in sym/privsubsasgn (line 1124)

L_tilde2 = builtin('subsasgn',L_tilde,struct('type','()','subs',{varargin}),R_tilde);

Error in sym/subsasgn (line 961)

C = privsubsasgn(L,R,inds{:});

Error in normcdf>localnormcdf (line 109)

p(sigma==0 & x<mu) = 0;

Error in normcdf (line 50)

[varargout{1:max(1,nargout)}] = localnormcdf(uflag,x,varargin{:});

Error in testing_convolution (line 75)

nCk =

@(n,kVec,z)arrayfun(@(k)lam_f*exp(-lam_f*alpha).*(nchoosek(n,k)*(1*ep).^(k).*((1-ep)*normcdf(z,u_f,sigma)).^(n-(k)).*(normcdf(alpha-u_f-z,sigma))^k),kVec);

Error in testing_convolution (line 75)

nCk =

@(n,kVec,z)arrayfun(@(k)lam_f*exp(-lam_f*alpha).*(nchoosek(n,k)*(1*ep).^(k).*((1-ep)*normcdf(z,u_f,sigma)).^(n-(k)).*(normcdf(alpha-u_f-z,sigma))^k),kVec);

Error in testing_convolution (line 81)

F(z,u_f,lam_f,ep)=int(sum(nCk(n,0:n,z)),alpha,0, z);%'IgnoreAnalyticConstraints');