# Multiple linear regression with restricted ranges in R

I'm trying to do a simple least-squares fit on some data. The formula is simply:

$Y = C_0X_1 + C_0X_2 + C_0X_3 + C_0X_4 + C_0X_5$

I have 24 rows of $Y$'s and $X$'s, and I'm trying to fit the $C$'s.

The problem, however, is that I need to constrain the ranges of the $C$'s to be non-zero, regardless of the quality of the fit. I've read the docs on lm in R, but I cannot find the option to constrain the variable ranges. Any ideas?

-
What is the reason to constrain the ranges? If it is just a technical question you should maybe ask on stackoverflow, since this is more of a stats theoretical site. If the reason for constraint is a inheret to your model, I'd advise you to consider a (semi)-logarithmic model depending on what you are doing. Cutting out data point entirely tends to be a really bad idea for inference. – IMA Jan 21 at 8:07
I'm trying to fit some biological data, and need to constrain the coefficients such that they are not nonsensical... – user19869 Jan 21 at 8:17
Oh I completely misunderstood you then. Sorry. May I ask: Is the linear model necessarily appropriate? What results do you get when estimating right now? – IMA Jan 21 at 8:21
no problem....linear model's a start...right now I get negative values for one of the coefficients, which I can't have. – user19869 Jan 21 at 8:26
I fixed up the formatting, but do you really mean to have $C_0$ five times? This would be most unusual; I think you want $C_1$, $C_2$ etc. Also, you didn't include an intercept; was that deliberate? – Peter Flom Jan 21 at 11:48
show 1 more comment

I had the same problem and found two solutions.

1) First option is that you need to obtain the likelihood by hand and then optimize in R by using the optim command. When you use this command, you need to specify starting values for the parameters to be estimated. As far as I know, R uses 0 for all the coefficients. You need to try different initial values. Of course, you will highly unlikely get the positives coefficients if you determine these initial values by hand. My suggestion is to generate random numbers using rnorm (or some other distribution) for the initials values and then use repeat until break command. You need to break the repetition when all of the coefficients are positive. Then, you need to save these numbers so that you no longer need to run the repeat until break command. For more info, see http://www.r-bloggers.com/logistic-regression/

Another solution would be to use a Bayesian approach. If you specify uniform priors whose lower limits are 0 for your coefficients, you will get positive coefficients. You can conduct this analysis in WinBugs in minutes. Here is an example:

model{
# modelâ€™s likelihood
for (i in 1:n){
time[i] ~ dnorm( mu[i], tau ) # stochastic componenent
mu[i] <- beta1 * cases[i] + beta2 * distance[i]
}
# prior distributions
tau ~ dgamma( 0.01, 0.01 )
beta1 ~ dunif( 0, 1000)
beta2 ~ dunif( 0, 1000)
# definition of sigma
s2<-1/tau
s <-sqrt(s2)
}

-

The issue could be not including an intercept, as mentioned above. You can add another $X_i$ and set it to 1 throughout if you need it in variable form.
Not including an intercept can make your results invalid if the mean of $Y$ truly does deviate from 0.

If that is not the case then I strongly suspect misspecification to be the issue here. There are two caveats that are of special interest, but I don't know on what level your knowledge is so I will keep this brief initially since it is pretty basic stuff.

First: Missing variable bias. Any variable with a strong influence (or in this case, correlation) on the dependent variable which is also correlated with the exogenous variable in question can not be left out - otherwise your estimator $C_i$ will be invalid.
This could explain an unexpectedly negative $C_i$ if there is another variable which you left out of the model.
In this case one would suspect a variable that influences $Y$ negatively, yet is correlated with your $C_i$. Its influence will be attributed to $C_i$, maybe pulling it just into negative range.

Second: The underlying effect is not linear.
I would like you to plot the exogenous variables against the endogenous for a first hint. Then please go ahead and run a RESET-Test. The R function is this here
Please don't do this by hand.

This tests your model against several more general Taylor approximations. If the relationship is not linear, this could be a good indicator. The test will generally "fail" your model if those other models do a better job at explaining the variation of $Y$ for whatever reason.

Of course in the end any estimation problem which makes your estimates biased could be the problem. In fact the RESET test may trigger because of other things than just the linearity of the model - but I suspect those two.

Edit: Because there was another approach posted above let me reiterate that you do NOT want to use some quick fix to force your coefficients into beeing positive. All this does is make your inference nonsense.
Using this model, especially using lm, there is a REASON why that coefficient is estimated as negative. If that can not be in reality, then it is almost certain that the whole model is errornous. Switching to a different method or trying to fix the "error" until you reach an acceptable value is very dangerous - statistical inference absolutely mandates that if you get bogus results then you have to take a step back, not sideways.

Also consider that you literally have 19 or 18 degrees of freedom on this regression. To use OLS you HAVE to have the complete set of small-samples assumptions in check, otherwise it comes as no surprise your regression goes astray.

-
 Some negative coeff. estimates? Maybe some of the X's are highly correlated and a multicollinearity problem shows up. Try PCA or exclude some of the highly correlated variables - if any. – Pantera Jan 21 at 22:53