Services & ResourcesWolfram ForumsMathGroup Archive

[Date Index] [Thread Index] [Author Index]

Re: Integration Program of Burr Distribution

• To: mathgroup at smc.vnet.net
• Subject: [mg91565] Re: [mg91545] Integration Program of Burr Distribution
• From: Daniel Lichtblau <danl at wolfram.com>
• Date: Wed, 27 Aug 2008 06:44:09 -0400 (EDT)
• References: <200808260730.DAA23794@smc.vnet.net>

Dewi Anggraini wrote:
> Dear All
> 
> I'd like to calculate the upper specification limit of Burr Distribution from my data below:
> 
> n = 69;
> x1 = {3, 5, 6, 14, 11, 7, 7, 10, 11, 5, 4, 19, 9, 2, 8, 5, 6, 6, 5,
>    5,
>        4, 5, 5, 6, 8, 5, 8, 6, 5, 16, 5, 18, 16, 21, 6, 3, 16, 8, 3,
>    8,
>        11, 2, 3, 4, 8, 7, 9, 10, 8, 11, 8, 10, 9, 12, 12, 9, 6, 12,
>    3, 9,
>        14, 7, 4, 13, 8, 14, 5, 8, 2};
> 
> with the initialisation pdf of Burr as follows:
> 
> BurrDistribution[x1_, c_,
>   k_] := (c*k)*(x1^(c - 1)/(1 + x1^c)^(k + 1))
> pdf = BurrDistribution[x1, c, k]
> 
> I already got the parameter value of c and k through MLE:
> 
> logl = Plus @@ Log[pdf]
> maxlogl = FindMinimum[{-logl, c > 0 && k > 0}, {c, 1}, {k, 2},
>      MaxIterations -> 1000]
> mle = maxlogl[[2]]
> 
> However, when I wanna calculate the integration of Burr with the value of c=37.8115, k=0.0135614 and x1 limit = [8, Infinity] to get the probability of data fall outside the upper specification limit given {8} with this program below:
> 
> In[2]:= Integrate[pdf[37.81151579009424, 0.01356141249769735, t], 
>    {t, 8, Infinity}]
> 
> Out[2]= Integrate[pdf[37.81151579009424, 0.01356141249769735, t], 
>    {t, 8, Infinity}]
> 
> That gives me the repetation not a real number.
> 
> Why is this happen? Do I make mistake somewhere? Because when I did the similar program towards Gamma Distribution, it works well.
 >
> I look forward for your reply and assistance. I highly appreciate the suggestions you contribute to my minor thesis.
> 
> Thank You.
> Regards,
> Dewi

Among the obvious problems:

(1) You defined pdf to be a (very large) list, not a function of three 
parameters.

(2) Your integration occurs at In[2]. This means you could not have 
(re)defined all the necessary information for pdf. I would surmise you 
used Quit[] or otherwise restarted your kernel, and did not redefine things.

More subtle:

(3) I doubt your MLE is what you intend. You minimized a sum of logs of 
functions of data values (ordinates, or y values). They are in no way 
related to the abscissae, or x values. These, I presume, are meant to be 
values 1,2,...,69. Also, you will require a sum of discrepancie measures 
(absolute values, or maybe squares). So I think your MLE optimization 
should involve an expression such as logtot below.

pdflist = BurrDistribution[Range[n], c, k];
logtot = Total[Abs[PowerExpand[Log[pdflist]] - Log[x1]]];

Or possibly it could just use
tot = Total[Abs[pdflist - x1]];

And you could get the optimal values via

minlogl =
   FindMinimum[{logtot, c > 0 && k > 0}, {c, 1}, {k, 2},
    MaxIterations -> 1000];
mle = minlogl[[2]]

If you now define

pdf[t_, c_, k_] := BurrDistribution[t, c, k]

(note order of arguments) then you can do quadrature as below. Observe 
that one can use mle to substitute numeric values for the parameters; 
you do not (and should not) cut/paste from the earlier computation.

NIntegrate[pdf[t, c, k] /. mle, {t, 8, Infinity}]


Daniel Lichtblau
Wolfram Research

• Follow-Ups:
  • Re: Re: Integration Program of Burr Distribution
    • From: Darren Glosemeyer <darreng@wolfram.com>

• References:
  • Integration Program of Burr Distribution
    • From: "Dewi Anggraini" <dewi_anggraini@student.rmit.edu.au>