solving a equation with integral

3 Ansichten (letzte 30 Tage)
Louis Liu
Louis Liu am 22 Aug. 2017
Bearbeitet: Louis Liu am 29 Aug. 2017
Hello,
I write some code like this,
>> C = 1.00;
>> n = 38;
>> alpha = 0.05;
>> Cp = C+0.33;
>> fun=@(y)chi2cdf((n-1)*(3*Cp*sqrt(n)-y).^2/(9*n*czero^2),n-1).*(normpdf(y+3*(Cp-C)*sqrt(n))+normpdf(y-3*(Cp-C)*sqrt(n)));
>> solve(integral(fun,0,3*Cp*sqrt(n))-alpha,czero)
and get the main error message : Undefined function or variable 'czero'.
Does anyone can tell me how to solve this equation? Using symbolic variable or other method?
Thanks!

Akzeptierte Antwort

Walter Roberson
Walter Roberson am 22 Aug. 2017
What you would like to do is roughly
C = 1.00;
n = 38;
alpha = 0.05;
Cp = C+0.33;
syms y czero real
fun = chi2cdf((n-1)*(3*Cp*sqrt(n)-y).^2/(9*n*czero^2),n-1).*(normpdf(y+3*(Cp-C)*sqrt(n))+normpdf(y-3*(Cp-C)*sqrt(n)));
solve( int(fun, y, 0,3*Cp*sqrt(n))-alpha, czero)
Unfortunately, chi2cdf does not expect symbolic values, and naively does a test on whether the input < 0 (because input values below 0 give an output of 0), but that test does not work on symbolic values. You can avoid the test by coding chi2pdf in terms of gampdf() but that only postpones the problem, as gampdf() has the same difficulty.
So, instead you have to work numerically:
C = 1.00;
n = 38;
alpha = 0.05;
Cp = C+0.33;
fun = @(y, czero) chi2cdf((n-1)*(3*Cp*sqrt(n)-y).^2/(9*n*czero^2),n-1).*(normpdf(y+3*(Cp-C)*sqrt(n))+normpdf(y-3*(Cp-C)*sqrt(n)));
czero_guess = 0;
czero_sol = fzero( @(czero) integral( @(y) fun(y, czero), 0, 3*Cp*sqrt(n) ) - alpha, czero_guess )
The output I am getting is -1.26021065783808
  3 Kommentare
Walter Roberson
Walter Roberson am 22 Aug. 2017
I had to choose something for the guess, and I had no idea of the value that would be found.
The only place your equation has czero is in the form czero^2 so your equation cannot tell the difference between negative and positive czero
Louis Liu
Louis Liu am 23 Aug. 2017
Oh, you are right! Thanks!

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Kategorien

Mehr zu MATLAB finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!