# optimization (using python) how to find the most suitable solver?

I have a fitting routine set up. It works, but pretty slow. I was wondering if there is a better method to use. I checked my (forward) code against some literature data and at least I do have no bugs and mistakes in my problem.

I fit the steric mass action model to some some isotherms like the ones shown below. Additionally, i assume some arbitrary numbers of unknown components and fit just the measured concentrations of my known components (it is a multicomponent system). Additionally, the measured compound itself consists of some charge variants with different adsorption characteristics. Therefore, I am playing around right know with the number of components.

The model for a three component system takes following form. (the first component (the salt concentration) has to be calculated numerically)

\begin{align} qs &=& \frac{\Lambda \, c_1}{c_1 + (\sigma_2 + \nu_2)* \{ k_2 \, ( \frac{qs}{c_1} )\}^{\nu_2 - 1} + (\sigma_3 + \nu_3)* \{ k_3 \, ( \frac{qs}{c_1} )\}^{\nu_3 - 1} } \\ \end{align}

\begin{align} q2 &=& \frac{\Lambda \, {k \, ( \frac{qs}{c_1} )^{\nu_2 - 1} } }{c_1 + (\sigma_2 + \nu_2)* \{ k_2 \, ( \frac{qs}{c_1} )\}^{\nu_2 - 1} + (\sigma_3 + \nu_3)* \{ k_3 \, ( \frac{qs}{c_1} )\}^{\nu_3 - 1} } \\ \end{align}

\begin{align} q3 &=& \frac{\Lambda \, {k \, ( \frac{qs}{c_1} )^{\nu_3 - 1} } }{c_1 + (\sigma_2 + \nu_2)* \{ k_2 \, ( \frac{qs}{c_1} )\}^{\nu_2 - 1} + (\sigma_3 + \nu_3)* \{ k_3 \, ( \frac{qs}{c_1} )\}^{\nu_3 - 1} } \\ \end{align}

qs is the stationary concentration of free ligands [mM] typically between 50 and 1500

$\Lambda$ is initial concentration of ligands; can be a fit parameter; typically between 100 and 1700

$k$ is a fit parameter typically between 0.1 and 40

$\sigma$ is a fit parameter typically between 1 and 50

$\nu$ is a fit parameter typically between 1 and 15

Since the solver fails on some sets of input parameters and i do not want to fiddle around, i use scipys basin-hopping solver together with the method L-BFGS-B (quite slow). The parameters are bound so they make physical sense. Starting points are generated randomly within the bounded domain.

I was wondering how someone used to optimization would set this problem up ? Would you try to provide an analytical gradient ?

Edit:

I do minimize the sum of squares. I thought about scipys leastsq, but it does not accept bounds. From the local solvers, fmin (using L-BFGS-B) was the most robust and fastest one.

• scipy.optimize.curve_fit? – kηives Feb 3 '15 at 16:30
• As far as i understood, curve_fit is a wrapper for lstsq and does not handle bounds. I tried it once, but i get the error Improper input: N=7 must not exceed M=1 – Moritz Feb 3 '15 at 21:25