д.т.н.
Южно-Российский государственный политехнический университет
профессор
Шестопал О.В., аспирант кафедры прикладной математики Южно-Российского государственного политехнического университета имени М.В.Платова
УДК 519.6 : 621.316
Неизвестную нам связь между выходной величиной Y и факторами Xl обозначим в виде Y=f(Xl, …, Xn), где n – объем выборки Xl, …, Xn.
Робастная множественная регрессия
Обзор методов построения робастной регрессии можно найти в [1-6].
В значительной степени эти подходы используют функции Хубера, основанные на оценке медианы абсолютного отклонения [1-4].
Одним из зарекомендовавших себя робастных методов статистики является метод Тейла-Сена. Как правило, его применяют для построения одномерной линейной регрессии. В литературе сообщается о нескольких вариантах его обобщения на случай многомерной линейной регрессии. Наиболее зарекомендовал себя медианный метод на основе метода Гаусса-Зейделя.
Сформулируем его в следующем виде. Предположим, что требуется построить робастную многомерную линейную регрессию в Rn: y=a0 +a1x1+…+anxn. Для отыскания ее коэффициентов применим медианный метод на основе процедуры Гаусса-Зейделя с выделением диагональных элементов:
ak(s+1)= Me{[yi – yj – a1(s+1)(x1,i – x1,j) –…– ak–1(s+1)(xk–1,i – x k–1,j) –
ak+1(s)(xk+1,i – x k+1,j)–…– an(s)(xn,i – x n,j)]/ (xk,i – x k,j)}, k=1,…,n;
a0(s+1)= Me(yi – a1(s+1)x1,i –…– an(s+1)xn,i);
где Me – статистическая оценка медианы по парам точек выборки i,j=1,…,n(i не равно j); s– номер итерации по методу Зейделя.
Понятно, что должно быть задано начальное приближение, которое может быть выбрано разнообразными способами, например, часто нулевым.
Записанная система уравнений относительно коэффициентов a0,a1,…,an является, вообще говоря, нелинейной, так как медиана не всегда обладает линейным свойством. Часто можно решать данную систему методом простой итерации или более эффективным методом Зейделя для нелинейной системы уравнений. Если у системы нет диагонального преобладания, то метод Зейделя может расходиться. Для преодоления этой проблемы автором предлагается метод решения на основе метода релаксации:
ak(s+1)= (1– rk)Me{[yi – yj – a1(s+1)(x1,i – x1,j) –…– ak–1(s+1)(xk–1,i – x k–1,j) –
ak+1(s)(xk+1,i – x k+1,j)–…– an(s)(xn,i – x n,j)]/ (xk,i – x k,j)} + rkak, k=1,…,n;
a0(s+1)= Me(yi – a1(s+1)x1,i –…– an(s+1)xn,i);
rk– коэффициент релаксации.
Пример. Построить робастную линейную регрессию в R2: y=a0 +a1x1 +a2x2 по данным статистической выборки из 10-100 точек, равномерно распределенных в единичном квадрате. Точное решение – уравнение плоскости y= 1+ x1 + 10x2. Выборочные значения функции yiполучались суммированием точных значений с равномерно распределенным шумом с амплитудой 0.25. Примерно в одном из 10 случаев моделировался выброс, равный 1.
Сходимость достигалась менее чем за 10 итераций. Коэффициенты регрессии вычислялись с точностью от 20% до 7% .
Обобщение метода Тейла-Сена на случай нелинейной регрессии
Наиболее употребительной является полиномиальная регрессия [1].
Сформулируем метод в следующем виде. Предположим, что требуется построить робастную квадратическую регрессию: y=a0 +a1x+a2x2. Если применить для отыскания ее коэффициентов медианный метод на основе процедуры Гаусса-Зейделя с выделением диагональных элементов по аналогии с линейной регрессией, то очень часто имеет место расходимость итерационного процесса. Система уравнений и соотношения метода Зейделя при этом имеют вид:
a1= Me[yi – yj – a2(xi 2– xj2)]/ (xi – x j)],
a2= Me[yi – yj – a1(xi– xj)]/ (xi 2– xj2)],
a0= Me(yi – a1xi – a2xi 2); (1)
a1(s+1)= Me[yi – yj – a2(s)(xi 2– xj2)]/ (xi – x j)],
a2(s+1)= Me[yi – yj – a1(s+1)(xi– xj)]/ (xi 2– xj2)],
a0(s+1)= Me(yi – a1(s+1)xi – a2(s+1)xi 2);
где Me – статистическая оценка медианы по парам точек выборки: i,j=1,…,n (i не равно j); s – номер итерации по методу Зейделя.
Для сходимости часто требуется задание очень точного начального приближения.
Требуемое количество итераций часто является непомерным. Автором предложена эффективная вычислительная процедура решения системы медианных уравнений. Перепишем первое уравнение системы в виде:
a1= Me[(yi– yj)/(xi–xj)– a2(xi+xj)].
Учтем, что корреляция двух слагаемых под знаком медианы относительно мала. Уравнение запишем в виде:
a1= Me[(yi– yj)/(xi–xj)] – a2Me[xi+xj] + r,
где невязка
r = Me[(yi – yj)/(xi – x j)– a2(xi + x j)] – Me[(yi – yj)/(xi – x j)+ a2(xi + x j)] (2)
часто, если не как правило, принимает относительно малые значения. Выразим второй коэффициент через первый:
a2= {Me[(yi – yj)/(xi – x j)] – a1 + r }/ Me[xi + x j]; (3)
и далее подставим это выражение в общее уравнение регрессии y=a0 +a1x+a2x2:
y=a0 + a1x + x2{Me[(yi – yj)/(xi – x j)] – a1 + r }/ Me[xi + x j],
откуда
y1(x) = y(x) – x2{Me[(yi – yj)/(xi – x j)] + r}/ Me[xi + x j],
y1(x) = a0 + a1(x– x2/Me[xi + x j]) = a0 + a1x1,
откуда находим медианную оценку:
a1= Me[(y1i– y1j)/(x1i–x1j)]. (4)
Далее система (1)-(4) решается методом Зейделя в следующем порядке: построение вариационного ряда, вычисление Me[xi+xj], выборочных значений yi. Следующий шаг – вычисление Me[(yi– yj)/(xi–xj)] и далее в цикле итерационного метода Зейделя поправка r, значения x1i,y1i, соответственно (4) значение a1, согласно (1) значениеa0, и наконец значение a2, в соответствии с (3).
Пример. Построить робастную нелинейную регрессию: y= a0 +a1x+a2x2 по данным статистической выборки из 10-100 точек, равномерно распределенных на интервале (0,1). Точное решение – уравнение параболы y= 1+ x+x2. Выборочные значения функции yiполучались суммированием точных значений с равномерно распределенным шумом с амплитудой 0.5. Примерно в одном из 10 случаев моделировался выброс, равный 5.
Сходимость достигалась за 1-2 итерации. Коэффициенты регрессии вычислялись с точностью от 20% до 7% (Рисунок 1).
Рисунок 1. Вычисление коэффициентов нелинейной регрессии
На рисунке МТС2 обозначает нелинейный вариант метода Тейла-Сена, MNK – метод наименьших квадратов, 100 точек, выбросы через 10 точек, равны 5. Пять итераций. Точность МТС2 много выше.
На рисунке 2 точки и черная линия – это выборочные значения и точное решение, синяя линия – робастный нелинейный вариант метода Тейла-Сена (высокая точность), красная линия – метод наименьших квадратов (точность плохая).
Рисунок 2. Кривые нелинейной регрессии
Рецензии:
5.11.2018, 14:41 Мирмович Эдуард Григорьевич
Рецензия: Работа хорошая, есть элементы (хоть и небольшие) новизны, использование медианного распределения при оценке отклонений рецензенту нравятся. Грамматику и синтаксис подправить, пробелы и др. аспекты формативрования подправить.Рецензенту хотелось бы, чтобы автор посмотрел работы, которые есть в Интернете:
1.Мирмович Э.Г. Некрасов О.Г. Интерполирование и аппроксимация данных полиномами степенного, экспоненциального и тригонометрического вида // Научные и образовательные проблемы гражданской защиты. Научный журнал. Химки: АГЗ МЧС России. № 4. 2010. С. 23–32. http://www.amchs.ru/images/stories/ucheba/ng/4_2010.pdf
2.Некрасов О.Н. Аппроксимация экспериментальных данных на основе линеаризуемых функций // Научные и образовательные проблемы гражданской защиты. Научный журнал. Химки АГЗ МЧС России. №1. 2008. – С. 81–88.
Может окажется возможным ссылку дать в обзорной части. Рекомендetncz к публикации после учёта замечаний рецензента.
С уважением к автору.
6.11.2018, 13:54 Соловьёв Виктор Григорьевич Отзыв: Уважаемый Сергей Александрович! С большим интересом прочитал Вашу публикацию решения задачи аппроксимации данных итерационным методом при использовании медианного распределения при оценке отклонений. Единственное, на что обращаю Ваше внимание: в работе показано сравнение метода МТС2 с методом MNK с аналитическим и графическим доказательством преимущества первого над вторым. Однако, даже по приведенному Вами графику видно явное несоответствие красной кривой линии методу MNK, так как сумма квадратов отклонений никак не минимальна, вследствие чего красная линия приподнята по оси ординат как минимум на 0.2. Тогда коэффициент k_0 должен быть существенно меньше 1.74. Рецензент не обратил на это внимание. По-всей видимости, расчет по методу MNK выполнен с ошибками. Другой вопрос, что Ваш метод дает определенные преимущества над MNK, особенно на концах отрезка. Поэтому есть смысл скорректировать данные по методу MNK и произвести соответствующие сравнения. С уважением, Виктор Григорьевич Соловьёв |
6.11.2018, 21:15 Некрасов Сергей Александрович Отзыв: Уважаемый Виктор Григорьевич! Вы не обратили внимание на условия проведения вычислительного эксперимента с МНК. Моделировались большие выбросы (10%). На графике они из-за масштаба не видны. Они-то (эти большие положительные выбросы) и потянули кривую МНК вверх. А кривая МТС2 по причине отличной устойчивости к выбросам медианной оценки почти не сдвинулась от нормального положения. С уважением, Некрасов С.А. |