Solving a NLLSQ Problem by Pieces?
The context here is astronomy and optics. The real point is in the last
paragraph.
I'm looking at a paper that deals with 5 NL (nonlinear) equations and 8
unknown parameters.
A. a=a0+arctan((y-y0)/(x-x0)
B. z=V*r+S*e**(D*r)
r=sqrt((x-x0)**2+(y-y0)**2)
and
C. cos(z)=cos(u)*cos(z)-sin(u)*sin(ep)*cos(b)
sin(a-E) = sin(b)*sin(u)/sin(z)
He's trying to estimate parameters of a fisheye lens which has taken
star images on a photo plate. For example, x0,y0 is the offset of the
center of projection from the zenith (camera not pointing straight up in
the sky.) Eq. 2 expresses some nonlinearity in the lens.
a0, xo, y0, V, S, D, ep, and E are the parameters. It looks like he uses
gradient descent (NLLSQ is nonlinear least squares in Subject.), and
takes each equation in turn using the parameter values from the
preceding one in the next, B. He provides reasonable initial estimates.
A final step uses all eight parameters. He re-examines ep and E, and
assigns new estimates. For all (star positions) on the photo plate, he
minimizes SUM (Fi**2*Gi) using values from the step for A and B, except
for x0,y0. He then does some more dithering, which I'll skip.
What I've presented is probably a bit difficult to understand without a
good optics understanding, but my question is something like is this
commonly done to solve a system of NLLSQ? It looks a bit wild. I guess
if one knows his subject well, then bringing some "extra" knowledge to
the process helps. As I understand it, he solves parameters in A, then
uses them in B, and so on. I guess that's a reasonable way to do it.
--
Wayne Watson (Watson Adventures, Prop., Nevada City, CA)
(121.015 Deg. W, 39.262 Deg. N) GMT-8 hr std. time)
Obz Site: 39° 15' 7" N, 121° 2' 32" W, 2700 feet
Air travel safety Plus Three/Minus 8 rule. Eighty %
of crashes take place in the first 3 minutes and
last 8 minutes. Survival chances are good in you're
paying attention. No hard drinks prior to those
periods. No sleeping pills either. Sit within 5 rows
of an exit door. --The Survivors Club, Ben Sherwood
Web Page:
The basic problem with nonlinear least squares fitting, as with other
nonlinear minimization problems, is that the standard algorithms find
only a local minimum. It's easy to miss the global minimum and instead
settle on a local minimum that is in fact a horrible fit.
To deal with this, there are global optimization techniques -
simulated annealing, genetic algorithms, and what have you (including
the simplest, explore the whole space with a "dense enough" grid then
fine-tune the best one with a local optimizer) - but they're
computationally very expensive and not always reliable. So when it's
possible to use domain-specific knowledge to make sure what you're
settling on is the real minimum, this is a better solution.
In this specific case, as in many optimization problems, the problem
can be viewed as a series of nested approximations. The crudest
approximation might even be linear in some cases; but the point is,
you solve the first approximation first because it has fewer
parameters and a solution space you understand better (e.g. maybe you
can be sure it only has one local minimum). Then because your
approximations are by assumption nested, adding more parameters
complicates the space you're solving over, but you can be reasonably
confident that you're "close" to the right solution already. (If your
approximations are "orthogonal" in the right sense, you can even fix
the parameters from previous stages and optimize only the new ones; be
careful with this, though.)
This approach is a very good way to incorporate domain-specific
knowledge into your code, but you need to parameterize your problem as
a series of nested approximations, and if it turns out the corrections
are not small you can still get yourself into trouble. (Or, for that
matter, if the initial solution space is complex enough you can still
get yourself in trouble. Or if you're not careful your solver can take
your sensible initial guess at some stage and zip off into never-never
land instead of exploring "nearby".)
If you're interested in how other people solve this particular
problem, you could take a look at the open-source panorama stitcher
"hugin", which fits for a potentially very large number of parameters,
including a camera model.
To bring this back nearer to on-topic, you will naturally not find
domain-specific knowledge built into scipy or numpy, but you will find
various local and global optimizers, some of which are specialized for
the case of least-squares. So if you wanted to implement this sort of
thing with scipy, you could write the domain-specific code yourself
and simply call into one of scipy's optimizers. You could also look at
OpenOpt, a scikit containing a number of global optimizers.
Anne
P.S. This question would be better suited to scipy-user or astropy
than numpy-discussion. -A
On 26 June 2010 13:12, Wayne Watson
The context here is astronomy and optics. The real point is in the last paragraph.
I'm looking at a paper that deals with 5 NL (nonlinear) equations and 8 unknown parameters. A. a=a0+arctan((y-y0)/(x-x0) B. z=V*r+S*e**(D*r) r=sqrt((x-x0)**2+(y-y0)**2) and C. cos(z)=cos(u)*cos(z)-sin(u)*sin(ep)*cos(b) sin(a-E) = sin(b)*sin(u)/sin(z)
He's trying to estimate parameters of a fisheye lens which has taken star images on a photo plate. For example, x0,y0 is the offset of the center of projection from the zenith (camera not pointing straight up in the sky.) Eq. 2 expresses some nonlinearity in the lens.
a0, xo, y0, V, S, D, ep, and E are the parameters. It looks like he uses gradient descent (NLLSQ is nonlinear least squares in Subject.), and takes each equation in turn using the parameter values from the preceding one in the next, B. He provides reasonable initial estimates.
A final step uses all eight parameters. He re-examines ep and E, and assigns new estimates. For all (star positions) on the photo plate, he minimizes SUM (Fi**2*Gi) using values from the step for A and B, except for x0,y0. He then does some more dithering, which I'll skip.
What I've presented is probably a bit difficult to understand without a good optics understanding, but my question is something like is this commonly done to solve a system of NLLSQ? It looks a bit wild. I guess if one knows his subject well, then bringing some "extra" knowledge to the process helps. As I understand it, he solves parameters in A, then uses them in B, and so on. I guess that's a reasonable way to do it.
-- Wayne Watson (Watson Adventures, Prop., Nevada City, CA)
(121.015 Deg. W, 39.262 Deg. N) GMT-8 hr std. time) Obz Site: 39° 15' 7" N, 121° 2' 32" W, 2700 feet
Air travel safety Plus Three/Minus 8 rule. Eighty % of crashes take place in the first 3 minutes and last 8 minutes. Survival chances are good in you're paying attention. No hard drinks prior to those periods. No sleeping pills either. Sit within 5 rows of an exit door. --The Survivors Club, Ben Sherwood
Web Page:
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Thanks for the info. Where do I find hugin? The paper I'm looking is from 1987. It is not clear what actual method he used. Probably something popular of those times. Things have changed. :-) Yes, scipy might be a better place to post this, and I will. Hmmm, I think I have had problems finding a connection to a scipy mail list. Astropy is lightly traveled, but I'll give it a shot. On 6/26/2010 10:32 AM, Anne Archibald wrote:
The basic problem with nonlinear least squares fitting, as with other nonlinear minimization problems, is that the standard algorithms find only a local minimum. It's easy to miss the global minimum and instead settle on a local minimum that is in fact a horrible fit.
To deal with this, there are global optimization techniques - simulated annealing, genetic algorithms, and what have you (including the simplest, explore the whole space with a "dense enough" grid then fine-tune the best one with a local optimizer) - but they're computationally very expensive and not always reliable. So when it's possible to use domain-specific knowledge to make sure what you're settling on is the real minimum, this is a better solution.
In this specific case, as in many optimization problems, the problem can be viewed as a series of nested approximations. The crudest approximation might even be linear in some cases; but the point is, you solve the first approximation first because it has fewer parameters and a solution space you understand better (e.g. maybe you can be sure it only has one local minimum). Then because your approximations are by assumption nested, adding more parameters complicates the space you're solving over, but you can be reasonably confident that you're "close" to the right solution already. (If your approximations are "orthogonal" in the right sense, you can even fix the parameters from previous stages and optimize only the new ones; be careful with this, though.)
This approach is a very good way to incorporate domain-specific knowledge into your code, but you need to parameterize your problem as a series of nested approximations, and if it turns out the corrections are not small you can still get yourself into trouble. (Or, for that matter, if the initial solution space is complex enough you can still get yourself in trouble. Or if you're not careful your solver can take your sensible initial guess at some stage and zip off into never-never land instead of exploring "nearby".)
If you're interested in how other people solve this particular problem, you could take a look at the open-source panorama stitcher "hugin", which fits for a potentially very large number of parameters, including a camera model.
To bring this back nearer to on-topic, you will naturally not find domain-specific knowledge built into scipy or numpy, but you will find various local and global optimizers, some of which are specialized for the case of least-squares. So if you wanted to implement this sort of thing with scipy, you could write the domain-specific code yourself and simply call into one of scipy's optimizers. You could also look at OpenOpt, a scikit containing a number of global optimizers.
Anne P.S. This question would be better suited to scipy-user or astropy than numpy-discussion. -A
On 26 June 2010 13:12, Wayne Watson
wrote: The context here is astronomy and optics. The real point is in the last paragraph.
I'm looking at a paper that deals with 5 NL (nonlinear) equations and 8 unknown parameters. A. a=a0+arctan((y-y0)/(x-x0) B. z=V*r+S*e**(D*r) r=sqrt((x-x0)**2+(y-y0)**2) and C. cos(z)=cos(u)*cos(z)-sin(u)*sin(ep)*cos(b) sin(a-E) = sin(b)*sin(u)/sin(z)
He's trying to estimate parameters of a fisheye lens which has taken star images on a photo plate. For example, x0,y0 is the offset of the center of projection from the zenith (camera not pointing straight up in the sky.) Eq. 2 expresses some nonlinearity in the lens.
a0, xo, y0, V, S, D, ep, and E are the parameters. It looks like he uses gradient descent (NLLSQ is nonlinear least squares in Subject.), and takes each equation in turn using the parameter values from the preceding one in the next, B. He provides reasonable initial estimates.
A final step uses all eight parameters. He re-examines ep and E, and assigns new estimates. For all (star positions) on the photo plate, he minimizes SUM (Fi**2*Gi) using values from the step for A and B, except for x0,y0. He then does some more dithering, which I'll skip.
What I've presented is probably a bit difficult to understand without a good optics understanding, but my question is something like is this commonly done to solve a system of NLLSQ? It looks a bit wild. I guess if one knows his subject well, then bringing some "extra" knowledge to the process helps. As I understand it, he solves parameters in A, then uses them in B, and so on. I guess that's a reasonable way to do it.
-- Wayne Watson (Watson Adventures, Prop., Nevada City, CA)
(121.015 Deg. W, 39.262 Deg. N) GMT-8 hr std. time) Obz Site: 39° 15' 7" N, 121° 2' 32" W, 2700 feet
Air travel safety Plus Three/Minus 8 rule. Eighty % of crashes take place in the first 3 minutes and last 8 minutes. Survival chances are good in you're paying attention. No hard drinks prior to those periods. No sleeping pills either. Sit within 5 rows of an exit door. --The Survivors Club, Ben Sherwood
Web Page:
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
--
Wayne Watson (Watson Adventures, Prop., Nevada City, CA)
(121.015 Deg. W, 39.262 Deg. N) GMT-8 hr std. time)
Obz Site: 39° 15' 7" N, 121° 2' 32" W, 2700 feet
Air travel safety Plus Three/Minus 8 rule. Eighty %
of crashes take place in the first 3 minutes and
last 8 minutes. Survival chances are good in you're
paying attention. No hard drinks prior to those
periods. No sleeping pills either. Sit within 5 rows
of an exit door. --The Survivors Club, Ben Sherwood
Web Page:
Found hugin (google), but I think it's a bit too general for my current
interests. You might enjoy an optimization video course from Stanford
Univ on the web.
http://see.stanford.edu/see/lecturelist.aspx?coll=17005383-19c6-49ed-9497-2b...
--
Wayne Watson (Watson Adventures, Prop., Nevada City, CA)
(121.015 Deg. W, 39.262 Deg. N) GMT-8 hr std. time)
Obz Site: 39° 15' 7" N, 121° 2' 32" W, 2700 feet
Air travel safety Plus Three/Minus 8 rule. Eighty %
of crashes take place in the first 3 minutes and
last 8 minutes. Survival chances are good in you're
paying attention. No hard drinks prior to those
periods. No sleeping pills either. Sit within 5 rows
of an exit door. --The Survivors Club, Ben Sherwood
Web Page:
participants (2)
-
Anne Archibald
-
Wayne Watson