[R] constrained regression

[Mike Lonergan]

Giles,

While I hesitate to suggest mutilating any of Prof. Ripley's code, one
approach to forcing coefficients to be positive might be to exploit his
step() function after hacking his extractAIC.lm() to heavily penalise
negative coefficients (since such terms are likely to need to become zero/be
dropped rather than turn positive). Something like:

 extractAIC.lm<-function (fit, scale = 0, k = 2, ...)
{
...
    c(edf, (dev + k * edf)*(1+1000000*length(fit$coef[fit$coef<0]))
}

should do the trick. (I guess. I haven't really tested it.) It will, of
course, need to be changed back if you then want to use it properly.

I believe the other constraints can also be converted to positive ones by
rearranging the inequalities and substituting.

eg:

 y= ax + bz + error
 a + b < 1

 is equivalent to

 y =a.(x-z) + z +d.(-z) +error

 where d =1-a-b, and so is positive.

 In this case the lm formula then becomes:

y~I(x-z)+offset(z) +I(-z)

which is a bit of a mess, but better than working for a living. I haven't
really had a need for this, and so have little idea how inconvenient it is
to use. Hopefully wiser heads will warn you off if they consider it all too
daft.

My one serious warning would be about using the resultant models for
prediction. In any cases where strongly negative coefficients are suppressed
this is likely to result in obvious trends in the residuals, and the
theoretical foundations of any confidence intervals round the fitted values
may be badly dented. Whether this matters will depend on what you're
intending to do with the models, but I'd certainly not bet money on the
results.

Cheers,

Mike.

     > -----Original Message-----
     > From: owner-r-help at stat.math.ethz.ch
     > [mailto:owner-r-help at stat.math.ethz.ch]On Behalf Of Heywood, Giles
     > Sent: 28 May 2002 15:12
     > To: r-help at stat.math.ethz.ch
     > Subject: [R] constrained regression
     >
     >
     > I want to do a linear regression where the coefficients
     > obey two linear constraints, and also are all non-negative.
     > What is the best way to do this? Computational speed is a
     > consideration as I must do it many times.
     >
     > When this question was asked previously on the list, quadprog
     > was suggested - is this the best solution?
     >
     > (I may have missed something obvious in the documentation,
     > but I have looked).
     >
     > thanks
     >
     > Giles
     >
     >
     > *************************************************************
     > *********
     > This communication is confidential and is intended only for
     > the person to whom it is addressed.  If you are not that person you
     > are not permitted to make use of the information and you are
     > requested
     > to notify <mailto:LONIB.Postmaster at commerzbankib.com>
     > immediately that
     > you have received it and then destroy the copy in your possession.
     > Commerzbank AG is regulated by the FSA for the conduct of investment
     > business in the UK.
     > *************************************************************
     > *********
     >
     > -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
     > .-.-.-.-.-.-.-.-.-
     > 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
     > _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
     > ._._._._._._._._._

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._