4
$\begingroup$

I am trying to fit some raw data using a function of the form

$f(r) = \sum_{i=1}^{K} d_kS_k(n_k,\alpha_k,r)$

where

$S_k(n_k,\alpha_k,r) = \frac{\alpha_k ^{n_k+3}}{(n_k+2)!}r^{n_k}\exp(-\alpha_kr)$

(note that the actual fit function is far more complex).

Now, for each $k$, this has two non-linear parameters ($\alpha_k$ and $n_k$) and one linear parameters ($d_k$). I am considering various methods for non-linear optimization (genetic algorithms, conjugate gradient method, etc.) but I am wondering how I should work with the linear functions.

Should I include them in the non-linear optimization (use them as a variables when I take the gradient for CGM, etc)? Should I optimize the non-linear parameters completely on $\alpha$ and $n$ and then perform a linear fit with $d$? Should I mix and match - perform one sequence of conjugate gradient steps and then perform a linear fit, rinse and repeat?

NOTE: If it helps at all, $d_k$ is a real number, $n_k$ is a positive integer, $\alpha_k$ is a positive real.

|cite|improve this question
$\endgroup$
  • $\begingroup$ In what sense do you want to fit? Least-squares? $\endgroup$ – David Ketcheson Dec 15 '11 at 15:53
  • $\begingroup$ you could kind of think of it as least squares. there is a least squares part for comparing the domain on which the data points are defined, then a squared integral part for comparing the contrived function over the rest of the space (final data point to infinity). even this is a simplification, but you can think of it as a least squares fit. $\endgroup$ – drjrm3 Dec 15 '11 at 16:36
3
$\begingroup$

For the least-squares case anyway, one standard method is Linda Kaufman's variable projection algorithm (see this as well). This is based on earlier work by Gene Golub and Victor Pereyra for the separable nonlinear least squares problem, in which they show how the structure of the Moore-Penrose pseudoinverse can be exploited to yield algorithms that are more efficient when there are linear parameters present. This is one FORTRAN implementation of the method.

(But see also this paper by Ruhe and Wedin for a different method.)

|cite|improve this answer
$\endgroup$
2
$\begingroup$

You should optimize all parameters simultaneously if you would like to obtain a (local epsilon-)optimal solution to the optimization problem corresponding to the fit you'd like to perform. J. M.'s suggestion is a good one for the least-squares case, but if the specifics of your fit aren't quite least-squares, the general strategy of fitting all parameters simultaneously will yield a better fit than optimizing groups of parameters separately.

If $n_k$ is an integer, you'll want to use a mixed-integer nonlinear programming (MINLP) solver (or algorithm). Including an integral in your objective in an MINLP solver will be difficult unless you can evaluate it numerically in closed form, which will limit you to fixed quadrature. MINLP solvers require that problems be posed explicitly as equations, not as algorithms. (It is possible to optimize a problem posed as part of an algorithm using the algorithms in this paper, but the class of algorithms supported does not include algorithms that use adaptive quadrature.)

|cite|improve this answer
$\endgroup$
  • $\begingroup$ At the very least, if the actual norm being considered isn't the 2-norm, the least-squares result can be used as a seed for optimizing in some other norm... $\endgroup$ – J. M. Dec 15 '11 at 19:39
  • $\begingroup$ @J.M.: What do you mean by a "seed"? As an initial guess? $\endgroup$ – Geoff Oxberry Dec 15 '11 at 20:18
  • $\begingroup$ i had actually been hoping to treat $n_k$ as a continuous variable during the calculation, but restrict it to an integer (by simply rounding a double precision variable) every once in a while. in this case the integration (which i can transform to a closed integral on 0,1) will be easier, as will the derivatives if i choose the conjugant gradient route. $\endgroup$ – drjrm3 Dec 15 '11 at 20:45
  • $\begingroup$ @Geoff: Yes. Iterative methods always need a good starting point, don't they? $\endgroup$ – J. M. Dec 15 '11 at 21:28
  • 1
    $\begingroup$ @Laurbert515: Rounding $n_{k}$ may be fine for your purposes. In general, when solving MINLPs, relaxing the integer variables to continuous variables, solving the resulting nonlinear program, and then rounding the relaxed integer variables back to the nearest integer (i.e., the procedure you're using now) does not necessarily yield an optimal solution. It just depends on how well you want to fit your model, and how much wall clock time you want to spend doing so. I presume that since you're treating $n_{k}$ as continuous that you're also substituting the gamma function for the factorial. $\endgroup$ – Geoff Oxberry Dec 15 '11 at 23:27

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.