Index: /trunk/doc/pslib/psLibADD.tex
===================================================================
--- /trunk/doc/pslib/psLibADD.tex	(revision 4176)
+++ /trunk/doc/pslib/psLibADD.tex	(revision 4177)
@@ -1,3 +1,3 @@
-%%% $Id: psLibADD.tex,v 1.74 2005-06-03 02:34:13 price Exp $
+%%% $Id: psLibADD.tex,v 1.75 2005-06-09 05:16:23 eugene Exp $
 \documentclass[panstarrs]{panstarrs}
 
@@ -433,53 +433,123 @@
 \subsection{Non-linear Minimization}
 
+Non-linear minimization techniques use an iterative approach to find a
+minimization since an analytical inversion is impractical or not
+possible.  These techniques use a starting guess for the parameters of
+interest, and make a sequence of new guess parameters based on the
+properties of the function at the previous position.  If the new
+parameters yield a reduced function value, the new position is used as
+the starting position for the next iteration.  Otherwise, the guess
+must be modified and another attempt is made.  Convergence may be
+determined based on the absolute amount of change in the function
+value, or by comparison with the expectation for a linear system. 
+
+The two common techniques used to construct a guess parameter set are
+the 'steepest descent method' and the 'gauss newton method'.  In the
+first case, the guess is selected some distance along the local
+gradient.  In the second case, a local Taylor expansion of the
+function is used to construct a linear model for the function, and the
+new guess is chosen to minimize that model.  The methods discussed
+below make use of combinations of these two methods.  Aside from the
+differences in their guess steps, the two methods differ in using
+first derivatives of the function supplied by the user, or by locally
+calculating the first derivatives.  
+
+Mathematically, we would like to choose the parameter set $a_m$ to
+minimize a function of those parameters $F(a_m)$.  We iterate by
+choosing a new parameter set $a^\prime_m = a_m + \delta_m$ based on
+the behavior of the function at $a_m$.  
+
+The steepest descent method chooses a step direction of $\bar{\delta}
+= -\nabla F$, or defining $g_m$ as a component of the gradient,
+$\delta_m = -g_m$.  The Gauss-Newton method uses a Taylor expansion of
+the function to solve for the step: $(\nabla^2 F) \bar{\delta} =
+-\nabla F$.  Defining $H_{m,n}$ as a component of the second
+derivatives, we can write the Gauss-Newton step as $\delta_m = -
+H_{m,n}^{-1} g_n$.
+
 \subsubsection{Levenberg-Marquardt Method}
 
