Лассо проходит через LARS (итерационный процесс, который начинается с некоторой начальной оценки ). По умолчанию но вы можете изменить это в большинстве реализаций (и заменить его на оптимальный вас уже есть). Ближайшее - , тем меньше число итераций LARS, которое вам придется пройти, чтобы добраться до . $\beta^0$ $\beta^0=0_p$ $\beta^*_{old}$ $\beta^*_{old}$ $\beta_{new}^*$ $\beta_{new}^*$

РЕДАКТИРОВАТЬ:

Из-за комментариев user2763361я добавляю больше деталей к моему первоначальному ответу.

Из комментариев ниже я понимаю, что пользователь 2763361 предлагает дополнить мой первоначальный ответ, чтобы превратить его в тот, который можно использовать напрямую (с полок), но при этом он очень эффективен.

Чтобы сделать первую часть, я проиллюстрирую решение, которое я предлагаю, шаг за шагом на игрушечном примере. Чтобы удовлетворить вторую часть, я сделаю это с помощью недавнего высококачественного решателя точек интерьера. Это связано с тем, что легче получить высокопроизводительную реализацию решения, которое я предлагаю, используя библиотеку, которая может решить проблему лассо с помощью подхода внутренней точки зрения, а не пытаться взломать LARS или симплексный алгоритм, чтобы начать оптимизацию с не стандартная отправная точка (хотя это второе место также возможно).

Обратите внимание, что иногда (в более старых книгах) утверждается, что внутренний подход к решению линейных программ медленнее, чем симплексный подход, и это, возможно, было верно давно, но в целом это не так сегодня и, конечно, не верно для крупномасштабных задач. (именно поэтому большинство профессиональных библиотек любят cplexиспользовать алгоритм внутренней точки), и вопрос, по крайней мере, неявно касается крупномасштабных проблем. Также обратите внимание, что решатель внутренних точек, который я использую, полностью обрабатывает разреженные матрицы, поэтому я не думаю, что с LARS будет большой разрыв в производительности (оригинальная мотивация для использования LARS заключалась в том, что многие популярные решатели LP в то время плохо обрабатывали разреженные матрицы и это характерные черты проблемы ЛАССО).

(Очень) хорошая реализация алгоритма внутренней точки с открытым исходным кодом находится ipoptв COIN-ORбиблиотеке. Еще одна причина , которую я буду использовать в ipoptтом , что она имеет имеет интерфейс R, ipoptr. Более подробное руководство по установке вы найдете здесь , ниже я даю стандартные команды для его установки ubuntu.

в bash, сделайте:

sudo apt-get install gcc g++ gfortran subversion patch wget
svn co https://projects.coin-or.org/svn/Ipopt/stable/3.11 CoinIpopt
cd ~/CoinIpopt
./configure
make 
make install

Затем, от имени пользователя root, в Rdo (я предполагаю svn, что файл subversion скопирован в файл ~/по умолчанию):

install.packages("~/CoinIpopt/Ipopt/contrib/RInterface",repos=NULL,type="source")

Отсюда я привожу небольшой пример (в основном из игрушечного примера, приведенного Джелмером Ипмой в качестве составной части его Rобертки ipopt):

library('ipoptr')
# Experiment parameters.
lambda <- 1                                # Level of L1 regularization.
n      <- 100                              # Number of training examples.
e      <- 1                                # Std. dev. in noise of outputs.
beta   <- c( 0, 0, 2, -4, 0, 0, -1, 3 )    # "True" regression coefficients.
# Set the random number generator seed.
ranseed <- 7
set.seed( ranseed )
# CREATE DATA SET.
# Generate the input vectors from the standard normal, and generate the
# responses from the regression with some additional noise. The variable 
# "beta" is the set of true regression coefficients.
m     <- length(beta)                           # Number of features.
A     <- matrix( rnorm(n*m), nrow=n, ncol=m )   # The n x m matrix of examples.
noise <- rnorm(n, sd=e)                         # Noise in outputs.
y     <- A %*% beta + noise                     # The outputs.
# DEFINE LASSO FUNCTIONS
# m, lambda, y, A are all defined in the ipoptr_environment
eval_f <- function(x) {
    # separate x in two parts
    w <- x[  1:m ]          # parameters
    u <- x[ (m+1):(2*m) ]

    return( sum( (y - A %*% w)^2 )/2 + lambda*sum(u) )
}
# ------------------------------------------------------------------
eval_grad_f <- function(x) {
    w <- x[ 1:m ]
    return( c( -t(A) %*% (y - A %*% w),  
               rep(lambda,m) ) )
}
# ------------------------------------------------------------------
eval_g <- function(x) {
    # separate x in two parts
    w <- x[  1:m ]          # parameters
    u <- x[ (m+1):(2*m) ]
    return( c( w + u, u - w ) )
}
eval_jac_g <- function(x) {
    # return a vector of 1 and minus 1, since those are the values of the non-zero elements
    return( c( rep( 1, 2*m ), rep( c(-1,1), m ) ) )
}
# ------------------------------------------------------------------
# rename lambda so it doesn't cause confusion with lambda in auxdata
eval_h <- function( x, obj_factor, hessian_lambda ) {
    H <- t(A) %*% A
    H <- unlist( lapply( 1:m, function(i) { H[i,1:i] } ) )
    return( obj_factor * H )
}
eval_h_structure <- c( lapply( 1:m, function(x) { return( c(1:x) ) } ),
                       lapply( 1:m, function(x) { return( c() ) } ) )
