2
$\begingroup$

I have the following simple task at hand: given the data

data={{0., 1.49021*10^-6}, {0.0104212, 0.}, {0.011298963693097342, 
  0.}, {0.012250704096577147, 0.0000300991}, {0.0132826120110092413, 
  0.}, {0.014401440149411573, 0.}, {0.0156145100229668453, 
  0.000047227}, {0.0169297598523363013, 
  0.000021779}, {0.018355796514665096, 
  0.000020087}, {0.019901951866215888, 0.}, {0.0215783, 
  0.}, {0.023395942502427995, 0.}, {0.025366642149038433, 
  0.0000145421}, {0.0275033388311069493, 
  0.0000804508}, {0.029820014900448666, 
  0.0000865697}, {0.0323318304778776877, 0.0000798299}, {0.0350552, 
  0.00011572}, {0.038008, 0.0000970346}, {0.0412095, 
  0.000107386}, {0.0446807074321427938, 
  0.0000825284}, {0.048444276877472492, 
  0.000182684}, {0.052524861334064740, 
  0.000133403}, {0.056949163781318632, 
  0.000239611}, {0.061746136458395343, 
  0.000244875}, {0.066947170325079339, 0.000275425}, {0.0725863, 
  0.000396276}, {0.0787004, 0.000548252}, {0.085329565860813736, 
  0.000886005}, {0.092517091871782364, 
  0.00110423}, {0.100310041450036191, 
  0.0014339}, {0.10875941096000785, 0.00176335}, {0.11792049231940185,
   0.00207987}, {0.127853234824554518, 
  0.00252122850098223242}, {0.13862263745326262, 
  0.0030170796711709233}, {0.150299174212275455, 
  0.0033715876826087561}, {0.16295925531288655, 
  0.0041892010005696428}, {0.176685727192516073, 
  0.00414758}, {0.191568, 
  0.00432007698126049211}, {0.20770470867296703, 
  0.00490949019252781349}, {0.2252, 
  0.00510288435189528265}, {0.244169388724396586, 
  0.0051927526915368576}, {0.26473639635925972, 
  0.0052866}, {0.28703581527328598, 0.0054873353922293725}, {0.311214,
   0.0055598252377202850}, {0.33742788021063624, 
  0.00580761254792859567}, {0.36585028739592496, 
  0.00602367890695240577}, {0.396667, 
  0.00671594187572154046}, {0.43007903676732617, 
  0.00722826511875487534}, {0.46630568613971068, 
  0.0082064615556007357}, {0.505584, 0.0090722}, {0.54817, 
  0.0100235944466727454}, {0.594344190553855456, 
  0.0111068159544501042}, {0.644407310950390055, 
  0.01220849026653985611}, {0.69868737510387846, 
  0.0136699881071071293}, {0.757539586895110695, 
  0.0154469611857013811}, {0.821349069929730446, 
  0.0172370134135830941}, {0.89053338775263302, 
  0.01902048386724843813}, {0.965545, 
  0.0209349937333521787}, {1.04687560679516989, 
  0.0231107424931955134}, {1.13506, 
  0.0249087477731436135}, {1.2306652969727638, 
  0.027179116847683504}, {1.3343273599816916, 
  0.028790797370158439}, {1.44672114178833011, 
  0.0308187}, {1.5685821372395816, 0.0330423}, {1.7007077937809543, 
  0.035765810727818562}, {1.84396, 0.0385458646152422482}, {1.99928, 
  0.0418738}, {2.16768919425212747, 
  0.0453100275928077234}, {2.35027915952626021, 
  0.0490147988668498693}, {2.5482491412285828, 
  0.052714906696813697}, {2.7628946371975260, 
  0.056930486441867330}, {2.9956202683441231, 
  0.060554032423127421}, {3.2479489703655182, 
  0.0655495947236072779}, {3.5215319596998365, 
  0.0714955617488584674}, {3.81815953893874127, 
  0.0777785}, {4.1397728124073332, 
  0.084507518747138657}, {4.4884763885770909, 
  0.091266254006932421}, {4.86655, 0.097710114642383951}, {5.27647, 
  0.103619629618812047}, {5.72092501825008881, 
  0.108076130674904863}, {6.2028130597590074, 
  0.111359951413410502}, {6.7252917546689925, 
  0.1138569293950436445}, {7.2917801567884117, 
  0.111399689696608317}, {7.9059853155099589, 
  0.106300879052989579}, {8.57192653441002861, 
  0.0982251}, {9.2939616732115304, 
  0.087549107149667511}, {10.076815665227805, 
  0.0735944}, {10.925611436903257, 0.0592198}, {11.845903431785432, 
  0.045930099273712776}, {12.8437139583063598, 
  0.032735604064644489}, {13.925572599229811, 
  0.022157895599934801}, {15.0985589406564191, 
  0.0141406780873452052}, {16.370348900201357, 
  0.0082515394244985029}, {17.749264957512075, 
  0.0046942188705470563}, {19.244330615830133, 
  0.0025202677958886972}, {20.8653294509886322, 
  0.001309568579533558067}, {22.622869134256639, 
  0.00067969169862283607}, {24.5284508479904275, 
  0.000329474452633988493}, {26.5945, 
  0.000155956667282389328}, {28.8346705675320649, 
  0.000072887856658869226}};

