Re: constraining parameters in multi-Gaussian 1D fitting [message #45360] |
Mon, 05 September 2005 19:35 |
Craig Markwardt
Messages: 1869 Registered: November 1996
|
Senior Member |
|
|
JD Smith <jdsmith@as.arizona.edu> writes:
> On Mon, 05 Sep 2005 12:07:42 -0500, Craig Markwardt wrote:
>
>>
>> "Jess" <jobrien@mso.anu.edu.au> writes:
>>> One constraint I am unable yet to do is: I = would like to be able to
>>> tie the peak flux of the Gaussians such that the peak flux of last
>>> Gaussian is always greater than that of the first Gaussian.
>>> I tried using
>>> parinfo((n_gauss-1)*3).tied = 'GT P[0]'
>>>
>>> However the tied structure of parinfo doesn't seem to be meant to
>>> accept operators like GT,LT, etc. ...
>>
>> True. MPFIT's TIED fields are limited to equality constraints only.
>
> What if you availed yourself of the ITERPROC procedure to enforce the
> constraint, dragging the fitter (kicking and screaming if necessary)
> back into line if it attempts to step out? Any reason this wouldn't
> work?
It might work, it might not. I suspect that in general the fitter
might get stuck. For example, if we *started* out by thinking the
tallest peak was on the left part of the curve - and made a model
function to match that - but the truth is that a peak on the right is
truly he tallest, then we can enforce all the constraints in the world
and won't come out with the best fit.
Now, with the additional info that the original poster provided, this
will probably not be the case, so s/he might be alright.
Craig
--
------------------------------------------------------------ --------------
Craig B. Markwardt, Ph.D. EMAIL: craigmnet@REMOVEcow.physics.wisc.edu
Astrophysics, IDL, Finance, Derivatives | Remove "net" for better response
------------------------------------------------------------ --------------
|
|
|
Re: constraining parameters in multi-Gaussian 1D fitting [message #45361 is a reply to message #45360] |
Mon, 05 September 2005 18:58  |
Craig Markwardt
Messages: 1869 Registered: November 1996
|
Senior Member |
|
|
"Jess" <jobrien@mso.anu.edu.au> writes:
> Hi Craig (and everyone else),
>
> Thanks for clarifying the use of "tied". In that case I will have to
> use limits. Perhaps I can re-fit a velocity profile altering my parinfo
> limits, if the first multi-Gaussian fit of that profile doesnt yield
> the right relationship between the peak fluxes of each Gaussian.
>
> Unfortunately I can't sequentially fit single Gaussians, as the
> velcoity profile is a sum of overlapping Gaussians. Each profile
> comprises emission from a line-of-sight through a edge-on galaxy disk
> in cylindrical rotation V(R), with gas distribution assumed to be
> azimuthally symmetric, monotonicly decreasing radial flux profile,
> isothermal and with isotropic velocity dispersion.
>
... extensive description deleted ...
>
> What do you think Craig? Is the solution too degenerate to be solved
> with curvefit algorithms like mpfit, despite the heavy constraints I am
> trying to place on most of the parameters in the multi-gaussian curve
> fit?
Sorry, I got a bit lost in that description.
If possible, I would try to fit the "new" peak while holding the other
peaks absolutely constant, and then only after the "new" peak fit has
converged, allow more of the parameters to vary. You want to be as
close as possible to the optimal solution before you allow the most
parameters to vary.
Happy fitting,
Craig
--
------------------------------------------------------------ --------------
Craig B. Markwardt, Ph.D. EMAIL: craigmnet@REMOVEcow.physics.wisc.edu
Astrophysics, IDL, Finance, Derivatives | Remove "net" for better response
------------------------------------------------------------ --------------
|
|
|
Re: constraining parameters in multi-Gaussian 1D fitting [message #45363 is a reply to message #45361] |
Mon, 05 September 2005 17:33  |
JD Smith
Messages: 850 Registered: December 1999
|
Senior Member |
|
|
On Mon, 05 Sep 2005 12:07:42 -0500, Craig Markwardt wrote:
>
> "Jess" <jobrien@mso.anu.edu.au> writes:
>> One constraint I am unable yet to do is: I = would like to be able to
>> tie the peak flux of the Gaussians such that the peak flux of last
>> Gaussian is always greater than that of the first Gaussian.
>> I tried using
>> parinfo((n_gauss-1)*3).tied = 'GT P[0]'
>>
>> However the tied structure of parinfo doesn't seem to be meant to
>> accept operators like GT,LT, etc. ...
>
> True. MPFIT's TIED fields are limited to equality constraints only.
What if you availed yourself of the ITERPROC procedure to enforce the
constraint, dragging the fitter (kicking and screaming if necessary)
back into line if it attempts to step out? Any reason this wouldn't
work?
JD
|
|
|
Re: constraining parameters in multi-Gaussian 1D fitting [message #45364 is a reply to message #45363] |
Mon, 05 September 2005 17:02  |
Jess
Messages: 11 Registered: June 2005
|
Junior Member |
|
|
Hi Craig (and everyone else),
Thanks for clarifying the use of "tied". In that case I will have to
use limits. Perhaps I can re-fit a velocity profile altering my parinfo
limits, if the first multi-Gaussian fit of that profile doesnt yield
the right relationship between the peak fluxes of each Gaussian.
Unfortunately I can't sequentially fit single Gaussians, as the
velcoity profile is a sum of overlapping Gaussians. Each profile
comprises emission from a line-of-sight through a edge-on galaxy disk
in cylindrical rotation V(R), with gas distribution assumed to be
azimuthally symmetric, monotonicly decreasing radial flux profile,
isothermal and with isotropic velocity dispersion.
Thus I assume to know the following constraints:
- number of gaussians in any profile
- "good" parameter estimates of all but one gaussian in any profile.
- The centroid of the unknown gaussian must be higher than all the
others.
- The peak flux of the unknown gaussian should be equal or higher than
all the others.
- For the other Gaussians in a profile, each comprising flux from a
radial bin higher than R at the line of nodes, the parameters held
fixed or strongly limited from having already fitted flux at those
radii when analysing velocity profile of sightlines further out. I
start at the outermost sightline, where I am only fittin gas at Rmax,
then fit the next innermost using the kinematics at the higher R
sightline as parameter estimates.
- In practise, the centroids of higher R gaussians can be projected to
the new sightline and held fixed as they are found accurately in the
outer sightline, while the peak flux and velocity dispersion are
vulnerable to noise and biased by the telescope spatial beam at the
outermost sightlines, so must be bounded for the multi-gaussian fit of
the next few sightlines, until they stabilise.
I will try to generate reasonable bounding limits to supply for these
high R gaussians, such that the output obeys reasonable assumptions of
the input galaxy:
- that radial flux profile is flat or decreasing,
- the centroids of the gaussians at lower radii in line-of-sight
profile is always higher than the Gaussian centroids from gas at higher
radii projected onto the line of sight.
If this works I should find that the 3 parameter fits of higher radius
Gaussian gas components in a sightline profile should stabilise after a
few sightline profiles, as I fit that gas component in successive
sightlines. Thus the multi-gaussian fit to each profile should have few
unknowns, 3 due to the unknown gaussian at the line-of-nodes radius,
the dispersion and peak flux of the next couple of gaussians, say 3-4.
So I never have more than 8-11 unknowns, though less in the outermost
profiles with fewer gaussians. One advantage of an exponential
decreasing radial flux profile is that the poor fits to the highest r
components, while they affect the velocity profile fit of the next
innermost sightline, they become unimportant, providing its fitting the
line-of-node gaussian well.
What do you think Craig? Is the solution too degenerate to be solved
with curvefit algorithms like mpfit, despite the heavy constraints I am
trying to place on most of the parameters in the multi-gaussian curve
fit?
Thanks for your advice,
Jess
|
|
|
Re: constraining parameters in multi-Gaussian 1D fitting [message #45365 is a reply to message #45364] |
Mon, 05 September 2005 10:07  |
Craig Markwardt
Messages: 1869 Registered: November 1996
|
Senior Member |
|
|
"Jess" <jobrien@mso.anu.edu.au> writes:
> One constraint I am unable yet to do is: I = would like to be able to
> tie the peak flux of the Gaussians such that the peak flux of last
> Gaussian is always greater than that of the first Gaussian.
> I tried using
> parinfo((n_gauss-1)*3).tied = 'GT P[0]'
>
> However the tied structure of parinfo doesn't seem to be meant to
> accept operators like GT,LT, etc. ...
True. MPFIT's TIED fields are limited to equality constraints only.
...
> However this requires assigning rather tight bounds to P[0] which I
> really dont know well. Is there a smarter way I can do this using
> 'tied' or another structure in parinfo?
One approach is to fit one gaussian at a time starting with the
strongest one and working your way down to the fainter ones. This
assumes that the peaks are well enough separated that they can be fit
in succession. If the peaks are blended, well, that's a tough
situation.
Long ago I had somebody trying to do something similar, i.e. fit an
arbitrary number of gaussians with arbitrary centroids, widths and
amplitudes. I warned him that I thought it would turn into an
unconstrained mess, it did, and I think he eventually gave up. Sorry,
but least squares fitting is not mind reading. The more external
constraints you can apply, based on knowledge of the problem, the
better.
Craig
--
------------------------------------------------------------ --------------
Craig B. Markwardt, Ph.D. EMAIL: craigmnet@REMOVEcow.physics.wisc.edu
Astrophysics, IDL, Finance, Derivatives | Remove "net" for better response
------------------------------------------------------------ --------------
|
|
|