[R] Fitting generalized additive models with constraints?

[David Reiss]

>> I am trying to fit a GAM for a simple model, a simple model, y ~ s(x0) +
>> s(x1) ; with a constraint that the fitted smooth functions s(x0) and s(x1)
>> have to each always be >0.
>>
>> From the library documentation and a search of the R-site and R-help
>> archives I have not been able to decipher whether the following is possible
>> using this, or other GAM libraries, or whether I will have to try to "roll
>> my own". I see from the mgcv docs that GAMs need to be constrained such
>> that the smooth functions have zero mean. Is there a way around this?
>>
>> Is such a constraint possible?

> It is possible to estimate a GAM subject to this constraint, but be aware that
> the mean levels of your component smooths are not identifiable, so there is
> an unavoidable abitrariness in the estimate....

> You have to have some sort of constraint on the smooths in a GAM to ensure
> identifiability, and a convenient way to set the model up is to write it as
> e.g.
>
> E(y) = a + f0(x0) + f1(x1)
>
> where `a' is the intercept and f0 and f1 are smooth functions which sum to
> zero over their respective covariate values. In this parameterization your
> constraint implies that
>
> a + f0(x0) + f1(x1)  > 0
>
> for all x0, x1. If this constraint is met then you can find constants b and c
> such that b+c=a such that f0(x0)+b>0 and f1(x1)+c>0 for all x0,x1. i.e. you
> redefine f0 as f0+b and f1 as f1+c, and you have a fitted model meeting the
> constraints.
>
> To fit the GAM subject to the constraints you can use mgcv:::pcls... ?pcls has
> some examples, but it does involve moderately low level programming. It's
> hard to impose the constraint exactly, so the usual approach would be to
> impose the constraint over a fairly fine grid of x0, x1 values. Also, you'll
> need to figure out how to select smoothing parameters. For many problems it
> suffices to estimate smoothing parameters on the unconstrained fit, and then
> use these to fit subject to constraints, but it depends on the problem....
>
> Hope that's some use.
>
> Simon

Hi Simon,
thanks very much for the advice. I will try to parse your response and
the pcls docs and see if I can get this to work. In the meantime, I
found a paper that tries to achieve a similar thing with the same
constraints as I am working with, using quadprog:
http://www.esajournals.org/esaonline/?request=get-abstract&issn=0012-9658&volume=083&issue=08&page=2256
-David