-In the Levenberg-Marquardt Method (LMM; see NR \S 15.5), we make a
-guess at the input parameters, evaluate the function of interest, vary
-the parameters by a particular choice based on the gradient, evaluate
-the function again, and adjust the parameters and the parameter
-varient based on the results.
-
-The LMM only works if the second derivative of the function can be
-considered negligible, as in the case of minimizing $\chi^2$.
-
-Given some ordinates, $x_i$, we would like to find the parameters,
-$a_k$, of the function $f(x_i; a_k)$ which minimize $\chi^2$ for some
-measurements, $y_i$ and associated errors, $\sigma_i$.  We start with
-a set of parameter guesses, $a_k$.  We calculate the gradient
-$\beta_k$ and the Hessian matrix $\alpha_{j,k}$ at this parameter
-selection as follows:
-\begin{eqnarray}
-\beta_k & = & \frac{\partial \chi^2}{\partial a_k} \\
-\alpha_{j,k} & = & \sum_i \frac{1}{\sigma_i^2} \frac{\partial f(c_i)}{\partial a_k} \frac{\partial f(c_i)}{\partial a_j}
-\end{eqnarray}
-
-We now define the new parameter guess for $a_k$ based on the gradient
-and Hessian by defining $A_{j,k}$ as a variant on $\alpha_{j,k}$ as
-follows:
-
-\begin{eqnarray}
-A_{j,k} & = & \alpha_{j,k} ~ (j \ne k) \\
-A_{j,k} & = & (1 + \lambda) \alpha_{j,k} ~ (j = k)
+In the Levenberg-Marquardt Method (LMM; see NR \S 15.5, Madsen et al),
+we make a guess at the input parameters, evaluate the function of
+interest, vary the parameters by a particular choice based on the
+gradient, evaluate the function again, and adjust the parameters and
+the parameter varient based on the results.  The LMM only works if the
+second derivative of the function can be considered negligible, as in
+the case of minimizing $\chi^2$.
+
+Consider the chi-square function, given some ordinates, $x_i$, we
+would like to find the parameters, $a_m$, of the function which
+minimize $\chi^2$ for some measurements, $y_i$ and associated errors,
+$\sigma_i$:
+\begin{eqnarray}
+\chi^2 (\bar{a})  & = & \sum_i \frac{1}{\sigma_i^2} \left( y_i - p(x_i;a_m) \right)^2  \\
+\end{eqnarray}
+
+We simplify this as:
+\begin{eqnarray}
+p_i (a_m)         & = & p(x_i;a_m) \\
+f_i (a_m)         & = & \frac{1}{\sigma_i} (y_i - p_i) \\
+\chi^2 (\bar{a})  & = & \sum_i f_i^2 \\
+\end{eqnarray}
+
+We write the minimization function $F(a_m) = \frac{1}{2} \chi^2$ to
+avoid various extra factors of 2.  We can now write out the needed
+derivatives in terms of $f_i$:
+\begin{eqnarray}
+F(a_m) = \frac{1}{2} \sum_i f_i^2 \\
+\nabla F(a_m) = \sum_i f_i \frac{\partial f_i}{\partial a_m} \\
+\nabla^2 F(a_m) = \sum_i \frac{\partial f_i}{\partial a_m} \frac{\partial f_i}{\partial a_n} \\
+\end{eqnarray}
+
+where we have dropped the second-derivatives of the function in the
+representation of $\nabla^2 F(a_m)$.  Since $\frac{\partial
+f_i}{\partial a_m} = -\frac{1}{\sigma_i}\frac{\partial p_i}{\partial
+a_m}$, we can write these in terms of the derivatives of $p_i$ only:
+\begin{eqnarray}
+\nabla F(a_m)   = -\sum_i \frac{f_i}{\sigma_i} \frac{\partial p_i}{\partial a_m} \\
+\nabla^2 F(a_m) = \sum_i \frac{1}{\sigma_i^2} \frac{\partial p_i}{\partial a_m} \frac{\partial p_i}{\partial a_n} \\
+\end{eqnarray}
+
+Writing these in matrix representation, and replacing $f_i$, we have:
+\begin{eqnarray}
+-g_m    = \sum_i \frac{(y_i - p_i)}{\sigma_i^2} \frac{\partial p_i}{\partial a_m} \\
+H_{m,n} = \sum_i \frac{1}{\sigma_i^2} \frac{\partial p_i}{\partial a_m} \frac{\partial p_i}{\partial a_n} \\
+\end{eqnarray}
+
+In the Levenberg-Marquart Method, we define a new guess using a
+combination of the Steepest Descent and Gauss-Newton methods discussed
+above.  We replace the Hessian matrix above with $A_{m,n}$ as a
+variant on $H_{m,n}$ as follows:
+
+\begin{eqnarray}
+A_{m,n} & = & H_{m,n} if (j \ne k) \\
+A_{m,n} & = & H_{m,n}(1 + \lambda) if (j = k)
 \end{eqnarray}
 %
 and solve the system of equations represented by:
 \begin{equation}
-A_{j,k} a^\prime_k = \beta_j
+A_{m,n} \delta_n = -g_m
 \end{equation}
 %
-where $a^\prime_k$ represents our new attempt at a parameter guess. We
-use this parameter set to evaluate the function.  If the new value of
-the function is lower than the previous guess, we accept this new set
-of parameters and decrease $\lambda$ by a factor of 10, otherwise we
-keep the old set, and increase the value of $\lambda$ by a factor of
-10.  We repeat this process until the value of the function changes by
-much less than the tolerance.  The resulting values of $a_k$ are the
-best-fit parameters for the system.
+where $alpha^\prime_n = \alpha_n + \delta_n$ represents our new
+attempt at a parameter guess. We use this parameter set to evaluate
+the function.  
+
+To evaluate the new parameter set, we compare the change in $\chi^2$
+with the change expected from the linear model (the Taylor expansion).
+If the linear model were correct, we would have expected a change (a
+reduction) in $\chi^2$ of $\Delta = \frac{\lambda}{2}\sum\delta_m^2 +
+\frac{1}{2}\sum\delta_m g_m$.  We use the 'gain ratio' $\rho =
+\frac{\chi^2_{\rm old} - \chi^2_{\rm new}}{\Delta}$ to judge the new
+step. If $\rho > 0$, we accept this new set of parameters and decrease
+$\lambda$ by a factor of 10, otherwise we keep the old set, and
+increase the value of $\lambda$ by a factor of 10.  We repeat this
+process until the value of the function changes by much less than the
+tolerance.  The resulting values of $a_m$ are the best-fit parameters
+for the system.
 
 The covariance matrix, $C_{i,j}$, which is the inverse of the matrix
-$\alpha_{j,k}$ allows simple calculation of the confidence limits of
-the parameters.
-
+$H_{m,n}$ allows simple calculation of the confidence limits of the
+parameters.
 
 %If the errors are normally distributed, the formal errors on the
@@ -510,5 +580,4 @@
 %\end{tabular}
 %\end{center}
-
 
 \subsubsection{Powell's method}