Быстрый неотрицательный квантиль и регрессия Хьюбера в R

Question

Я ищу быстрый способ сделать неотрицательную квантиль и регрессию Губер в R (то есть с ограничением, что все коэффициенты>0). Я пытался использовать CVXR пакет для квантильной и Huber регрессии и quantreg пакет для квантильной регрессии, но CVXR очень медленно и quantreg кажется, глючит, когда я использую ограничения неотрицательности. Кто-нибудь знает хорошее и быстрое решение в R, например, используя Rcplex пакет или R gurobi API, тем самым используя более быстрые оптимизаторы CPLEX или gurobi?

Обратите внимание, что мне нужно запустить размер проблемы, например, ниже 80 000 раз, при этом мне нужно только обновить y вектор в каждой итерации, но по-прежнему использовать ту же матрицу предиктора X, В этом смысле я чувствую, что это неэффективно, что в CVXR Теперь я должен сделать obj <- sum(quant_loss(y - X %*% beta, tau=0.01)); prob <- Problem(Minimize(obj), constraints = list(beta >= 0)) в каждой итерации, когда проблема остается неизменной, и все, что я хочу обновить, это y, Любые мысли, чтобы сделать все это лучше / быстрее?

Минимальный пример:

## Generate problem data
n <- 7 # n predictor vars
m <- 518 # n cases 
set.seed(1289)
beta_true <- 5 * matrix(stats::rnorm(n), nrow = n)+20
X <- matrix(stats::rnorm(m * n), nrow = m, ncol = n)
y_true <- X %*% beta_true
eps <- matrix(stats::rnorm(m), nrow = m)
y <- y_true + eps

Неотрицательная квантильная регрессия с использованием CVXR:

## Solve nonnegative quantile regression problem using CVX
require(CVXR)
beta <- Variable(n)
quant_loss <- function(u, tau) { 0.5*abs(u) + (tau - 0.5)*u }
obj <- sum(quant_loss(y - X %*% beta, tau=0.01))
prob <- Problem(Minimize(obj), constraints = list(beta >= 0))
system.time(beta_cvx <- pmax(solve(prob, solver="SCS")$getValue(beta), 0)) # estimated coefficients, note that they ocasionally can go - though and I had to clip at 0
# 0.47s
cor(beta_true,beta_cvx) # correlation=0.99985, OK but very slow

Синтаксис неотрицательной регрессии Губера такой же, но будет использоваться

M <- 1      ## Huber threshold
obj <- sum(CVXR::huber(y - X %*% beta, M))

Неотрицательная квантильная регрессия с использованием quantreg пакет:

### Solve nonnegative quantile regression problem using quantreg package with method="fnc"
require(quantreg)
R <- rbind(diag(n),-diag(n))
r <- c(rep(0,n),-rep(1E10,n)) # specify bounds of coefficients, I want them to be nonnegative, and 1E10 should ideally be Inf
system.time(beta_rq <- coef(rq(y~0+X, R=R, r=r, tau=0.5, method="fnc"))) # estimated coefficients
# 0.12s
cor(beta_true,beta_rq) # correlation=-0.477, no good, and even worse with tau=0.01...

Answer 1

Чтобы ускорить CVXR, вы можете получить данные о проблеме один раз в начале, затем изменить их в цикле и передать их непосредственно в интерфейс R решателя. Код для этого

prob_data <- get_problem_data(prob, solver = "SCS")

Затем разберите аргументы и передайте их scs из библиотеки scs. (См. Solver.solve в solver.R). Вам придется углубиться в детали канонизации, но я ожидаю, что если вы просто меняете y на каждой итерации, это должно быть простым изменением.