CODEWITHSUNDEEP
pythonscipycurve-fitting
physics_for_all
physics_for_all

Reputation: 2263

Power law with a constant factor using curve_fitting

I want to fit my x and y data using power law with a constant factor. My power law model is y(r) = F0 + F*(r)**alpha where F0 is a constant. My code is,

x = [0.015000000000000001, 0.024999999999999998, 0.034999999999999996, 0.044999999999999998, 0.055, 0.065000000000000002, 0.075000000000000011, 0.085000000000000006, 0.094999999999999987, 0.125, 0.17500000000000002, 0.22500000000000003, 0.27500000000000002]

y= [5.6283727993522774, 4.6240796612752799, 3.7366642904247769, 3.0668203445969828, 2.5751865553847577, 2.0815063796430979, 1.7152655187581032, 1.4686235817532258, 1.2501921057958358, 0.80178306738561222, 0.43372429238424598, 0.26012305284446235, 0.19396186239328625]


def func(x,m,c,c0):
      return c0 + x**m * c

coeff,var=curve_fit(func,x,y)

print coeff


Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/usr/local/lib/python2.6/dist-packages/scipy/optimize/minpack.py", line 511, in curve_fit
    raise RuntimeError(msg)
RuntimeError: Optimal parameters not found: Number of calls to function has reached maxfev = 800.

Then I changed maxfev=2000 then it gives me wrong coeff values. If I change, my slope m to (-m) in func then it gives me right answer but my slope will be negative. Is there any way to overcome this problem?

Upvotes: 3

Views: 3746

Answers (2)

kiriloff
kiriloff

Reputation: 26333

Looking at your data, an idea would be to rescale your problem, given the order of magnitude of difference between X and Y.

An example is given in details in first post.

Upvotes: 0

David Robinson
David Robinson

Reputation: 78620

The issue is that curve_fit is starting with default guesses for the parameters that are too poor (namely, it starts them all at 1).

Instead, use what you know about the data to make a very rough guess: at the very least, you know m must be negative (since it's a power law). Therefore, try starting m at -1. (You can start the intercept term at 0 and the slope term at 1 since those are reasonable defaults).

def func(x,m,c,c0):
    return c0 + x**m * c

coeff,var=curve_fit(func,x,y, p0=[-1, 1, 0])

This gives you the correct output:

[-0.34815029  2.16546037 -3.4650323 ]

(Note that you can start m with any number between 0 and -9 and it still converges to this result).

Upvotes: 6

Related Questions