Because of the constraint that you impose on the coefficients, the model you specify is not a linear model and so lm
can not be used to fit it. You would need to use a non-linear regression, such as nls
.
> summary(nls(y ~ b0 + b1*x1 + b2*x2 + b1*b2*x1*x2, start=list(b0=0, b1=1, b2=1)))
Formula: y ~ b0 + b1 * x1 + b2 * x2 + b1 * b2 * x1 * x2
Parameters:
Estimate Std. Error t value Pr(>|t|)
b0 0.987203 0.049713 19.86 <2e-16 ***
b1 0.494438 0.007803 63.37 <2e-16 ***
b2 0.202396 0.003359 60.25 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1121 on 47 degrees of freedom
Number of iterations to convergence: 5
Achieved convergence tolerance: 2.545e-06
You can really see that the model is non-linear when you re-write it as
> summary(nls(y ~ b0+(1+b1*x1)*(1+b2*x2)-1, start=list(b0=0, b1=1, b2=1)))
Formula: y ~ b0 + (1 + b1 * x1) * (1 + b2 * x2) - 1
Parameters:
Estimate Std. Error t value Pr(>|t|)
b0 0.987203 0.049713 19.86 <2e-16 ***
b1 0.494438 0.007803 63.37 <2e-16 ***
b2 0.202396 0.003359 60.25 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1121 on 47 degrees of freedom
Number of iterations to convergence: 5
Achieved convergence tolerance: 2.25e-06