[R] Fit with constraints (fwd)

[Douglas Bates]

Guido Masarotto <guido at hal.stat.unipd.it> writes:

> On Mon, Nov 22, 1999 at 08:50:31PM +0100, Agustin Lobo wrote:
> > 
> > I would like to fit a model of the type
> > 
> > Y = Xf + e
> > 
> > where f are cover fractions (or percent cover). Therefore,
> > I must constrain the fit to
> > 1. sum(f) = 1
> > 2. 0<= f <=1
...
>   Since R doesn't have 'constrained optimization' (yet) what I have tipically
>   done with similar  problems is 'proper' re-parametrization, e.g.,
>   f1 = exp(g1)/s f2=exp(g2)/s f3=1/s with s=1+exp(g1)+exp(g2) then
>   likelihood is maximized in (g1,g2). Of course, if you want
>   asymptotic standard errors, you have to transform the Fisher matrix
>   for (g1,g2) to the one of (f1,f2) using the jacobian of the transformation.
>   Anyway, if possible, I will try to use a profile likelihood approach.
...
> > The problem is that I get very odd f, which
> > are very dependent on what f I choose to be
> > the one to be calculated by difference (the f3 above).
> > 
> 
>  This can be a problem of the algorithm but also of the data,
>  or better, of your model with your data.
>  Hence, if you really have only three f's, I suppose that a
>  graphical study of the likelihood can be usefull.

If you use the nls function in R to fit the model in the
parameterization that Guido described, you can graphically examine the
profile likelihood with 
 pr <- profile(nlsObject)
 plot(pr)

If you have installed R-0.90.0 you can see an example of the types of
plots available with
 library(nls)
 example(plot.profile.nls)

P.S. The main title and subtitle on those plots are in the wrong
 positions.  Either I misunderstood the mtext function documentation or
 there was a bug introduced in the graphics code sometime around
 Nov. 20.  I'll file a bug report if I can confirm the bug.

P.P.S. The plots of the absolute value of tau versus the parameter
 being profiled should use yaxs = "i".  At the time I wrote the
 example I didn't realize that this worked.

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._