R:强制回归系数之和为1

0vvn1miw  于 2023-02-06  发布在  其他
关注(0)|答案(2)|浏览(346)

我试着运行一个简单的OLS回归,其限制条件是两个变量的系数之和等于1。
我想要:

Y = α + β1 * x1 + β2 * x2 + β3 * x3,
where β1 + β2 = 1

我发现如何使系数之间的关系如下:

β1 = 2* β2

但我还没有找到如何作出限制,如:

β1 = 1 - β2

在这个简单的例子中,我应该如何做呢?

data <- data.frame(
  A = c(1,2,3,4),
  B = c(3,2,2,3),
  C = c(3,3,2,3),
  D = c(5,3,3,4)
)

lm(formula = 'D ~ A + B + C', data = data)

谢谢!

euoag5mw

euoag5mw1#

β 1 + β 2 = 1

要获得β1 + β2 = 1,您必须拟合的模型为

fit <- lm(Y ~  offset(x1) + I(x2 - x1) + x3, data = df)

那是
Y = α + x1+ β 2 *(x2-x1)+ β 3 * x3
代之以 * β 1 = 1-β 2 *;x_new = x2 - x1,并且x1的系数为1。

β 1 + β 2 + β 3 = 1

fit <- lm(Y ~  offset(x1) + I(x2 - x1) + I(x3 - x1), data = df)

Y = α + x1+ β 2 *(x2-x1)+ β 3 *(x3-x1)
在替换β1 = 1 - β2 - β3之后

β 1 + β 2 + β 3 +...= 1

我认为模式很清楚......您只需从剩余变量(x2x3...)中减去一个变量x1,并使该变量的系数x1为1。

示例β 1 + β 2 = 1

# Data
df <- iris[, 1:4]
colnames(df) <- c("Y", paste0("x", 1:3, collaapse=""))

# β1 + β2 = 1
fit <- lm(Y ~  offset(x1) + I(x2 - x1) + x3, data = df)
coef_2 <- coef(fit)
beta_1 <- 1 - coef_2[2]
beta_2 <- coef_2[2]
ryhaxcpt

ryhaxcpt2#

    • 1)CVXR**我们可以通过指定目标和约束条件直接使用CVXR计算系数。我们假设D是响应,A和B的系数之和必须为1,b [1]是截距,b [2]、b [3]和b [4]分别是A、B和C的系数。
library(CVXR)

b <- Variable(4)
X <- cbind(1, as.matrix(data[-4]))
obj <- Minimize(sum((data$D - X %*% b)^2))
constraints <- list(b[2] + b[3] == 1)
problem <- Problem(obj, constraints)
soln <- solve(problem)

bval <- soln$getValue(b)
bval
##            [,1]
## [1,]  1.6428605
## [2,] -0.3571428
## [3,]  1.3571428
## [4,] -0.1428588

目标是残差平方和,它等于:

soln$value
## [1] 0.07142857
    • 2)pracma**我们也可以使用pracma软件包来计算系数,我们指定了X矩阵、响应向量、约束矩阵(在本例中,作为第三个参数给出的向量被视为一行矩阵)和约束的右侧。
library(pracma)

lsqlincon(X, data$D, Aeq = c(0, 1, 1, 0), beq = 1) # X is from above
## [1]  1.6428571 -0.3571429  1.3571429 -0.1428571
    • 3)limSolve**这个包也可以求解带约束的回归问题的系数。参数与(2)中的相同。
library(limSolve)
lsei(X, data$D, c(0, 1, 1, 0), 1)

给出:

$X
                    A          B          C 
 1.6428571 -0.3571429  1.3571429 -0.1428571 

$residualNorm
[1] 0

$solutionNorm
[1] 0.07142857

$IsError
[1] FALSE

$type
[1] "lsei"
    • 4)nls**这可以用公式表示为nls的问题,其中B系数等于1减去A系数。
nls(D ~ b0 + b1 * A + (1-b1) * B + b2 * C, data, 
  start = list(b0 = 1, b1 = 1, b2 = 1))
##     D ~ b0 + b1 * A + (1 - b1) * B + b2 * C
##     data: data
##      b0      b1      b2 
##  1.6429 -0.3571 -0.1429 
##  residual sum-of-squares: 0.07143
##
## Number of iterations to convergence: 1 
## Achieved convergence tolerance: 2.803e-08

检查

我们可以使用另一个答案中的lm方法来再次检查上述内容:

lm(D ~ I(A-B) + C + offset(B), data)

给出:

Call:
lm(formula = D ~ I(A - B) + C + offset(B), data = data)

Coefficients:
(Intercept)     I(A - B)            C  
     1.6429      -0.3571      -0.1429

I(A-B)系数等于原始公式中A的系数,减去1就是C的系数,我们看到所有方法都得到相同的系数。

相关问题