#13273: extension of sage.numerical.optimize.find_fit
-----------------------------------------------------------+----------------
Reporter: hackstein | Owner:
hackstein
Type: enhancement | Status:
new
Priority: major | Milestone:
sage-5.3
Component: numerical | Resolution:
Keywords: fitting, interpolation, two-dimensional | Work issues:
enhancement
Report Upstream: N/A | Reviewers:
tba
Authors: hackstein | Merged in:
Dependencies: | Stopgaps:
-----------------------------------------------------------+----------------
Comment (by hackstein):
A first attempt would be the following, but it still needs careful review:
def find_fit2dim(data, model, initial_guess = None, parameters = None,
variables = None, solution_dict = False):
r"""
Finds numerical estimates for the parameters of the two dimensional
function model to
give a best fit to data.
INPUT:
- ``data`` -- A four dimensional table of floating point numbers of the
form `[[x_{1,1}, x_{1,2},f_{1,1}, f_{1,2}],
[x_{2,1}, x_{2,2}, f_{2,1},f_{2,2}],
\ldots,
[x_{n,1}, x_{n,2},f_{n,1},f_{n,2}]]` given as a
numpy array.
- ``model`` -- Either a pair of symbolic expressions, of symbolic
functions, or of
Python functions. ``model`` has to be a pair of functions of the
variables
`(x_1, x_2, \ldots, x_k)` and free parameters
`(a_1, a_2, \ldots, a_l)`.
- ``initial_guess`` -- (default: ``None``) Initial estimate for the
parameters `(a_1, a_2, \ldots, a_l)`, given as either a list, tuple,
vector or numpy array. If ``None``, the default estimate for each
parameter is `1`.
- ``parameters`` -- (default: ``None``) A list of the parameters
`(a_1, a_2, \ldots, a_l)`. If model is a symbolic function it is
ignored, and the free parameters of the symbolic function are used.
- ``variables`` -- (default: ``None``) A list of the variables
`(x_1, x_2, \ldots, x_k)`. If model is a symbolic function it is
ignored, and the variables of the symbolic function are used.
- ``solution_dict`` -- (default: ``False``) if ``True``, return the
solution as a dictionary rather than an equation.
EXAMPLES:
ALGORITHM:
Uses ``scipy.optimize.leastsq`` which in turn uses MINPACK's lmdif and
lmder algorithms.
"""
import numpy
if not isinstance(data, numpy.ndarray):
try:
data = numpy.array(data, dtype = float)
except (ValueError, TypeError):
raise TypeError, "data has to be a numpy array"
elif data.dtype == object:
raise ValueError, "the entries of data have to be of type float"
if data.ndim != 4:
raise ValueError, "data has to be a four dimensional table of floating
point numbers"
from sage.symbolic.expression import Expression
if isinstance(model[1], Expression) and if
isinstance(model[2],Expression):
if variables is None:
variables = list(x_1,x_2)
if parameters is None:
parameters =
list(Set(list(model[1].variables())+list(model[2].variables())))
for v in variables:
parameters.remove(v)
if data.shape[1] != len(variables) + 2:
raise ValueError, "each row of data needs %d entries, only %d entries
given" % (len(variables) + 2, data.shape[1])
if parameters is None or len(parameters) == 0 or \
variables is None or len(variables) == 0:
raise ValueError, "no variables given"
if initial_guess == None:
initial_guess = len(parameters) * [1]
if not isinstance(initial_guess, numpy.ndarray):
try:
initial_guess = numpy.array(initial_guess, dtype = float)
except (ValueError, TypeError):
raise TypeError, "initial_guess has to be a list, tuple, or numpy array"
elif initial_guess.dtype == object:
raise ValueError, "the entries of initial_guess have to be of type float"
if len(initial_guess) != len(parameters):
raise ValueError, "length of initial_guess does not coincide with the
number of parameters"
if isinstance(model[1], Expression) and if isinstance(model[2],
Expression):
var_list = variables + parameters
var_names = map(str, var_list)
func = model._fast_float_(*var_names)
else:
func = model
def function(x_data, params):
result = numpy.zeros(len(x_data))
for row in xrange(len(x_data)):
fparams = numpy.hstack((x_data[row], params)).tolist()
result[row] = func(*fparams)
return result
def error_function(params, x_data, y1_data, y2_data):
result = numpy.zeros(len(x_data))
for row in xrange(len(x_data)):
fparams = x_data[row].tolist() + params.tolist()
result[row] = func(*fparams)
return np.array([result[:,-2] - y1_data,result[:,.1]-y2_data])
x_data = data[:, 0:len(variables)]
y1_data = data[:, -2]
y2_data=data[:,-1]
from scipy.optimize import leastsq
estimated_params, d = leastsq(error_function, initial_guess, args =
(x_data, y1_data,y2_data))
if isinstance(estimated_params, float):
estimated_params = [estimated_params]
else:
estimated_params = estimated_params.tolist()
if solution_dict:
dict = {}
for item in zip(parameters, estimated_params):
dict[item[0]] = item[1]
return dict
return [item[0] == item[1] for item in zip(parameters, estimated_params)]
--
Ticket URL: <http://trac.sagemath.org/sage_trac/ticket/13273#comment:2>
Sage <http://www.sagemath.org>
Sage: Creating a Viable Open Source Alternative to Magma, Maple, Mathematica,
and MATLAB
--
You received this message because you are subscribed to the Google Groups
"sage-trac" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to
[email protected].
For more options, visit this group at
http://groups.google.com/group/sage-trac?hl=en.
