6

对于这样的模型，我正在尝试执行线性回归：

Y = aX1 + bX2 + c

所以，Y ~ X1 + X2

假设我有以下响应向量：

set.seed(1)
Y <- runif(100, -1.0, 1.0)

以及以下预测变量矩阵：

X1 <- runif(100, 0.4, 1.0)
X2 <- sample(rep(0:1,each=50))
X <- cbind(X1, X2)

我想对系数使用以下约束：

a + c >= 0  
c >= 0

所以对 b 没有限制。

我知道 glmc 包可用于应用约束，但我无法确定如何将它应用于我的约束。例如，我也知道可以使用 contr.sum 使所有系数总和为 0，但这不是我想要做的。solve.QP() 似乎是另一种可能性，meq=0可以使用设置使所有系数都> = 0（同样，这不是我的目标）。

注意：解决方案必须能够处理响应向量 Y 中的 NA 值，例如：

Y <- runif(100, -1.0, 1.0)
Y[c(2,5,17,56,37,56,34,78)] <- NA
4

1 回答 1

4

solve.QP可以传递任意线性约束，因此它当然可以用于建模您的约束a+c >= 0和c >= 0.

首先，我们可以添加一列 1X来捕获截距项，然后我们可以复制标准线性回归solve.QP：

X2 <- cbind(X, 1)
library(quadprog)
solve.QP(t(X2) %*% X2, t(Y) %*% X2, matrix(0, 3, 0), c())$solution
# [1]  0.08614041  0.21433372 -0.13267403

使用来自问题的样本数据，使用标准线性回归没有满足任何约束。

通过修改Amat和bvec参数，我们可以添加两个约束：

solve.QP(t(X2) %*% X2, t(Y) %*% X2, cbind(c(1, 0, 1), c(0, 0, 1)), c(0, 0))$solution
# [1] 0.0000000 0.1422207 0.0000000

在这些约束条件下，通过将 a 和 c 系数都设置为 0 来最小化残差平方。

Y您可以通过删除有问题的观察结果来处理或X2像lm函数一样处理缺失值。作为预处理步骤，您可能会执行以下操作：

has.missing <- rowSums(is.na(cbind(Y, X2))) > 0
Y <- Y[!has.missing]
X2 <- X2[!has.missing,]
于 2017-08-08T21:57:15.397 回答