I am working on a multiple linear regression problem where I would like to constrain only some of the parameters to non-negative values. There have been discussions of how to solve for the parameters posted on SE here and on another site. Both solutions used optim() or constrOptim() to minimize the residual sum of squares, which worked very well for me and gives the same values as lm() for the unconstrained version of the problem.

But I would like now to take this a step further and estimate standard errors of those parameters, similar to what I'd get if I used lm(). However, there is nothing in the optim() object that would suggest a route to estimating errors, and I haven't found anything in the control settings that would suggest a route either. So I'm a little bit at a loss as to how to proceed. My hunch is that this approach - an optimization problem that assumes the data values are fixed and the parameter solutions represent a global minimum - does not allow for error. Is there any validity to that, or am I just missing something basic?

This is a small reproducible example, adapted from the example provided in the first link to work with the optim() function, rather than constrOptim().

    min.RSS <- function(data, par){
      with(data, sum((par[1]*x1 + par[2]*x2 - y)^2))
    dat = data.frame(x1=c(1,2), x2=c(2,3), y=c(5,6))
    result = optim(par=c(0,1), min.RSS, data=dat, method="L-BFGS-B", lower=c(-Inf,0), upper=c(0,Inf))

Thanks for any guidance you can provide.


I am assuming that you have really good reasons for constraining the coefficients and have thought hard about the potential problems with doing so. Note that a linear model is almost always an approximation to reality, and as a result coefficients that "ought to be non-negative" may end up negative in a linear model based on data with correlated predictors. Forcing coefficients to be non-negative poses a serious risk of generating a suboptimal model.

One way to approach your problem would be to repeat your analysis on multiple bootstrap samples from your data set, and use the distributions of the coefficients to provide the desired estimates. The boot package in R is one tool for this approach, although it might take some effort to write the needed function properly.

If you look carefully at the results of your bootstrap analyses, you will almost certainly find that the mean values, over those bootstrapped samples, of coefficients you constrained will differ from their values in your original analysis. That's because constraining non-negativity on a coefficient that might be negative will probably introduce bias in the technical statistical sense. You will need to consider how you will deal with that eventuality, based on your understanding of the subject matter.

  • 1
    $\begingroup$ EdM, thanks for your response. Yes, I have good theoretical reasons for constraining the coefficients. But you're right that I have correlation in the predictors, and so the usefulness of the coefficients for their individual predictive value is minimal. But I think you're probably right that the bootstrapping approach is best. $\endgroup$ – phalteman Jun 29 '15 at 12:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.