Модель полиномиальной регрессии16.06.2018 10:02

Выражаясь простым языком, модель регрессии в математической статистике строится на основе известных данных, в роли которых выступают пары чисел. Количество таких пар заранее определено. Если представить себе, что первое число в паре — это значение координаты $x$ , а второе — $y$ , то множество таких пар чисел можно представить на плоскости в декартовой системе координат в виде множества точек. Данные пары чисел берутся не случайно. На практике, как правило, второе число зависит от первого. Построить регрессию — это значит подобрать такую линию (точнее, функцию), которая как можно точнее приближает к себе (аппроксимирует) множество вышесказанных точек.

Для чего это всё нужно? Прежде всего, это нужно для составления т.н. прогнозов. Часто нужно узнать $y$ , зная только $x$ , если он отличается от тех иксов, на основе которых строилась регрессия. Приведу простой пример. Имеется статистика зависимости роста человека от его возраста на основе 100 исследованных разных человек. Таким образом, мы имеем 100 пар чисел {возраст; рост}. При этом «рост» — зависимая величина, а «возраст» — независимая. Грамотно построив модель регрессии, мы можем с определённой достоверностью «предсказать» рост по любому значению возраста.

На практике, в зависимости от ситуации, в построении моделей регрессии используют линейные, параболические, степенные и другие типы функций. В курсе математической статистики чаще всего рассматривается линейная модель регрессии. Иногда затрагивают случай посложнее — параболическую модель. Делая обобщение, несложно догадаться, что линейная и параболическая модели являются частными случаями более сложной модели — полиномиальной. Построить модель регрессии — это значит найти параметры той функции, которая будет в ней фигурировать. Для линейной регрессии — два параметра: коэффициент и свободный член.

Полиномиальная регрессия может применяться в математической статистике при моделировании трендовых составляющих временных рядов. Временной ряд — это, по сути, ряд чисел, которые зависят от времени. Например, средние значения температуры воздуха по дням за прошедший год, или доход предприятия по месяцам. Порядок моделируемого полинома оценивается специальными методами, например, критерием серий. Цель построения модели полиномиальной регрессии в области временных рядов всё та же — прогнозирование.

Для начала рассмотрим задачу полиномиальной регрессии в общем виде. Все рассуждения основаны на обобщении рассуждений в задачах линейной и параболической регрессии. После этих рассуждений я перейду к частному случаю — рассмотрению данной модели для временных рядов.

Пусть даны два ряда наблюдений $x_i$ (независимая переменная) и $y_i$ (зависимая переменная), $i=\overline{1,n}$ . Уравнение полинома имеет вид

$y=\sum\limits_{j=0}^k b_jx^j,\ \ \ \ \ (1)$

где $b_j$ — параметры данного полинома, $j=\overline{0,k}$ . Среди них $b_0$ — свободный член. Найдём по методу наименьших квадратов (МНК) параметры $b_j$ данной регрессии.

По аналогии с линейной регрессией, МНК также основан на минимизации следующего выражения:

$S=\sum\limits_{i=1}^n\left(\hat y_i-y_i\right)^2\to\min\ \ \ \ \ (2)$

Здесь $\hat y_i$ — теоретические значения, являющиеся значениями полинома (1) в точках $x_i$ . Подставляя (1) в (2), получаем

$S=\sum\limits_{i=1}^n\left(\sum_{j=0}^kb_jx_i^j-y_i\right)^2\to\min.$

На основании необходимого условия экстремума функции $(k+1)$ переменных $S=S(b_0, b_1,\dots,b_k)$ приравняем к нулю её частные производные, т.е.

$S'_{b_p}=2\sum\limits_{i=1}^nx_i^p\left(\sum\limits_{j=0}^kb_jx_i^j-y_i\right)=0,\ \ \ p=\overline{0,k}.$

Поделив левую и правую часть каждого равенства на 2, раскроем вторую сумму:

$\sum\limits_{i=1}^nx_i^p\left(b_0+b_1x_i+b_2x_i^2+\dots+b_kx_i^k\right)-\sum\limits_{i=1}^nx_i^py_i=0,\ \ \ p=\overline{0,k}.$

