I'm fumbling around trying to figure out how to use minpack (lmstr)
for a weighted non-linear least squares fit and I know it shouldn't be
as tough as I'm making it. For weights, I'm using a full covariance
matrix, S, so in Bx = y the linear problem is:
(K^T S^(-1) K) x = K^T S^(-1) y
where K is the Jacobian matrix. My confusion is that lmstr (and most
similar routines) all expect you to provide:
m, n = rows, cols of K
x = initial guess (overwritten with solution)
y = data to fit
fcn = function to evaluate F(x) and Jacobian matrix (dF/dx)
It isn't obvious to me how to weight the data or write the model
function 'fcn' when weights are used?
Does anyone have a basic fortran example showing how to set this up?
Thanks,
Dave
If S = L*L' (cholesky) and you define A = inv(L)*B and z = inv(L)*y
then the unweighted ls problem with A and z will have the same solutions
as the weighted one.
in the function routine fcn, before exit multiply the function vector and
the Jacobian (if you do it analytically) by inv(L) where
S = L*L^T is the cCholesky decomposition of S
hth
peter