Уменьшайте эффекты выброса Используя устойчивую регрессию

Можно уменьшать эффекты выброса в моделях линейной регрессии при помощи устойчивой линейной регрессии. Эта тема задает устойчивую регрессию, показывает, как использовать его, чтобы подбирать линейную модель и сравнивает результаты со стандартной подгонкой. Можно использовать fitlm с 'RobustOpts' аргумент пары "имя-значение", чтобы подбирать устойчивую модель регрессии. Или можно использовать robustfit просто вычислить устойчивые параметры коэффициента регрессии.

Почему использование устойчивая регрессия?

Устойчивая линейная регрессия менее чувствительна к выбросам, чем стандартная линейная регрессия. Стандартная линейная регрессия использует обычные наименьшие квадраты, соответствующие, чтобы вычислить параметры модели, которые связывают данные об ответе с данными о предикторе с одним или несколькими коэффициентами. (Дополнительную информацию см. в Оценке Многомерных Моделей Регрессии.) В результате выбросы имеют большое влияние на подгонку, потому что обработка на квадрат остаточным значениям увеличивает эффекты этих экстремальных точек данных. Модели, которые используют стандартную линейную регрессию, описанную в том, Что такое Модель Линейной регрессии?, основаны на определенных предположениях, таких как нормальное распределение ошибок в наблюдаемых ответах. Если распределение ошибок асимметрично или подвержено выбросам, предположения модели делаются недействительным, и оценки параметра, доверительные интервалы, и другие вычисленные статистические данные становятся ненадежными.

Устойчивое использование регрессии метод вызвало итеративно повторно взвешенные наименьшие квадраты, чтобы присвоить вес каждой точке данных. Этот метод менее чувствителен к большим изменениям в мелких деталях данных. В результате устойчивая линейная регрессия менее чувствительна к выбросам, чем стандартная линейная регрессия.

Итеративно повторно взвешенные наименьшие квадраты

В методе взвешенных наименьших квадратов подходящий процесс включает вес как дополнительный масштабный коэффициент, который улучшает подгонку. Веса определяют, насколько каждое значение отклика влияет на итоговые оценки параметра. Низкокачественная точка данных (например, выброс) должна иметь меньше влияния на подгонку. Чтобы вычислить веса w _i, можно использовать предопределенные функции веса, такие как функция bisquare Туки (см. аргумент пары "имя-значение" 'RobustOpts' \in fitlm для большего количества опций).

Итеративно перевзвешенный алгоритм наименьших квадратов автоматически и итеративно вычисляет веса. При инициализации алгоритм присваивает равный вес каждой точке данных и оценки коэффициенты модели с помощью обычных наименьших квадратов. В каждой итерации алгоритм вычисляет веса w _i, давая более низкий вес точкам дальше от предсказаний модели в предыдущей итерации. Алгоритм затем вычисляет коэффициенты модели b с помощью метода взвешенных наименьших квадратов. Итерация останавливается, когда значения содействующих оценок сходятся в заданном допуске. Этот алгоритм одновременно стремится найти кривую, которая соответствует объему данных с помощью подхода наименьших квадратов, и минимизировать эффекты выбросов.

Для получения дополнительной информации смотрите Шаги для Итеративно Перевзвешенных Наименьших квадратов.

Сравните результаты стандартного и устойчивого метода наименьших квадратов

Скрипт Open Live Script

В этом примере показано, как использовать устойчивую регрессию с fitlm функция, и сравнивает результаты устойчивой подгонки к стандартному методу наименьших квадратов.

Загрузите moore данные. Данные о предикторе находятся в первых пяти столбцах, и данные об ответе находятся в шестом.

load moore
X = moore(:,1:5);
y = moore(:,6);

Соответствуйте наименьшим квадратам линейная модель к данным.

mdl = fitlm(X,y)

mdl = 
Linear regression model:
    y ~ 1 + x1 + x2 + x3 + x4 + x5

Estimated Coefficients:
                    Estimate          SE          tStat       pValue 
                   ___________    __________    _________    ________

    (Intercept)        -2.1561       0.91349      -2.3603      0.0333
    x1             -9.0116e-06    0.00051835    -0.017385     0.98637
    x2               0.0013159     0.0012635       1.0415     0.31531
    x3               0.0001278    7.6902e-05       1.6618     0.11876
    x4               0.0078989         0.014      0.56421     0.58154
    x5              0.00014165    7.3749e-05       1.9208    0.075365


Number of observations: 20, Error degrees of freedom: 14
Root Mean Squared Error: 0.262
R-squared: 0.811,  Adjusted R-Squared: 0.743
F-statistic vs. constant model: 12, p-value = 0.000118

Подбирайте устойчивую линейную модель к данным при помощи 'RobustOps' аргумент пары "имя-значение".

mdlr = fitlm(X,y,'RobustOpts','on')

mdlr = 
Linear regression model (robust fit):
    y ~ 1 + x1 + x2 + x3 + x4 + x5

Estimated Coefficients:
                    Estimate         SE         tStat       pValue 
                   __________    __________    ________    ________

    (Intercept)       -1.7516       0.86953     -2.0144    0.063595
    x1             1.7006e-05    0.00049341    0.034467     0.97299
    x2             0.00088843     0.0012027      0.7387     0.47229
    x3             0.00015729    7.3202e-05      2.1487    0.049639
    x4              0.0060468      0.013326     0.45375     0.65696
    x5             6.8807e-05    7.0201e-05     0.98015     0.34365


Number of observations: 20, Error degrees of freedom: 14
Root Mean Squared Error: 0.249
R-squared: 0.775,  Adjusted R-Squared: 0.694
F-statistic vs. constant model: 9.64, p-value = 0.000376

Визуально исследуйте остаточные значения этих двух моделей.

tiledlayout(1,2)
nexttile
plotResiduals(mdl,'probability')
title('Linear Fit')
nexttile
plotResiduals(mdlr,'probability')
title('Robust Fit')

Figure contains 2 axes. Axes 1 with title Linear Fit contains 2 objects of type line. Axes 2 with title Robust Fit contains 2 objects of type line.

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

Найдите индекс выброса.

outlier = find(isoutlier(mdlr.Residuals.Raw))

outlier = 1

Постройте веса наблюдений в устойчивой подгонке.

figure
b = bar(mdlr.Robust.Weights);
b.FaceColor = 'flat';
b.CData(outlier,:) = [.5 0 .5];
xticks(1:length(mdlr.Residuals.Raw))
xlabel('Observations')
ylabel('Weights')
title('Robust Fit Weights')

Figure contains an axes. The axes with title Robust Fit Weights contains an object of type bar.

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

Шаги для итеративно перевзвешенных наименьших квадратов

Итеративно перевзвешенный алгоритм наименьших квадратов выполняет эту процедуру:

Начните с первоначальной оценки весов и подбирайте модель методом взвешенных наименьших квадратов.
Вычислите настроенные остаточные значения. Настроенными остаточными значениями дают
$r_{adj} = \frac{r_{i}}{\sqrt{1 - h_{i}}}$
где r_{, i} является обычными остаточными значениями наименьших квадратов и h _i, является значениями рычагов метода наименьших квадратов. Рычаги настраивают остаточные значения путем сокращения веса точек данных высоких рычагов, которые оказывают большое влияние на метод наименьших квадратов (см. Матрицу Шляпы и Рычаги).
Стандартизируйте остаточные значения. Стандартизированными настроенными остаточными значениями дают
$u = \frac{r_{adj}}{K s} = \frac{r_{i}}{K s \sqrt{1 - h_{i}}}$
где K является постоянной настройкой, и s является оценкой стандартного отклонения остаточного члена, данного s = MAD/0.6745.
MAD среднее абсолютное отклонение остаточных значений их медианы. Постоянные 0.6745 делают оценку несмещенной для нормального распределения. Если матрица данных предиктора X имеет столбцы p, программное обеспечение исключает самый маленький p абсолютные отклонения при вычислении медианы.
Вычислите устойчивые веса w _i в зависимости от u. Например, bisquare весами дают
$w_{i} = {\begin{array}{l} {(1 - u_{i}^{2})}^{2} & , | u_{i} | < 1 \\ 0 & , | u_{i} | \geq 1 \end{array}$
Оцените устойчивые коэффициенты регрессии b. Веса изменяют выражение для оценок параметра b можно следующим образом
$b = \hat{β} = {(X^{T} W T)}^{- 1} X^{T} W y$
где W является диагональной матрицей веса, X является матрицей данных предиктора, и y является вектором отклика.
Оцените ошибку метода взвешенных наименьших квадратов
$e = \sum_{1}^{n} w_{i} {(y_{i} - {\hat{y}}_{i})}^{2} = \sum_{1}^{n} w_{i} r_{i}^{2}$
где w_{, i} является весами, y _i, является наблюдаемыми ответами, ŷ_{, i} является подходящими ответами и r_{, i} является остаточными значениями.
Итерация останавливается, если подгонка сходится, или максимальное количество итераций достигнуто. В противном случае выполните следующую итерацию наименьших квадратов, соответствующих путем возврата к второму шагу.

Документация