Раскрывая скобки, перенесём в каждом $p$ -ом выражении последнее слагаемое с $y_i$ вправо и поделим обе части на $n$ . В результате у нас получилось $(k+1)$ выражений, образующие систему линейных нормальных уравнений относительно $b_p$ . Она имеет следующий вид:

$\left\{ \begin{array}{l} b_0+b_1\overline x+b_2\overline{x^2}+\dots+b_k\overline{x^k}=\overline y\\ b_0\overline x+b_1\overline{x^2}+b_2\overline{x^3}+\dots+b_k\overline{x^{k+1}}=\overline{xy}\\ b_0\overline{x^2}+b_1\overline{x^3}+b_2\overline{x^4}+\dots+b_k\overline{x^{k+2}}=\overline{x^2y}\\ \ldots\ldots\ldots\ldots\ldots\ldots\ldots\ldots\ldots\ldots\ldots\ldots\ldots\\ b_0\overline{x^k}+b_1\overline{x^{k+1}}+b_2\overline{x^{k+2}}+\dots+b_k\overline{x^{2k}}=\overline{x^ky} \end{array} \right.\ \ \ \ \ (3)$

Можно переписать систему (3) в матричном виде: $AB=C$ , где

$A=\left( \begin{array}{ccccc} 1 & \overline x & \overline{x^2} & \ldots & \overline{x^k}\\ \overline x & \overline{x^2} & \overline{x^3} & \ldots & \overline{x^{k+1}}\\ \overline{x^2} & \overline{x^3} & \overline{x^4} & \ldots & \overline{x^{k+2}}\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ \overline{x^k} & \overline{x^{k+1}} & \overline{x^{k+2}} & \ldots & \overline{x^{2k}} \end{array} \right),\ \ B=\left(\begin{array}{c} b_0\\b_1\\b_2\\\vdots\\b_k \end{array} \right),\ \ C=\left(\begin{array}{c} \overline y\\\overline{xy}\\\overline{x^2y}\\\vdots\\\overline{x^ky} \end{array} \right).$

Теперь перейдём к применению вышеизложенных фактов в случае с временными рядами. Пусть дан временной ряд $x_t$ , где $t=\overline{1,n}$ . Требуется построить полиномиальный тренд порядка $k$ , который аппроксимирует данный временной ряд как можно точнее. В качестве независимой переменной $x$ будем брать $t$ , исходя из определения временного ряда. Эти иксы представляют собой ряд натуральных чисел, обозначающие период времени. В качестве $y$ берутся значения временного ряда $x_t$ . При этом видно, что значения элементов $a_{ij}$ матрицы системы $A$ не зависят от $x_t$ . Так как в общем случае, очевидно,

$a_{ij}=\overline{x^{i+j-2}}=\frac1n\sum\limits_{r=1}^nx_r^{i+j-2},$

то в случае с временными рядами

$a_{ij}=\frac1n\sum\limits_{r=1}^nr^{i+j-2},$

где $i,j=\overline{1,(k+1)}.$

Элементы $c_j$ матрицы-вектора свободных членов $C$ в общем случае получаются, как

$c_j=\overline{x^{j-1}y}=\frac1n\sum\limits_{r=1}^nx_r^{j-1}y_r.$

А в случае с временными рядами

$c_j=\frac1n\sum\limits_{r=1}^nr^{j-1}x_r,$

где $j=\overline{1,(k+1)}.$

Таким образом, решив систему (3), мы сможем найти искомые параметры полиномиального тренда $b_0,\dots,b_k.$

Для заполнения матриц системы и для её решения можно применять один из численных методов при моделировании тренда на компьютере. При этом результат вычисления будет достаточно точным.

В результате трендовая составляющая примет вид:

$T_t=\sum\limits_{i=0}^kb_it^i,\ \ \ t=0,1,2,\dots.$

Стоит также отметить, что смоделированная трендовая составляющая $T_t$ , определена не только на текущие периоды $[1;n]$ , но и на будущие периоды © Habrahabr.ru