python - exponential decay fitting -

i trying fit data distributed in time following exponential decay. tried follow fitting examples on web, code doesn't fit data. straight line results fit. maybe there wrong initial parameters? until have used gaussian , line fits, using same method, maybe not correct case. code take data web, directly executable. question: why doesn't code result in fit? many in advance.

#!/usr/bin/env python  import pyfits, os, re, glob, sys scipy.optimize import leastsq numpy import * pylab import * scipy import *  rc('font',**{'family':'serif','serif':['helvetica']}) rc('ps',usedistiller='xpdf') rc('text', usetex=true) #------------------------------------------------------  tmin = 56200 tmax = 56249  data=pyfits.open('http://heasarc.gsfc.nasa.gov/docs/swift/results/transients/weak/gx304-1.orbit.lc.fits') time  = data[1].data.field(0)/86400. + data[1].header['mjdreff'] + data[1].header['mjdrefi'] rate  = data[1].data.field(1) error = data[1].data.field(2) data.close()  cond = ((time > 56210) & (time < 56225)) time = time[cond] rate = rate[cond] error = error[cond]  right_exp = lambda p, x: p[0]*exp(-p[1]*x) err = lambda p, x, y:(right_exp(p, x) -y) v0= [0.20, 56210.0, 1] out = leastsq(err, v0[:], args = (time, rate), maxfev=100000, full_output=1) v = out[0] #fit parameters out xxx = arange(min(time), max(time), time[1] - time[0]) ccc = right_exp(v, xxx) fig = figure(figsize = (9, 9)) #make plot ax1 = fig.add_subplot(111) ax1.plot(time, rate, 'g.') #spectrum ax1.plot(xxx, ccc, 'b-') #fitted spectrum savefig("right exp.png")  axis([tmin-10, tmax, -0.00, 0.45]) 

your problem ill conditioned because array times contains big numbers when used in exp(-a*time) giving values close 0., tricks err function because rate array contains small values close 0., leading small errors. in other words, high a in exponential function gives solution.

to fix can:

• change decay function include initial time:
  exp(-a*(time-time0))
• change input data start smaller number:
  time -= time.min()

for both options have change initial guess v0, e.g. v0=[0.,0.]. first solution seems more robust , not have manage changes in time array. initial guess time0 time.min():

right_exp = lambda p, x: p[0]*exp(-p[1]*(x-p[2])) err = lambda p, x, y:(right_exp(p, x) -y) v0= [0., 0., time.min() ] out = leastsq(err, v0, args = (time, rate)) v = out[0] #fit parameters out xxx = arange(min(time), max(time), time[1] - time[0]) ccc = right_exp(v, xxx) fig = figure(figsize = (9, 9)) #make plot ax1 = fig.add_subplot(111) ax1.plot(time, rate, 'g.') #spectrum ax1.plot(xxx, ccc, 'b-') #fitted spectrum fig.show() 

giving:

still, final results depending on v0, e.g. v0=[1.,1.,time.min()] decays fast , not find optimum.