How to constrain regression coefficients to be proportional

Question

Given the regression model $Y = \beta_1 X_1 + \beta_2 X_2 + \epsilon$ I would like to constrain $\beta_2$ to be proportional to $\beta_1$, that is $\beta_2 = \theta\cdot\beta_1$.

$Y = \beta_1 X_1 + \theta\cdot\beta_1 X_2 + \epsilon$

where $\beta_i$ are vectors and $X_i$ are matrices.

How can I get an estimate of $\beta_1$ and $\theta$?

Is $\theta>0$ ? Otherwise, you get $\theta$ by fitting regression per usual, and then $\theta=\frac{\beta_2}{\beta_1}$. — John Madden, Nov 03 '22 at 14:03
Do you mean proportional to a known constant $\theta$? Why not fit an unconstrained regression, divide the two and say, there! I have $\theta$. — AdamO, Nov 03 '22 at 14:03
$\theta$ can take any real number. The motivation comes from https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8984094/ where one set of coefficients (of a treatment) is a proportional constant of another (control). Some inference comes from seeing if $\theta$ is positive or negative. — Reylined, Nov 03 '22 at 14:06
To clarify, $X_1$ and $X_2$ are matrices of the same shape. The betas are a vector of corresponding coefficients. Therefore, I cannot simply divide the two. — Reylined, Nov 03 '22 at 14:16
Please edit your post to include new information in comments in the post! We want posts to be self-contained, motivated, and finally, comments might get deleted! — kjetil b halvorsen, Nov 04 '22 at 13:10
Note that the second problem is more difficult than the first, because one is linear and the other is not. — Federico Poloni, Nov 04 '22 at 13:13
See the very related Linear model with constraints on coefficients in terms of ratios — kjetil b halvorsen, Nov 04 '22 at 13:18

score 10 · Accepted Answer · answered Nov 03 '22 at 17:32

Why not just maximize the likelihood of

$$ Y = \vec{\beta_1} \mathbf{X}_1 + c \vec{\beta_1} \mathbf{X}_2 + \epsilon$$

set.seed(123)
n <- 100
p <- 3
x1 <- matrix(rnorm(n*p), n, p)
x2 <- matrix(rnorm(n*p), n, p)
b1 <- matrix(rnorm(p), p, 1)
theta <- 2
b2 <- theta*b1
y <- x1%%b1 + x2%%b2 + rnorm(n, 0, 3)
out <- nlm(function(b) {
  yhat <- x1%%b[1:3] + b[4]x2%*%b[1:3]
  err <- y-yhat
  -sum(dnorm(y, mean = yhat, sd = sd(err), log=T))
  }, hessian = T, p=rep(0,4))

gives

> out
$minimum
[1] 254.4209
$estimate
[1]  1.11534411  0.00171469 -0.06829704  2.40625637
$gradient
[1] -4.579198e-05 -5.371703e-06  1.440981e-05 -2.793441e-05
$hessian
          [,1]       [,2]      [,3]       [,4]
[1,] 71.713305 -4.9055719 -8.782636 29.8421070
[2,] -4.905572 62.8017176 13.926606  0.9367334
[3,] -8.782636 13.9266064 63.406844 -5.1652359
[4,] 29.842107  0.9367334 -5.165236 14.1889287
$code
[1] 1
$iterations
[1] 19

How to constrain regression coefficients to be proportional

1 Answers1