Plotting and Fitting

edit

Have a look at find_fit on the page http://www.sagemath.org/doc/reference...

You could also use the R interface in Sage to produce things, using their standard fitting and error bar tools.

Edit: R is probably the industry standard for statistical/modeling/fitting stuff. Of the statistical variety. http://www.r-project.org is the site. You can do

sage: r?

for basic info, and see the whole documentation at this page for more info (apparently it's not in the reference manual yet. I also have a sample at sagenb, but that site is so overloaded that it might be hard to access...

I'm not sure what you mean by "with considering the errorbars". Do you mean merely that you want the covariance of the fit (which leastsq returns), or that your data has errors in y (the usual model) and you want that to affect the fitting process, or that your data has errors in both x and y and you want that to affect the fit?

Anyway, here are two relatively straightforward ways to fit and get covariances.

sage: import numpy as np
sage: import scipy.optimize
sage: import scipy.odr
sage: import random
sage: random.seed(3)
sage: 
sage: #
sage: # generate the data
sage: #
sage: 
sage: x = np.linspace(0, 10, 100)
sage: y = np.sin(0.7*x)+2.4*x  # underlying function
sage: y += np.random.randn(len(x))/10 # add some noise
sage: x_err = [0.1]*len(x) # estimated errors
sage: y_err = [0.3]*len(x) # estimated errors
sage: 
sage: # first, we can use curve_fit
sage: 
sage: def fit1(x, a, b):
....:         return np.sin(x*a)+b*x
....: 
sage: # note: can't handle errors in x
sage: popt, pcov = scipy.optimize.curve_fit(fit1, x, y, sigma=y_err)
sage: popt
array([ 0.70124584,  2.39768117])
sage: pcov
array([[  6.97634001e-06,  -2.58153187e-06],
       [ -2.58153187e-06,   4.22745005e-06]])
sage: 
sage: # or use the more powerful odr approach, which unfortunately has
sage: # some different conventions
sage: 
sage: def fit2(B, x):
....:         return np.sin(x*B[0])+B[1]*x
....: 
sage: mymodel = scipy.odr.Model(fit2)
sage: mydata = scipy.odr.RealData(x, y, sx=x_err, sy=y_err)
sage: myodr = scipy.odr.ODR(mydata, mymodel, beta0=[1., 1.], maxit=1000)
sage: myoutput = myodr.run()
sage: myoutput.pprint()
Beta: [ 0.70139136  2.39769177]
Beta Std Error: [ 0.00258362  0.00198363]
Beta Covariance: [[  9.08706871e-05  -2.35028948e-05]
 [ -2.35028948e-05   5.35658971e-05]]
Residual Variance: 0.0734568446584
Inverse Condition #: 0.679929358499
Reason(s) for Halting:
  Numerical error detected

Here is an old thread from sage-devel about error bars on plots, with some suggestions. One of them suggests using matplotlib, and perhaps this would also work for curve fitting.

Sorry for the short answer, but I don't know much about this :(