I'd like to find the best fit of the form

$$ f (x) = \left( \frac{x}{a} \right)^{\alpha} e^{- (x/b)^{\beta}} $$ where $x$ is the variable and $a, b, \alpha, \beta$ are fit parameters.

So I run

fit = FindFit[
  data, {(x/a)^α E^(-(x/b)^β), α > 0.01, β > 
    0.01, a > 0, a < 100, b > 0, 
   b < 100}, {{a, 1}, {b, 2}, {α, 1}, {β, 3}}, x]

And all I get is a bunch of errors, like

Indeterminate expression 0.^0. encountered

and

The function value {Indeterminate,6.5666,Indeterminate,Indeterminate} is not
a list of real numbers with dimensions {4} at {a,b,α,β} = {1.`,2.`,1.`,3.`}

At first I thought I messed up the order within the FindFit function call, but I set it up like the help says: FindFit[data, model, parameters, variables]

I even removed the starting values from parameter list and constraints on parameters from model, but the errors are the same.

Even though the FindFit returns some parameter values, it's clearly bogus, so there must be some caveat in the minimalization routine that somehow encounters indeterminate numbers.

I find the 0^0 indet error pretty strange, as I explicitly limited both $\alpha$ and $\beta$ to be larger than 0.01, so there is no way to get anything of the form 0^0 anywhere in the process.

Can anyone shed some light on this issue please?

$\endgroup$
4
$\begingroup$

Getting rid of a degenerate data point and good starting values are what you need.

The first data with $x=0$ is the troublemaker that gives you the

0^0 error message

error message.

Using @Bill 's starting values and dropping the first data point:

fit = NonlinearModelFit[data[[2 ;; Length[data]]], (x/a)^α E^(-(x/b)^β),
  {{a, 35}, {b, 7}, {α, 1}, {β, 3}}, x];
fit["BestFitParameters"]
(* {a -> 35.0462, b -> 9.45751, α -> 1.1014, β -> 2.76364} *)

Show[ListPlot[data, PlotRange -> All],
  Plot[fit[x], {x, 0, 30}]]

data and fit

Data points 2, 3, 5, 6, 10, 11, and 12 are also problematic in that the response variable is zero when it is likely that those data points are censored in that the actual measurement is something below a measurement threshold. In such cases Tobit regression might be appropriate.

| improve this answer | |
$\endgroup$
  • $\begingroup$ Hi and thank you for your answer. I am curious as to why the point x = 0 causes trouble. I don't see why would this point be any special: since I limit the powers alpha and beta to be strictly larger than zero (> 0.01), the situation where expression 0^0 occurs, shouldn't even come up. So why? $\endgroup$ – user16320 Jan 2 at 19:23
2
$\begingroup$

Try this

fit = FindFit[data,{Abs[(x/a)^α] E^(-Abs[(x/b)^β]) ,α>0.5, β>0.5, 0<a<100, 0<b<100},
  {{a,35},{b,7},{α,1},{β,3}},x]

and if you get about the answer I get then use the fit values from that for a plot.

Show[Plot[(x/a)^α E^(-(x/b)^β)/.{a->35.0475,b->9.4576,α->1.10139,β->2.76371}, {x,0,28}],
  ListPlot[data]]

I would like a better away to avoid complex values coming from those fractional powers, but that was the best I could do at the moment.

| improve this answer | |
$\endgroup$
  • $\begingroup$ But how does MMA even arrive at complex values coming from fractional power? The input x is nonnegative (smallest value of x from the data is 0), I manually capped all parameters to be positive reals, so how does MMA arrive at any expression that evaluates to a complex number? :O $\endgroup$ – user16320 Jan 2 at 12:12
0
$\begingroup$

(Comment not an answer...)

Is it of interest the fit of multiple functions of the required form?

Here is a plot comparing the fit from Bill's answer and the this function:

 0.0028877 E^(-(x^3/1331)) x^0.8 + 
 0.00697135 E^(-(x^3/1728)) x^0.8 + 
 0.000302041 E^(-(x^3/5832)) x^0.8 + 
 0.00044096 E^(-0.00474737 x^2.5) x^1.1 + 
 0.00603381 E^(-0.00474737 x^2.5) x^1.2 + 
 0.00331148 E^(-0.000558514 x^3.5) x^1.5

enter image description here

| improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.