python - Scipy minimize, fmin, leastsq type problems (setting array element with sequence), bad fit -
i'm having trouble scipy.optimize.fmin , scipy.optimize.minimize functions. i've checked , confirmed arguments passed function of type numpy.array, return value of error function. also, carreau function returns scalar value.
the reason of arguments, such size, this: need fit data given model (carreau). data taken @ different temperatures, corrected shift factor (which fitted model), end several sets of data should used calculate same 4 constants (parameters p).
i read can't pass fmin function list of arrays, had concatenate data x_data_lin, keeping track of different sets size parameter. t holds different test temperatures, while t_0 one-element array holds reference temperature.
i positive (triple checked) arguments passed function, result, one-dimensional arrays. here's code aside that:
import numpy np import scipy.optimize scipy.optimize import fmin simplex def err_func2(p, x, y, t, t_0, size): result = array([]) temp = 0 in range(0, int(len(size)-1)): j in range(int(temp), int(temp+size[i])): result = np.append(result, (carreau(p, x[j], t[i], t_0[0])-y[i])) temp += size[i] return result p1 = simplex(err_func2, initial_guess, args=(x_data_lin, y_data_lin, t_list, t_0, size), full_output=0)
here's error:
traceback (most recent call last): file "c:\python27\scripts\projects\carreau - wlf\carreau_model_fit.py", line 146, in <module> main() file "c:\python27\scripts\projects\carreau - wlf\carreau_model_fit.py", line 105, in main args=(x_data_lin, y_data_lin, t_list, t_0, size), full_output=0) file "c:\python27\lib\site-packages\scipy\optimize\optimize.py", line 351, in fmin res = _minimize_neldermead(func, x0, args, callback=callback, **opts) file "c:\python27\lib\site-packages\scipy\optimize\optimize.py", line 415, in _minimize_neldermead fsim[0] = func(x0) valueerror: setting array element sequence.
it's worth noting got leastsq function working while passing lists of arrays. unfortunately, did poor job of fitting data. but, took me lot of time , research point, i'll post code follows. if interested in seeing of code, gladly post it, if can recommend me somewhere upload few files(as includes imported script , of course sample data):
##def error_function(p, x, y, t, t_0): ## result = array([]) ## index in range(len(x)): ## result = np.append(result,(carreau(p, x[index], ## t[index], t_0) - y[index])) ## return result ## p1, success = scipy.optimize.leastsq(error_function, initial_guess, ## args=(x_list, y_list, t_list, t_0), ## maxfev=10000)
:( going post picture of graphed data leastsq fit, don't have requisite 10 points.
late edit: have gotten optimize.curvefit , optimize.leastsq work (which not-so-coincidentally give same answer), curve bad. i've been trying figure out optimize.minimize, it's been bit of headache. simplex (fmin, nelder_mead, whatever want call it) run, produces crazy answer close. i've never worked nonlinear optimization problems before, , don't know direction head.
here's working curve_fit code:
def temp_shift(t_s, t, t_0): """ function calculates a_t temperature shift factor polymer viscosity curves. variable standard temperature, t_s """ c_1 = 8.86 c_2 = 101.6 return(np.exp( (c_1*(t_0-t_s) / (c_2+(t_0-t_s))) - (c_1*(t-t_s) / (c_2 + (t-t_s))) )) def pass_data(t, t_0): def carreau_2(x, p0, p1, p2, p3): visc_0 = p0 m = p1 n = p2 t_s = p3 a_t = temp_shift(p3, t, t_0) return (visc_0 * a_t / (1 + m * x * a_t)**n) return carreau_2 initial_guess = array([20000, 3, 0.94, -20]) p1, conv = scipy.optimize.curve_fit(pass_data(t_all, t_0), x_data_lin, y_data_lin, initial_guess)
here's sample data:
x_data_lin = array([0.01998, 0.04304, 0.2004, 0.43160, 0.92870, 2.0000, 4.30900, 9.28500, 15.51954, 21.94936, 37.52960, 90.41786, 204.35230, 331.58495, 811.92250, 1694.55309, 3464.27648, 8826.65738, 14008.00242]) y_data_lin = array([13520.00000, 13740.00000, 12540.00000, 9384.00000, 5201, 3232.00000, 2094.00000, 1484.00000, 999.00000, 1162.05088 942.56946, 705.62489, 429.47341, 254.15136, 185.22916, 122.07113, 76.46324, 47.85064, 25.74315, 18.84875]) t_all = array([190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190, 190]) t_0 = 80
here's picture of result of curve_fit (now have 10 points , can post!). note there 3 curves drawn because used 3 sets of data optimize curve, @ 3 different temperatures. polymers have property shear_rate - viscosity relationship stays same, shifted temperature factor a_t:
i'd appreciate suggestions how improve fit, or how define function optimize.minimize works, , method (nelder-mead, powel, bfgs) might work.
another edit add: got nelder-mead (optimize.fmin, , default of optimize.minimize) function work - i'll include revised error function below. before, summed result array , returned it. led extremely negative values (obviously, since function's goal minimize). squaring result before summing solved problem. note changed function take advantage of numpy's array broadcasting, suggested jaminsore (thanks jamin!)
def err_func2(p, x, y, t, t_0): return ((carreau(p, x, t, t_0)-y)**2).sum()
unfortunately, nelder-mead function gives me same result leastsq , curve_fit. can see in graph above it's not optimal fit; in fact, @ point, microsoft excel's solver function doing better job on data.
at least, hope thread can useful beginners scipy.optimize in future, since it's taken me quite awhile discover of this.
unlike leastsq
, fmin
can deal error functions return scalar if possible have rewrite error function returns scalar. here simple working example.
import necessary libraries
import numpy np scipy.optimize import fmin
define helper function (you'll see later)
def prob(a, b): return (1 + np.exp(b - a))**-1
simulate data
true_ = np.random.normal(size = 100) #parameters we're trying recover b = np.random.normal(size = 20) exp_ = prob(true_[:, none], b) #expected a_s, b_s = true_.shape[0], b.shape[0] noise = np.random.uniform(size = (a_s, b_s)) response = (noise > (1 - exp_)).astype(int)
define our error function (i'm using lambda
s not recommended in practice)
# sum of squared residuals err_func = lambda : ((prob(a[:, none], b) - response) ** 2).sum() result = fmin(err_func, np.zeros_like(true_)) #solve
if remove .sum()
@ end of error function definition, same error.
Comments
Post a Comment