# The starting point.
x0 = c( rep(0, m), 
        rep(1, m) )
# The constraint functions are bounded from below by zero.
constraint_lb = rep(   0, 2*m )
constraint_ub = rep( Inf, 2*m )
ipoptr_opts <- list( "jac_d_constant"   = 'yes',
                     "hessian_constant" = 'yes',
                     "mu_strategy"      = 'adaptive',
                     "max_iter"         = 100,
                     "tol"              = 1e-8 )
# Set up the auxiliary data.
auxdata <- new.env()
auxdata$m <- m
    auxdata$A <- A
auxdata$y <- y
    auxdata$lambda <- lambda
# COMPUTE SOLUTION WITH IPOPT.
# Compute the L1-regularized maximum likelihood estimator.
print( ipoptr( x0=x0, 
               eval_f=eval_f, 
               eval_grad_f=eval_grad_f, 
               eval_g=eval_g, 
               eval_jac_g=eval_jac_g,
               eval_jac_g_structure=eval_jac_g_structure,
               constraint_lb=constraint_lb, 
               constraint_ub=constraint_ub,
               eval_h=eval_h,
               eval_h_structure=eval_h_structure,
               opts=ipoptr_opts,
               ipoptr_environment=auxdata ) )

Я хочу сказать, что если у вас есть новые данные, вам просто нужно

обновить ( не заменить) матрицу ограничений и вектор целевой функции для учета новых наблюдений.
изменить начальные точки внутренней точки от

x0 = c (rep (0, m), rep (1, m))

к вектору решения, которое вы нашли ранее (до добавления новых данных). Логика здесь следующая. Обозначим новый вектор коэффициентов (те, которые соответствуют набору данных после обновления) и исходные. Также обозначим вектор в приведенном выше коде (это обычное начало для метода внутренней точки). Тогда идея заключается в том, что если: $\beta_{new}$ $\beta_{old}$ $\beta_{init}$ x0

| β_{i n i t} - β_{n e w} |_{1} > | β_{n e w} - β_{o l d} |_{1} (1)

$|\beta_{init}-\beta_{new}|_1>|\beta_{new}-\beta_{old}|_1\quad(1)$

тогда можно получить намного быстрее, начав внутреннюю точку с а не с наивного . Усиление будет тем более важным, когда размеры набора данных ( и ) больше. $\beta_{new}$ $\beta_{old}$ $\beta_{init}$ $n$ $p$

Что касается условий, при которых выполняется неравенство (1), это:

когда велика по сравнению с (это обычно имеет место, когда , число проектных переменных велико по сравнению с , числом наблюдений) $\lambda$ $|\beta_{OLS}|_1$ $p$ $n$
когда новые наблюдения не являются патологически влиятельными, например, когда они согласуются со случайным процессом, который породил существующие данные.
когда размер обновления небольшой по сравнению с размером существующих данных.

user603
источник

Спасибо, но я боюсь, что не пойду. LARS создает кусочно-линейный путь (с точками для наименьших углов и, возможно, больше точек для лассо.) У каждой точки есть свой набор . Когда мы добавим больше наблюдений, все бета-версии могут двигаться (кроме , который всегда .) Не могли бы вы расширить свой ответ? Благодарю.

p + 1

$p+1$

β

$\beta$

β^{0}

$\beta^0$

0_{p}

$0_p$

NPE

@aix: вы хотите обновить весь путь лассо или просто решение? (т.е. исправлен штраф за разреженность?).

user603

Я искал, чтобы обновить весь путь. Однако, если есть хороший способ сделать это для фиксированного штрафа ( в формуле ниже), это может быть хорошим началом. Это то, что вы предлагаете?

λ

$\lambda$

{\hat{β}}^{l a s s o} = \underset{β}{argmin} {\frac{1}{2} \sum_{i = 1}^{N} (y_{i} - β_{0} - \sum_{j = 1}^{p} x_{i j} β_{j})^{2} + λ \sum_{j = 1}^{p} | β_{j} |}

$\hat{\beta}^{lasso} = \underset{\beta}{\operatorname{argmin}} \left \{ {1 \over 2} \sum_{i=1}^N(y_i-\beta_0-\sum_{j=1}^p x_{ij} \beta_j)^2 + \lambda \sum_{j=1}^p |\beta_j| \right \}$

NPE

@aix. Да, все зависит от реализации, которую вы используете, и возможностей, к которым у вас есть доступ. Например: если у вас есть доступ к хорошему решению для lp, вы можете передать ему последние оптимальные значения и он будет очень эффективно переносить 1-2 шага к новому решению. Вы должны добавить эти детали к вашему вопросу.

β

$\beta$

user603

Какие библиотеки вы знаете, которые могут сделать это без необходимости редактировать исходный код?

user2763361

Обновление лассо соответствует новым наблюдениям

Ответы:

РЕДАКТИРОВАТЬ: