Уменьшите эффекты выбросов с помощью робастной регрессии

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

Зачем использовать робастную регрессию?

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

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

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

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

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

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

Сравнение результатов стандартной и робастной аппроксимации методом наименьших квадратов

Открыть 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 являются значениями рычага методом наименьших квадратов. Рычаги корректируют невязки путем уменьшения веса точек данных с высоким плечом, которые оказывают большой эффект на аппроксимацию методом наименьших квадратов (см. Hat Matrix и Leverage).
Стандартизируйте невязки. Стандартизированные скорректированные невязки даны
$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 являются невязками.
Итерация останавливается, если подгонка сходится или достигнуто максимальное количество итераций. В противном случае выполните следующую итерацию аппроксимации методом наименьших квадратов, вернувшись ко второму шагу.

См. также

fitlm | LinearModel | plotResiduals | robustfit

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