2
$\begingroup$

I have a set of data (42 points):

{{1.1*10^6, 1/21}, {600., 3/7}, {80., 5/7}, {100., 29/42}, {600., 3/
  7}, {600., 3/7}, {70., 31/42}, {500., 19/42}, {25., 37/42}, {25., 
  37/42}, {300., 4/7}, {25., 37/42}, {15000., 5/42}, {150., 2/
  3}, {1400., 13/42}, {10., 1}, {60., 16/21}, {200., 13/21}, {10000., 
  4/21}, {10000., 4/21}, {1000., 5/14}, {50., 17/21}, {1500., 2/
  7}, {200., 13/21}, {300., 4/7}, {5000., 3/14}, {10., 1}, {80000., 1/
  14}, {20., 13/14}, {1575., 11/42}, {400., 1/2}, {1000., 5/
  14}, {9.1671*10^6, 1/42}, {60000., 2/21}, {20., 13/14}, {50., 17/
  21}, {400., 1/2}, {150., 2/3}, {10., 1}, {10000., 4/21}, {350., 11/
  21}, {2200., 5/21}}

that's look like (Log-Log plot):

enter image description here where x-coordinate represent area (of forest fire, in m^2), and y-coordinate represent number of fires that are equal or greater than value on x-coordinate, divided by total numbers of fire (cumulative distribution).

"First part" of data fits very good with so-called "stretched exponential":

(y=exp -(x/A)^beta),

where A and beta are constants that I must find.

"Second part" of data fits better to Power law

(y=C*x^(-alpha)),

where C and alpha are also constants that I must find.

My question is: how to make Mathematica to find best fit for this data and for BOTH function at the same time? I want to get X*-value that is "border" between this two function, all constants that are present in functions and curve that is smooth but consist of two part: left part ("before" X*) should be a part of stretched exponential (curve in this graph) and right part ("after" X*) should be plot of Power law (straight line in LogLogPlot).

How to do all of that?

Thanks for any help!

Improve this question
$\endgroup$
11
  • 2
    $\begingroup$ You could try Piecewise. $\endgroup$
    – b.gates.you.know.what
    Aug 11, 2015 at 10:20
  • $\begingroup$ It is not an unreasonable question, but I am not sure it is a Mathematica question... $\endgroup$
    – Oleksandr R.
    Aug 11, 2015 at 10:56
  • $\begingroup$ I have tried a lot of "by hand" methods but I must, to avoid human mistake, left to Mathematica to find "border value" (X*). $\endgroup$
    – danijel
    Aug 11, 2015 at 11:50
  • $\begingroup$ To Oleksandr R.: I think it is Mathematica question because in forest fire (and mathematical - physicas similar problems, like spread of oil on sea for example), if fire behavior is like stretched exponential (where is "A" value, witch represent area) then fire is "under control" until area of that size. If it is similar to Power law, who have no scale, than small fire behavior is similar to big one (and vice versa) and fire behavior is not independent of human activities, it depends only on "outside" factors like (big) rain or natural border (like water, big slope bold rock ect.) $\endgroup$
    – danijel
    Aug 11, 2015 at 12:04
  • $\begingroup$ Sory, in last post "is not independent" should be "is independent" $\endgroup$
    – danijel
    Aug 11, 2015 at 12:05

1 Answer 1

Reset to default
4
$\begingroup$

see comment below by OleskandrR

If he posts I will up vote. Original post

This is not ideal (and done with little time) but may motivate better answers. Note this manually (visually tries to minimize the discontinuity at "breakpoint") using Manipulate. The NonlinearModelFit treats the breakpoint as fixed so you can play and chose better way. I look forward to better ways but am too time poor to think about.

Set up:

data = {{1.1*10^6, 1/21}, {600., 3/7}, {80., 5/7}, {100., 
    29/42}, {600., 3/7}, {600., 3/7}, {70., 31/42}, {500., 
    19/42}, {25., 37/42}, {25., 37/42}, {300., 4/7}, {25., 
    37/42}, {15000., 5/42}, {150., 2/3}, {1400., 13/42}, {10., 
    1}, {60., 16/21}, {200., 13/21}, {10000., 4/21}, {10000., 
    4/21}, {1000., 5/14}, {50., 17/21}, {1500., 2/7}, {200., 
    13/21}, {300., 4/7}, {5000., 3/14}, {10., 1}, {80000., 
    1/14}, {20., 13/14}, {1575., 11/42}, {400., 1/2}, {1000., 
    5/14}, {9.1671*10^6, 1/42}, {60000., 2/21}, {20., 13/14}, {50., 
    17/21}, {400., 1/2}, {150., 2/3}, {10., 1}, {10000., 4/21}, {350.,
     11/21}, {2200., 5/21}};

Putative model function:

f[x_, a_, b_, c_, d_, g_] := 
 Boole[x < d] Exp[-(x/a)^b] + (1 - Boole[x < d]) c x^(-g)

NonlinearModelFit with startiing values from playing with Manipulate:

nlm = NonlinearModelFit[data, 
  f[x, a, b, c, d, 
   g], {{a, 836}, {b, 0.41}, {c, 3.87}, {d, 1980}, {g, 0.36}}, x]

Note:

nlm["AdjustedRSquared"] yields: 0.996793

nlm["ParameterTable"]:

enter image description here

Visualizing:

Show[
 ListLogLogPlot[data], 
 LogLogPlot[nlm2[x], {x, 10, 10000000}, PlotRange -> {0, 1}]]

Normal@nlm:

(2.28415 (1 - Boole[x < 1980.]))/x^0.282831 + 
 E^(-0.0243126 x^0.551201) Boole[x < 1980.]

Good luck...note again this just uses the cut point from Manipulate look forward to better answers

enter image description here

Improve this answer
$\endgroup$
5
  • $\begingroup$ It is probably better to enforce continuity as a constraint and let NonlinearModelFit find the value of d appropriately to that. But the big problem with this model is that it is overdetermined and so d is not really treated as a parameter. nlm=NonlinearModelFit[data,{f[x,a,b,c,d,g],Exp[-(d/a)^b]==c d^(-g)&&a>1&&b>0&&c>0&&d>1},{{a,500},{b,1},{c,3},{d,200},{g,0.3}}, x,Method->{NMinimize,Method->"NelderMead"}] gives a different result, which looks better to my eye, although I would still have zero confidence that it is correct. $\endgroup$
    – Oleksandr R.
    Aug 11, 2015 at 12:26
  • $\begingroup$ @OleksandrR. Yes I agree with constraint and also 'zero confidence'...just thought might illustrate a start...without knowing what data is and why model plausible or what aim is...just play...thank you for your comment as I always learn... $\endgroup$
    – ubpdqn
    Aug 11, 2015 at 12:32
  • $\begingroup$ I want to note that my comment was not intended as a criticism of your answer, which I think shows basically the right approach. I just meant to point out that there are difficulties in applying the approach straightforwardly because of the mathematical structure of the model. $\endgroup$
    – Oleksandr R.
    Aug 11, 2015 at 12:47
  • $\begingroup$ @OleksandrR. I appreciated and agreed with your comment and just wanted to direct the answer to your comment which I think gets closer to the aim...very time poor at present... $\endgroup$
    – ubpdqn
    Aug 11, 2015 at 12:50
  • $\begingroup$ Thanks for both of you for all comments and suggestions. I was hoping that there is some specific (and "simply" :) order fort his type of a problem. I've tried to find a solution on many web pages; on this also, I will search further... And thanks for all future suggestions! $\endgroup$
    – danijel
    Aug 11, 2015 at 12:56

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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