Каковы различия между регрессией Риджа с использованием R glmnet и Python scikit-learn?

11

Я изучаю раздел LAB §6.6, посвященный регрессии Риджа / Лассо, в книге Джеймса Виттена «Hastie», Tibshirani (2013) «Введение в статистическое обучение с приложениями в R» .

Более конкретно, я пытаюсь применить модель scikit-learn Ridgeк набору данных 'Hitters' из пакета R 'ISLR'. Я создал такой же набор функций, как показано в коде R. Тем не менее, я не могу приблизиться к результатам glmnet()модели. Я выбрал один параметр настройки L2 для сравнения. Аргумент «альфа» в scikit-learn.

Python:

regr = Ridge(alpha=11498)
regr.fit(X, y)

http://nbviewer.ipython.org/github/JWarmenhoven/ISL-python/blob/master/Notebooks/Chapter%206.ipynb

Р:

Обратите внимание, что аргумент alpha=0в glmnet()означает, что следует применить штраф L2 (регрессия Риджа). Документация предупреждает, что нельзя вводить одно значение для lambda, но результат такой же, как в ISL, где используется вектор.

ridge.mod <- glmnet(x,y,alpha=0,lambda=11498)

Что вызывает различия?

Редактировать:
при использовании penalized()из штрафного пакета в R, коэффициенты такие же, как с scikit-learn.

ridge.mod2 <- penalized(y,x,lambda2=11498)

Возможно, тогда может возникнуть вопрос: «В чем разница между регрессией Риджа glmnet()и penalized()во время нее?

Новая оболочка Python для реального кода на Фортране, используемая в пакете R glmnet
https://github.com/civisanalytics/python-glmnet

Jordi
источник
5
Абсолютно незнаком с регрессией RGB. Но по умолчанию, sklearn.linear_model.Ridgeоценка непересекаемого перехвата (стандартная) и штраф такой, что ||Xb - y - intercept||^2 + alpha ||b||^2минимизируется для b. Там могут быть факторы 1/2или 1/n_samplesили оба перед штрафом, что сразу же меняет результаты. Чтобы устранить проблему масштабирования штрафа, установите штраф в обоих случаях на 0, устраните все расхождения и проверьте, что делает добавление штрафа. И между прочим, ИМХО, вот правильное место, чтобы задать этот вопрос.

Ответы:

9

В моем ответе отсутствует коэффициент , см. Ответ @visitors ниже для правильного сравнения.1N


Вот две ссылки, которые должны прояснить отношения.

Документация sklearn говорит, что linear_model.Ridgeоптимизирует следующую целевую функцию

|Xβy|22+α|β|22

Бумага glmnet говорит, что эластичная сеть оптимизирует следующую целевую функцию

|Xβy|22+λ(12(1α)|β|22+α|β|1)

Обратите внимание, что две реализации используют совершенно по-разному, sklearn использует для общего уровня регуляризации, в то время как glmnet использует для этой цели, резервируя для торговли между риджем и лассо-регуляризацией. ααλα

Сравнивая формулы, похоже, что установка и в glmnet должна восстановить решение .α=0λ=2αsklearnlinear_model.Ridge

Мэтью Друри
источник
И я полностью пропустил это в комментарии @eickenberg. Я должен использовать standardize = FALSEв , glmnet()чтобы получить одинаковые результаты.
Jordi
@Jordi Вы должны определенно стандартизировать, если используете linear_model.Ridgeдля анализа реального мира.
Мэтью Друри
Я понимаю, что linear_model.Ridgeмодель sklearn стандартизирует функции автоматически. Нормализация не обязательна. Я удивляюсь, почему мне нужно деактивировать стандартизацию, glmnet()чтобы модели давали идентичные результаты.
Джорди
10

Ответ Мэтью Друри должен иметь коэффициент 1 / N. Точнее...

В документации glmnet говорится, что эластичная сеть минимизирует функцию потерь

1NXβy22+λ(12(1α)β22+αβ1)

Документация Sklearn говорит, что linear_model.Ridgeминимизирует функцию потерь

Xβy22+αβ22

что эквивалентно минимизации

1NXβy22+αNβ22

Чтобы получить одно и то же решение от glmnet и sklearn, обе их функции потерь должны быть равны. Это означает установку и в glmnet.α=0λ=2Nαsklearn

library(glmnet)
X = matrix(c(1, 1, 2, 3, 4, 2, 6, 5, 2, 5, 5, 3), byrow = TRUE, ncol = 3)
y = c(1, 0, 0, 1)
reg = glmnet(X, y, alpha = 0, lambda = 2 / nrow(X))
coef(reg)

Выход glmnet: –0.03862100, –0.03997036, –0.07276511, 0.42727955

import numpy as np
from sklearn.linear_model import Ridge
X = np.array([[1, 1, 2], [3, 4, 2], [6, 5, 2], [5, 5, 3]])
y = np.array([1, 0, 0, 1])
reg = Ridge(alpha = 1, fit_intercept = True, normalize = True)
reg.fit(X, y)
np.hstack((reg.intercept_, reg.coef_))

выход склеарна: –0,03862178, –0,0399697, –0,07276535, 0,42727921

посетитель
источник
4
Различные определения параметров и их масштабирование, используемые в разных библиотеках, являются распространенным источником путаницы.
AaronDefazio
1
Я не ожидал бы, что и Гун, и я ошиблись бы.
Майкл Р. Черник
2
Да, вы оба ошиблись. Ваши причины отклонения моего редактирования дают понять, что вы оба не видели мой комментарий "Отсутствует коэффициент 1 / N" на stats.stackexchange.com/review/suggested-edits/139985
посетитель
Ваше редактирование, вероятно, было отклонено, потому что оно изменилось гораздо больше, чем то, что вы заявляете. Если вы хотите отредактировать мой пост и изменить только отсутствующий фактор, сделайте это, но изменение моих ссылок, формулировок и кода также является излишним. Комментарии о вашем несправедливом обращении в вашем ответе неуместны и не имеют отношения к содержанию вопроса, пожалуйста, удалите их. Ваша формулировка также омрачила мой ответ, это неправильный способ ответить на отклоненное редактирование. Мы хотели бы, чтобы вы внесли ценный вклад в наше сообщество, но, пожалуйста, ознакомьтесь с нашими нормами, прежде чем нас покинуть.
Мэтью Друри
1
@visitor Извините, если я немного грубоват. Я действительно должен просто попытаться сообщить, что вы, кажется, хороший потенциальный участник сайта, и я хочу, чтобы у вас был хороший опыт. У нас есть некоторые социальные нормы, как и у любой другой группы, и у вас будет лучший опыт, если вы будете их знать. Я все еще думаю, что «ответ Мэтью Друри неправильный» довольно резкий, есть, конечно, лучшие способы сообщить, что в моем ответе ошибочно отсутствует фактор . «X ответ неверен» читается как личная атака. 1N
Мэтью Друри