[Date Index]
[Thread Index]
[Author Index]
Re: nonLinearFit and nonLinearRegress
*To*: mathgroup at smc.vnet.net
*Subject*: [mg55319] Re: [mg55298] nonLinearFit and nonLinearRegress
*From*: Daniel Lichtblau <danl at wolfram.com>
*Date*: Sat, 19 Mar 2005 04:45:21 -0500 (EST)
*References*: <200503181034.FAA14747@smc.vnet.net>
*Sender*: owner-wri-mathgroup at wolfram.com
John Hawkin wrote:
> I'm trying to do a nonLinear fit on some data, fitting a parameter
> which appears in the limits of the integral that I'm trying to fit.
> When I do a nonLinearFit, I get about a page of errors, but it does
> work. It takes about 10 minutes to do it for 20 points with a
> reasonable initial estimate. However when I use the nonLinearRegress
> command, even with the best fit parameter from nonLinearFit given as
> the initial guess, it takes hours to work. Is it possible there is a
> way to speed this up, and is this normal?
>
> The commands that I'm using to integrate this function are:
>
> h[B]:=2/(B pi) NIntegrate[x Sin[x] Exp[-(x/B)^3/2], {x, 0, inf},
> Method-Oscillatory];
> holts[c_, d_, Q_]:= NIntegrate[h[B], {B, c/Q, d/Q}];
>
> where holts is the function I am fitting (the h[B] integral is known
> as the Holtsmark integral). The fitting command I'm using is:
>
> NonlinearFit[fitData, holts[c,d,Q], {c,d}, {Q, estimateQ}];
>
> where fitData contains sets of size 3, with the values of c and d as
> the first two values, and the data point that I'm fitting the integral
> to as the 3rd value. If anyone has any ideas on how I can perform the
> nonlinear regression fast enough so that I can do it many times (maybe
> 50), I would greatly appreciate it. Thanks very much,
>
> John Hawkin
One possibility for speeding it would be to use the analytic form for h.
n[9]:= InputForm[i1 = Integrate[x*Sin[x]* Exp[-(x/B)^3/2], {x, 0,
Infinity}, Assumptions->B>0]]
Out[9]//InputForm=
(2*B^3*((-1)^(1/3)*Sqrt[3]*Pi*
((-1)^(1/3)*BesselI[-2/3, ((2/3 + (2*I)/3)*B^(3/2))/Sqrt[3]] -
BesselI[2/3, ((2/3 + (2*I)/3)*B^(3/2))/Sqrt[3]] -
(-1)^(1/3)*BesselJ[-2/3, ((2/3 + (2*I)/3)*B^(3/2))/Sqrt[3]] +
BesselJ[2/3, ((2/3 + (2*I)/3)*B^(3/2))/Sqrt[3]]) +
9*HypergeometricPFQ[{1}, {1/3, 2/3, 5/6, 7/6}, -B^6/2916]))/27
This can be evaluated at machine precision for B not very large. If your
upper bound d/Q is large for plausible values of Q you might need to
increase WorkingPrecision in the outer integration holts[...].
I don't know if this approach will help but it might be worth trying.
Obviously it relies on the fact that we can find an analytic form for
h[]. Taking it a step further, here is what I get for
Integrate[i1, {B, c/Q, d/Q}, Assumptions->{Q>0,0<c<d}]
(((-2 + 2*I)*(-1)^(2/3)*Pi*
(c^(5/2)*BesselI[-5/3, ((2/3 + (2*I)/3)*(c/Q)^(3/2))/Sqrt[3]] -
d^(5/2)*BesselI[-5/3, ((2/3 +
(2*I)/3)*(d/Q)^(3/2))/Sqrt[3]]))/Q^(5/2) +
(2*2^(2/3)*Pi*((-1)^(1/3)*c^(5/2)*BesselI[5/3,
((-2/3 + (2*I)/3)*(c/Q)^(3/2))/Sqrt[3]] -
d^(5/2)*BesselI[5/3, ((2/3 - (2*I)/3)*(d/Q)^(3/2))/Sqrt[3]]))/
((1 + I)^(1/3)*Q^(5/2)) - (2*(-1)^(11/12)*(1 + I)^(2/3)*2^(1/6)*Pi*
(c^(5/2)*BesselI[5/3, ((2/3 + (2*I)/3)*(c/Q)^(3/2))/Sqrt[3]] -
d^(5/2)*BesselI[5/3, ((2/3 +
(2*I)/3)*(d/Q)^(3/2))/Sqrt[3]]))/Q^(5/2) +
9*2^(2/3)*Sqrt[3]*Pi*(-Hypergeometric0F1Regularized[-2/3,
(((-2*I)/27)*c^3)/Q^3] + Hypergeometric0F1Regularized[-2/3,
(((-2*I)/27)*d^3)/Q^3]) +
(3*(-(c^4*HypergeometricPFQ[{1}, {1/3, 5/6, 7/6, 5/3},
-c^6/(2916*Q^6)]) +
d^4*HypergeometricPFQ[{1}, {1/3, 5/6, 7/6, 5/3},
-d^6/(2916*Q^6)]))/Q^4)/18
As for verification of correctness, you're on your own. Assuming it is
basically sound, I have no idea offhand whether this will suffer from
branch cut problems when evaluated at specific {d,c,Q}.
Another possibility for approaching a problem such as this might be to
use a series or similar such approximation for the first integrand,
integrated termwise. This will not rely on vagaries of NIntegrate error
control, but you'd need to control it yourself in choosing how many
terms in the approximation to use.
Daniel Lichtblau
Wolfram Research
Prev by Date:
**Re: Add new option to Notebook[]**
Next by Date:
**findroot**
Previous by thread:
**Re: nonLinearFit and nonLinearRegress**
Next by thread:
**Re: nonLinearFit and nonLinearRegress**
| |