Кинетические уравнения для прогнозирования нестационарных временных рядов

( Kinetic equations for forecasting of non-stationary time-series
Preprint, Inst. Appl. Math., the Russian Academy of Science)

Введение

В настоящей работе исследуется возможность применения кинетических уравнений к задаче прогнозирования нестационарных временных рядов. Такие ряды встречаются во многих практических задачах, связанных с обработкой данных и прогнозированием в медицине, геологии, биологии, экономике и других областях. Случайные временные ряды (но, как правило, стационарные) возникают также и при численном решении некоторых систем обыкновенных дифференциальных уравнений, когда при определенных значениях параметров дискретизации получается разбегание траекторий для сколь угодно близких начальных условий. Использование кинетических уравнений для описания прогнозируемого изменения выборочной функции распределения временного ряда является новым подходом, который для нестационарных случайных процессов имеет определенные преимущества по сравнению с традиционными адаптивными методами анализа данных.

Основная проблема при исследовании нестационарных рядов состоит в том, что любая конечная выборка не принадлежит генеральной совокупности, т.е. не существует стационарного распределения вероятности, реализацией которого является данный случайный процесс. В силу этого обстоятельства нет и теорем об асимптотической нормальности отклонений выборочных моментов от моментов генеральной совокупности. Но тогда с увеличением объема выборки не происходит уточнения статистических свойств временного ряда. Отсутствует также и такое традиционное (для стационарных процессов) свойство как асимптотическая несмещенность оценок выборочных моментов.

Напомним некоторые определения из теории анализа выборочных статистик, важные для понимания отличий подхода, предлагаемого в настоящей работе, от традиционных методов, развитых для стационарных временных рядов [1, 2].

Оценкой величины , получаемой в результате n наблюдений над некоторой случайной величиной , называется статистика, т.е. некоторая функция от результатов наблюдений, по значениям которой судят о значениях величины . Такой величиной, например, являются моменты гипотетического стационарного распределения, которым, как предполагается, определяются вероятности осуществления изучаемых событий.

Оценка  величины  называется несмещенной, если ее математическое ожидание равно оцениваемой величине: . Это равенство означает, что если бы нам было известно распределение вероятностей, то, усреднив по этому распределению наблюдаемую эмпирическую статистику, мы получили бы среднее значение параметра . Оценка называется асимптотически несмещенной, если ее математическое ожидание стремится к оцениваемой величине при увеличении объема выборки.

Оценка называется состоятельной, если она сходится по вероятности к оцениваемой величине: . Если оценка состоятельна, то увеличение объема выборки приводит к уточнению эмпирического значения изучаемой величины. Оценка называется суперсостоятельной, если указанная сходимость по вероятности имеет асимптотику, главный член которой убывает быстрее, чем , т.е. такая сходимость более быстрая, чем сходимость к нормальному распределению в рамках центральной предельной теоремы.

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

Доверительным интервалом с доверительной вероятностью p для, например, математического ожидания , оценкой которого по выборке объема n является выборочное среднее , называется интервал с границами

                                            ,

где

                               ,

а величина  является -квантилем стандартного нормального распределения, т.е. такой величиной , что

                                     ,

причем

есть вероятность того, что случайная величина  не превосходит по модулю значение , а  есть плотность нормального распределения с нулевым средним и единичной дисперсией. Такой подход к оценке средних величин основан на том, что для стационарных процессов отклонения выборочных средних от среднего по генеральной совокупности распределены асимптотически нормально.

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

В математическом плане задача прогноза вероятности того, что значение временного ряда будет лежать в заданных пределах, сводится к построению модели временной эволюции выборочной функции распределения (ВФР), построенной по выборке конечного объема, один или оба конца которой, возможно, скользящие.

В работах авторов [3, 4] предложена методика анализа нестационарных временных рядов, основанная на том факте, что для любого ряда можно подобрать выборку такого объема, что изменение построенной ВФР при сдвиге на любое конечное число шагов будет меньше любого наперед заданного числа. Это означает, что, задавая желательную точность анализа данных, можно определить промежуток времени (горизонт прогноза), на котором с указанной точностью справедливы методы, разработанные для стационарного случая. При этом важной частью методики является установление связи между объемом выборки, точностью прогноза и горизонтом прогноза. В настоящей работе показано, что этими величинами определяется также и тип моделей, которые предлагаются для описания эволюции ВФР.

1. Определение промежутка квазистационарности ВФР

Пусть  представляет значение случайной величины x в дискретный момент времени t. Обозначим  ВФР, полученную по выборке  объема T, где натуральное число n принадлежит скользящему окну . Символ геометрического отрезка  здесь означает набор соответствующих натуральных чисел. Среднее значение функции a(x) от случайной величины x по ВФР  будем обозначать угловыми скобками:

                                       .                                      (1)

В тексте для краткости изложения мы иногда будем опускать указание зависимости от объема выборки T и/или момента времени t, если это не будет препятствовать пониманию смысла предложения.

В каждый момент времени ВФР  представляет собой функцию одного аргумента x. Введем норму в пространстве функций распределения, являющихся суммируемыми по x функциями, т.е. в пространстве :

                                                    .                                           (2)

Определение 1. Расстоянием между двумя ВФР f и h называется величина

                     .                         (3)

В этом определении сравниваемые ВФР могут быть построены по выборкам разных объемов и/или быть отнесены к различным моментам времени.

Определение 2. ВФР  временного ряда  будем называть q-e-стационарной на временном промежутке , если

                    .                           (4)

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

Если неравенство (4) выполняется для всех q, то ВФР (и сам ряд) будем называть просто e-стационарным. Если неравенство (4) выполнено лишь для конкретного значения  (или для некоторого набора этих значений), то ряд и соответствующую ВФР будем называть ограниченно t-e-стационарными.

Если исследуемый ряд таков, что промежуток изменения величины x равномерно ограничен по времени, область интегрирования в (4) без ограничения общности можно считать отрезком . В этом случае из -стационарности распределения следует также и -стационарность его моментов. Функционал

                                 (5)

совместно с определением 2 играет роль критерия -стационарности распределения случайной величины x. Он показывает, насколько ВФР с одинаковым объемом статистической базы отличаются одна от другой в моменты времени, разделенные интервалом τ.

В [3] показано, что из неотрицательности ВФР и ее нормированности на единицу в любой момент времени следует оценка

                                        .                                             (6)

Неравенство (6) дает возможность сделать важный вывод о том, что при фиксированном τ функционал (5) равномерно ограничен по t. Поэтому  всегда можно подобрать такой объем выборки , что ВФР будет -стационарной. Таким образом, если нельзя сравнить выборочную функцию распределения со стационарным распределением, то можно добиться близости двух нестационарных выборочных распределений. Это позволит прогнозировать временной ряд с некоторой заданной точностью возможного изменения функции распределения его значений.

При увеличении точности в определении -стационарности, т.е. при уменьшении , объем выборки, при которой достигается условие , растет. В силу равномерной ограниченности для каждого момента времени t существует такое минимальное значение , что при всех  и при всех  значения функционала  не превосходят . Рассмотрим максимальное по периоду наблюдений значение статистики

                                               .                                           (7)

В силу оценки (6) имеем .

Статистику  будем называть горизонтным рядом для прогнозирования ряда . Ее анализ позволяет получить эмпирическую оценку максимального горизонта прогноза, внутри которого распределение остается -стационарным, а также оценку минимального объема выборки, необходимого для прогнозирования ВФР внутри заданного горизонта.

2. Точность прогнозирования ВФР и временного ряда

Вычисление текущих оценок статистических характеристик нестационарного процесса  приводит к появлению принципиально неустранимых погрешностей: погрешности, возникающей за счет конечности времени усреднения, т.е. вследствие недостаточной репрезентативности объема выборки, а также погрешности, возникающей за счет изменения статистики на интервале усреднения. Будем называть погрешность (в смысле среднего квадратичного) оценки статистики за счет конечности промежутка  погрешностью первого типа и обозначать , а погрешность за счет нестационарности статистики – погрешностью второго типа и обозначать . Задача оптимизации прогнозного алгоритма состоит в минимизации функционала полной ошибки  в оценке статистики. Предполагая, что процессы, приводящие к ошибкам первого и второго типов, независимы, в качестве функционала полной ошибки выбираем

                                                    .                                                       (8)

При прогнозировании возникает задача оптимизации объема T выборки в терминах функционала ошибки (8). Эта задача состоит в определении такого объема , при котором , и построении соответствующего численного алгоритма нахождения  для заданного нестационарного временного ряда. Оптимизационная задача возникает в силу того, что ошибки  и  как функции объема выборки имеют различное поведение. Именно, чтобы уменьшить ошибку вследствие неполной статистической репрезентативности, следует увеличивать объем , а для уменьшения влияния нестационарности на статистические характеристики временного ряда следует уменьшать объем выборки. Алгоритмические методы нахождения  описаны одним из авторов в [5].

Таким образом, формулировка условия близости двух ВФР объединяет две статистики данного временного ряда : статистику объема выборки  и статистику горизонта прогноза . Именно, требуется определить  так, чтобы ВФР, построенные по выборкам этих оптимальных объемов на промежутке , различались бы меньше, чем на заданное число, характеризующее ошибку прогноза.

Определение 3. Ошибкой  прогноза для временного ряда  будем называть среднеквадратичное отклонение прогнозных значений  от фактических  на промежутке горизонта прогноза:

                                          .                                   (9)

Определение 4. Ошибкой  прогноза ВФР в момент времени t будем называть интегральное абсолютное отличие прогнозной ВФР  от фактической , построенных по выборкам равных объемов:

                                            .                                       (10)

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

                       .                     (11)

Различие между прогнозным и фактическим средними выборочными величинами можно оценить из неравенства

                                           .                                           (12)

Рассмотрим ошибку прогнозирования значения ряда в некоторый момент времени. Прогнозным значением ряда в момент времени t естественно считать среднее значение . Формальной ошибкой такого прогноза является . Однако при сравнении с фактом следует учесть, что фактическая ВФР изменилась по сравнению с прогнозной, так что оценка квадрата фактической ошибки прогноза составит

                                        (13)

Минимизация верхней оценки ошибки прогноза, представляемой формулой (13), и есть задача минимизации ошибки (8). Идея такой минимизации состоит в следующем.

Введем плотность распределения  значений горизонтного ряда , т.е. определим вероятность того, что расстояние между двумя ВФР, построенных по выборке объема T и сдвинутых по времени на фиксированный промежуток , не превосходит  для всех . Для краткости зависимость от момента времени t в аргументах  опущена.

Эта плотность  строится по имеющимся данным к текущему моменту времени t следующим образом. Для каждого момента времени , , строятся ВФР по выборкам объемов , после чего для каждой из этих ВФР вычисляется функционал . Промежуток значений , принимаемых функционалом , разбивается при этом на некоторое количество N отрезков, например, равномерно, так что правый конец k-го отрезка есть . Номер отрезка, фиксирующий заданную в (5) величину близости  между ВФР, есть . Промежуток  возможных значений объемов выборок содержит  целочисленных точек, так что плоскость  покрыта -сетью. Отдельную ячейку сети нумеруем двумя индексами i,k. Если значение  соответствует номерам i,k некоторой ячейки, помещаем в нее индикатор «1». Ячейки, номерам которых не отвечают значения функционала , заполняем нулями. Совокупность ячеек с фиксированными значениями i или k образует полосу. В каждый момент времени  k-ая полоса состоит, таким образом, из нулей и единиц. Значением  горизонтного ряда в k-ой полосе называется индекс i ячейки, следующей за ячейкой с наибольшим индексом, содержащей ноль. Если нулевая ячейка – последняя, то значение горизонтного ряда в этой полосе отсутствует. Плотность  представляет собой выборочную функцию распределения временного ряда , построенную по выборке объема , доставляющей исходному ряду условие ограниченной t-e-стационарности (4) в определении 2.

Если требуется, чтобы было выполнено условие полной t-e-стационарности, то процедуру формирования членов ряда  следует дополнить рассмотрением всех значений  от 1 до . Именно, при каждом значении  рассматриваются величины , которые определяются по вышеописанному алгоритму. Затем находятся величины .

Определим интегральную функцию распределения горизонтного ряда

                            .                                 (14)

Смысл функции (14) в том, что если вместо величины  из (7), во многих случаях оказывающейся равной своему теоретическому максимуму , взять некоторое меньшее значение , то с вероятностью  будут выполнены условия соответственно ограниченной или полной t-e-стационарности.

Введенные понятия позволяют окончательно сформулировать задачу об оптимизации объема выборки для прогнозирования временного ряда на заданный горизонт  с заданной точностью  в смысле определения (9). Поскольку для t-e-стационарных рядов справедлива оценка (13) , то оптимальный объем для прогнозирования на горизонт  при соблюдении условия t-e-стационарности ВФР определим следующим образом.

Определение 5. Оптимальным объемом  для прогнозирования временного ряда  на горизонт  называется такой объем выборки, при котором сумма относительной дисперсии ВФР, построенной по выборке этого объема на текущий момент начала прогноза, и квадрата расстояния между двумя ВФР, сдвинутыми на интервал , минимальна.

Чтобы провести оптимизацию объема выборки согласно определению 5, положим сначала . Рассмотрим зависимость эмпирической выборочной дисперсии от объема выборки, усредненной по промежутку времени, на котором распределение горизонтного ряда t-e-стационарно. Уменьшение оптимального объема выборки играет принципиальную роль при прогнозировании, поскольку выборочная дисперсия часто ведет себя немонотонно как функция объема выборки. Простая оценка оптимального объема дается тогда формулой

                                      .                                   (15)

Более тонкая оптимизация основана на анализе суммы выборочной дисперсии и квадрата ожидаемого расстояния между двумя ВФР, сдвинутыми на промежуток . Эта оптимизация использует локальные свойства плотности распределения горизонтного ряда. В частности, может оказаться, что максимальных значений объемов выборки, равных , среди значений  весьма мало, так что соответствующая вероятность  близка к нулю. Выберем тогда некоторое значение . Для него функционал  не превосходит  с вероятностью , определяемой формулой (14). Однако для выборок объемов от  до  функционал нормы может быть больше . Этим увеличением неточности в близости двух ВФР можно пренебречь, если  есть -квантиль распределения , т.е. такая величина, для которой

                                                 .                                                   (16)

В качестве объема выборки  берется ближайшее целое число к решению уравнения (16). Тогда ожидаемая неточность функционала (5) будет не больше, чем

          .                  (17)

Последнее равенство в (17) следует из того, что при  .

Итак, в результате получается уточненная оценка оптимального объема выборки при сохранении требования t-e-стационарности ВФР:

                                          .                                   (18)

Наконец, при минимизации ошибки прогноза исходного временного ряда можно еще расширить допустимую неточность в прогнозе ВФР с тем, чтобы уменьшить совокупную ошибку прогноза за счет возможного уменьшения выборочной дисперсии. Описанная методика определения оптимального объема выборки является пошаговой с шагом, равным горизонту прогноза t.

Определение 6. Временной ряд , для которого оптимальный объем выборки при прогнозировании на заданный промежуток времени , является квазистационарным (т.е. t-e-стационарным) временным рядом, будем называть устойчиво прогнозируемым рядом.

Из (6) следует, что в качестве максимальной оценки минимально допустимого объема выборки  можно взять , т.е. при желании каждый временной ряд можно считать устойчиво прогнозируемым в смысле определения 6. Оптимизация же прогноза ряда заключается в том, что в некоторых случаях суммарная ошибка (13) может быть меньше для меньших допустимых объемов выборки. Поэтому желательно, чтобы именно эти оптимальные выборки представляли собой квазистационарную статистику.

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

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

3. Выборочные моменты и их эволюция

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

Выборочный момент первого порядка по выборке на промежутке  определяем как

                                              .                                                (19)

Выборочные центральные моменты порядка  определяем по формуле

                                   .                                        (20)

Рассмотрим, например, временные ряды, образованные дневной ценой закрытия акций компаний General Motors (GM), General Electric (GE) и Microsoft (MS) за 2004-2005 гг. по данным [6]. Результаты усреднения (19-20) по объему выборки в 900 значений представлены на Рис. 1-2.

Рис. 1. Среднее скользящее значение (левая шкала) и дисперсия (правая шкала) цен на акции GE.

Здесь показан только фрагмент графика эволюции моментов ряда GE. Локальное уменьшение среднего значения и дисперсии (абсолютной, а не относительной) обусловлено снижением цен, а не стабилизацией случайного процесса. График приведен именно с целью демонстрации существенного непостоянства во времени моментов распределения. Похожим образом выглядит зависимость от времени средних выборочных значений ВФР для цен на акции компании MS (Рис. 2):

Рис. 2. Среднее скользящее значение (левая шкала) и дисперсия (правая шкала) цен на акции MS.

На Рис. 3 приведены для сравнения скользящие средние временного ряда цен на электроэнергию на НОРЭМ, исследованного в нашей предыдущей работе [3].

Рис. 3. Среднее скользящее значение (левая шкала, Ряд 1) и дисперсия (правая шкала, Ряд 2) цен на электроэнергию на НОРЭМ.

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

4. Кинетические уравнения для ВФР

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

Введем двумерную ВФР , случайных величин x и  исходного ряда и ряда его производной, полученного взятием разностей  в соседние моменты времени. Можно рассмотреть ВФР, зависящую и от большего числа производных, т.е. функцию вида , но надо помнить об ограниченности объема выборки заданной величиной T. В силу конечности выборки невозможно эмпирически определить производную более высокого порядка, чем T-1. Более того, чтобы оперировать одинаковыми массивами данных по каждой из r производных, необходимо иметь исходный ряд в количестве T+r элементов.

Следовательно, чтобы можно было содержательно обсуждать задачу об эволюции ВФР , число доступных для анализа элементов ряда должно быть на 1 больше, чем в выборке. Такое включение естественным образом содержится в отрезке : крайний левый элемент  окна  будем считать виртуальным, т.е. не принадлежащим анализируемой выборке, в соответствии с определениями (19-20), а сами выборочные средние будем определять по данным , , …, . Тем самым мы будем считать, что в окне  содержится не T+1, а ровно T элементов.

Формально из сохранения во времени нормировки функции распределения следует уравнение Лиувилля

                             .                              (21)

Обычно предполагается, что дивергенция в (21) берется в пространстве фазовых переменных, в котором введена динамически-инвариантная мера. Это означает, что в системе нет источников вероятности, т.е. вероятностное пространство является полным. Тогда

                                                          ,

после чего уравнение (21) может быть записано в виде

                                               .

Такая запись согласуется с представлениями статистической механики, в которой компоненты фазовой координаты  считаются независимыми. Однако следует заметить, что в случае выборочного распределения конечного объема такой независимости может и не быть, т.к. корреляция между выборочными значениями  и  для рядов, встречающихся, например, в эконометрике, как правило, отлична от нуля. В то же время «интеграл» от совместного эмпирического распределения двух величин по одной из них, т.е. сумма по всем соответствующим ее значениям, даст эмпирическое распределение по оставшейся переменной в силу самого правила построения совместной ВФР. Таким образом, на стадии вывода уравнения эволюции эмпирической ВФР методически правильно считать компоненты  независимыми. То, что при одной реализации случайного процесса (т.е. на некоторой фазовой траектории) фазовые координаты коррелируют, не препятствует выбору удобного базиса в фазовом пространстве.

Имея целью создание прогнозной модели для «однофазной» ВФР , мы ограничимся рассмотрением кинетических уравнений для распределения одной переменной, хотя на стадии вывода соответствующего уравнения понадобится совместное распределение .

Поскольку

                                              ,                                               (22)

то из (21) следует уравнение Лиувилля

                                                    ,                                                   (23)

где введена функция  по формуле

                                           .                                      (24)

Введем средние по x величины в момент времени t:

          .                     (25)

Поскольку уравнение (23) записано относительно эмпирической ВФР, то будем его называть «эмпирическим» уравнением Лиувилля. Это означает, что скорость изменения ВФР  является параметром, определяемым по данной конкретной выборке в момент t, а не выводится из каких-либо уравнений, как в традиционной статистической механике.

Найдем скорость изменения скользящего среднего  в силу уравнения (23). Поскольку, как отмечалось выше, фактическая область изменения наблюдаемых величин x конечна, то мы не имеем возможности строить эмпирическую ВФР с бесконечной областью определения, т.к. все интегралы, вычисляемые численно, – принципиально собственные. Тогда при реализации численной процедуры следует ввести ряд фиктивных промежутков изменения случайной величины x, количество которых на единицу больше порядка разностных операторов, фигурирующих в модели. В этих фиктивных промежутках значения ВФР и ее производных равны нулю. Удобно считать эти граничные точки без ограничения общности точками . Тогда из (24) и (25) получаем интегрированием по частям эмпирическое уравнение эволюции первого момента ВФР в виде

                                       .                                (26)

Это уравнение означает, что изменение со временем среднего значения ряда равно средней скорости изменения значений ряда по скользящей выборке, что совпадает с представлениями математической статистики в дискретном случае. Таким образом, уравнение Лиувилля (23) является адекватной моделью для описания эволюции первого момента ВФР.

Рассмотрим эволюцию центральных выборочных моментов, определяемых в соответствии с (20) как

                             .                                        (27)

Скорость изменения этих величин также может быть найдена по уравнению (23). Имеем

Преобразуем полученное выражение. Первые два слагаемых в нем, содержащие , в силу (26) взаимно уничтожаются, а оставшийся интеграл может быть выражен через ковариации скорости и степеней исходной случайной величины x:

                 .

Поскольку при  одна из величин под знаком ковариации является константой, то соответствующая ковариация равна нулю. Таким образом, окончательно получаем уравнение

                                                          (28)

В частности, для уравнения эволюции выборочной дисперсии получаем

                                                .                                                     (29)

С другой стороны, изменение тех же величин (27) может быть вычислено в момент времени t непосредственно по имеющимся выборкам (см. выше Рис. 1-3). Получающиеся при этом значения будут, вообще говоря, отличаться от (28), поскольку эмпирическое уравнение Лиувилля (23) без включения в рассмотрение высших разностей (ускорений и т.д.) представляет лишь приближенную эволюционную модель ВФР. Мы сознательно пошли на такие упрощения, желая сократить объем вычислений, т.к. величина массива данных растет как степень объема одномерной выборки. Возникает вопрос: можно ли, не увеличивая существенно объем вычислений, т.е. оставаясь в том же пространстве, уточнить модель (23), чтобы некоторые, по крайней мере первые центральные моменты ВФР, действительно эволюционировали бы в соответствии с уточненным уравнением?

Для ответа на этот вопрос используем эвристические соображения о том, что при надлежащих условиях на моменты некоторого стационарного распределения оно (распределение) определяется ими однозначно (теорема Гамбургера [7]). Именно, теорема Гамбургера утверждает, что если некоторая последовательность чисел  такова, что ряд по обратным степеням этих величин расходится

                                               ,

то существует единственная функция распределения , моменты которой равны этим числам .

Как было показано в п.1, с некоторой точностью можно говорить об -стационарных выборочных распределениях, которые хотя и не имеют пределов (по вероятности или в слабом смысле и т.п.), могут доставлять оценки тех или иных статистик, удовлетворяющие практические нужды.

Поэтому, хотя требование состоятельности оценки, например, асимптотически несмещенной дисперсии не выполняется, желательно, чтобы с точностью, контролируемой исследователем, соответствующим выборочным средним по нестационарным распределениям можно было придать практический смысл. Полезно, например, ввести понятие -несмещенной оценки по выборке объема T. Такой оценкой мы будем называть статистику, отклонение среднего по которой от среднего по предположительно стационарному распределению не превосходит . В этом случае отличие модулей разности смещенных и несмещенных оценок, построенных по двум близким (т.е. -стационарным) ВФР будет иметь порядок .

Определение 7. Оценка  величины a по -стационарному распределению называется -несмещенной, если отличие модулей разности этой оценки и оценки  той же величины, но являющейся несмещенной, если распределение стационарно, имеет порядок .

Например, пусть A есть оценка дисперсии  согласно определениям (19-20). Если бы распределение было стационарным, то величина  являлась бы несмещенной оценкой дисперсии этого распределения. Поскольку моменты -стационарного распределения тоже -стационарны, то . Так как оценки B и A связаны равенством  и, по оценке (6), , то для  получаем оценку . Таким образом, выборочная дисперсия (21) и аналогично высшие центральные моменты являются -несмещенными оценками -стационарного распределения.

Поскольку выборочный момент порядка k имеет k-1 степень свободы, то для выборки фиксированного объема T осмысленно можно ставить вопрос об исследовании только первых T-1 моментов. Более того, чтобы получающиеся оценки были -несмещенными, требуется выполнение условия . Это означает, что корректными оценками моментов нестационарного распределения конечной выборки являются только первые несколько выборочных моментов. Именно, если при прогнозировании на  шагов при заданном уровне -стационарности ВФР определен оптимальный объем выборки , то при построении модели эволюции ВФР мы можем аргументированно использовать оценки только для моментов порядка до .

Приведенные соображения позволяют уточнить эволюционную модель для ВФР. Например, пусть у нас есть эмпирическое значение производной по времени дисперсии данной ВФР, равное . Обозначим невязку этого значения с уравнением (29) через :

                                           .                                             (30)

Рассмотрим для ВФР вместо уравнения Лиувилля (23) уравнение Фоккера-Планка:

                                        .                                             (31)

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

Уравнение Фоккера-Планка применяется, как правило, для описания процесса с независимыми нормально распределенными приращениями. При независящих от времени величинах  и  оно описывает перенос и диффузию вероятности для гауссовского белого шума. Для нестационарного процесса оно играет роль модельного и является определенным уточнением уравнения Лиувилля, если есть основания считать, что в исследуемом процессе присутствуют случайные блуждания. Некоторым «оправданием» такой не вполне корректно обоснованной модели служит то, что уравнение Фоккера-Планка активно (и также не вполне обоснованно) применяется во многих задачах статистической физики. Необоснованность же уравнения связана с тем, что гауссов белый шум является удобной, но абстрактной математической моделью, на практике не реализующейся, поскольку не может реально наблюдаться процесс, имеющий автокорреляционную функцию в виде дельта-функции.

С целью обобщения уравнения (31) при построении модели эволюции ВФР, которое должно приводить к верным эмпирическим значениям производной по времени для выборочных моментов вплоть до порядка r, вводим эмпирически определяемые величины

                                            (32)

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

                            .                                  (33)

Легко проверяется, что в силу этого уравнения .

Важно подчеркнуть, что количество учитываемых эмпирических выборочных моментов в (33) зависит от горизонта и желаемой точности прогноза. Если оказывается, что при этом близость между ВФР должна быть не больше заданного числа , то, как указывалось выше, учитываемое количество моментов не должно превосходить , которое, в свою очередь, не превосходит 2t. Например, при e~(1/T; 2/T) корректные модели эволюции ВФР могут быть основаны на уравнениях Лиувилля (23) или Фоккера-Планка (31). Если e~(2/T; 4/T), то можно использовать уравнение (31) или включить в рассмотрение третью пространственную производную, и т.д.

4. Сравнительный анализ точности прогнозных моделей

При численном решении уравнений (23) и (31) производные по аргументу x аппроксимировались схемами с центральными разностями. Для модели, основанной на уравнении Лиувилля, численно решалось уравнение

         .           (34)

Для уравнения Фоккера-Планка дискретная схема определения ВФР на следующем временном слое имеет вид

                   (35)

При решении конечно-разностного уравнения Лиувилля (34) на каждом шаге по времени проводилась коррекция найденной функции распределения. Отрицательные значения «насильственно» занулялись, а затем распределение вновь нормировалось на единицу.

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

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

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

Аналогичные действия выполнялись и при численном решении уравнения Фоккера-Планка.

Прогноз ВФР на основе численного расчета по уравнениям (34) и (35) показал, что модель на основе уравнения Фоккера-Планка является более точной. Хотя в абсолютных величинах повышение точности незначительное – в среднем ошибка прогноза по уравнению (35) меньше ошибки прогноза по уравнению (34) на слабо меняющуюся от объема выборки величину, равную 0,21%, – в относительных величинах (т.е. по отношению к ошибке прогноза) это может быть заметным эффектом. Например, при прогнозировании ВФР с точностью 5% относительный эффект использования более точной модели составил приблизительно 4%, тогда как при прогнозировании с точностью 1,5% – уже 14%.

Фактическая ошибка прогноза ряда оказывается несколько больше, чем декларируемая точность прогноза ВФР. Это связано с тем, что кроме чисто статистической ошибки прогноза погрешность содержится и в самом методе. Так, при построении ВФР ошибка в определении принадлежности элементов ряда некоторому интервалу  имеет значение порядка , где  есть количество ячеек, на которое разбит интервал изменения величины . В наших расчетах использовалось разбиение на 100 ячеек, т.е. ВФР строилась с точностью 1%. Кроме того, ошибка возникает и при процедуре перенормировки ВФР после каждого шага расчета по времени. Оказалось, что из-за появления отрицательных значений нормировка не сохраняется примерно на 0,5%. Наконец, неточность вносится и тем, что скорость , дающая прогноз изменения ВФР, берется в один – предшествующий – момент времени. На самом деле эта скорость также является случайной величиной, которая должна прогнозироваться на основе анализа ускорений и т.д., однако увеличение размерности фазового пространства привело бы к существенному замедлению счета.

Отметим, что точность конструирования ВФР по данным наблюдения, которые считаются «абсолютно точными», может быть увеличена единственным образом – путем использования более мелкого разбиения. Естественно, эта мелкость является ограничивающим фактором для самой постановки вопроса о прогнозе ВФР с заданной точностью, поскольку точность прогноза не может превышать точность конструирования самой ВФР. Однако и сама мелкость разбиения не может быть выбрана произвольной. Дело в том, что эмпирическая вероятность попадания значения случайной величины в некоторый интервал является по своему существу величиной стационарной, т.к. она определяется по одной конкретной выборке. Для нестационарных процессов такая выборка единственна, т.к. у нас нет возможности построить ансамбль вместо данной реализации случайного процесса. Поэтому, чтобы можно было корректно говорить о частотах попадания значений элементов выборки в заданный интервал как о квазистационарных вероятностях, надо, чтобы число элементов в выборке было во всяком случае не меньше, чем число ячеек разбиения. Для выборок малого объема это условие заведомо нарушается, т.е. ВФР не всегда статистически репрезентативна. Поскольку закон распределения вероятностей не известен априори, то под него невозможно заранее подогнать разбиение области значений случайной величины, чтобы эти значения заполняли бы ячейки равномерно. Приемлемость такого подхода связана с тем, что для нестационарных процессов частоты попадания значений в некоторый интервал все равно нельзя трактовать как оценку соответствующей вероятности, а возможность использования их в этом качестве обеспечивается функционалом (5) как мерой близости двух распределений.

Подчеркнем, что кинетическое уравнение относительно ВФР как функции только одного фазового аргумента x является прогнозной моделью первого порядка, когда по «известной» скорости строится прогноз самого ряда. Поэтому корректно сравнивать этот прогноз с другими – стандартными – методами также первого порядка: прогнозами методом скользящих средних для регрессионной и авторегрессионной моделей. Для получения модели второго порядка в уравнении Лиувилля (21) следует рассмотреть фазовое пространство координат, скоростей и ускорений. Строго говоря, уравнение Фоккера-Планка неявно учитывает зависимость от высших производных, поскольку в нем используется некоторая «поправленная» локальная средняя скорость, однако эта модель является все же менее точной, чем уравнение Лиувилля в фазовом пространстве большей размерности.

В качестве тестовых примеров оценки точности метода были взяты ценовые ряды акций компаний GM, GE, MSFT, FORD по данным [6].

По результатам анализа можно сделать два основных вывода. Во-первых, среди стандартных прогнозных методов первого порядка наша модель имеет более высокую точность (хотя и не очень значительно). Во-вторых, ошибка прогнозов по стандартным моделям скользящих средних минимальна, если в качестве объема усреднения берется оптимальный объем, вычисленный по горизонтной статистике для данного ряда. Например, для ряда GE (Рис. 4) ошибка по оптимальному и неоптимальному объемам выборки в среднем составила для стандартных методов соответственно 10-11% (показаны на Рис. 4) и 17-19% (не показаны), а для прогноза по кинетическому уравнению – 9%. Корреляционные методы имеют более высокую точность, чем регрессионные. Аналогичное поведение демонстрируют и ряды для GM, MSFT, FORD.

Разумеется, регрессионные методы первого порядка не претендуют на высокую точность в описании даже стационарных процессов, что и демонстрируется графиками на Рис. 4. Однако следует отметить, что кинетический прогноз, представленный гладкой кривой, имеет более высокую адаптивность, хотя и построен по тем же самым скользящим выборкам. Кроме того, важным его преимуществом является то, что в регрессионных методах ошибка прогноза накапливается с ростом горизонта прогнозирования, тогда как точность кинетического прогноза остается в фиксированных пределах.

Рис. 4. Сравнение прогнозов ряда GE на 1 шаг вперед

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

5. Уравнение Лиувилля для динамических систем с особенностями

Используемое в настоящей работе эволюционное уравнение для ВФР является не вполне законным с той точки зрения, что формально у нас нет динамической системы, в силу которой эволюционирует распределение начальных условий, т.е. распределение ансамбля таких динамических систем. Кроме того, нет и самого ансамбля. В качестве псевдо-скоростей, псевдо-ускорений и т.п. мы использовали выборочные эмпирические зависимости по текущему распределению, которые носят лишь формальную аналогию с величинами, входящими в уравнение Лиувилля статистической механики. В данной работе мы не обсуждаем метод выделения т.н. главных компонент ряда [8], который с некоторой точностью позволяет подобрать динамическую систему, дающую похожую псевдо-случайную траекторию. Термин «уравнение Лиувилля» используется нами здесь в том смысле, что это уравнение эволюционного типа, которое по виду ВФР в предыдущий момент времени позволяет найти ее в следующий момент. В силу вероятностного характера динамики, порождаемой случайным временным рядом, эта эволюция должна пониматься в некотором усредненном аспекте. Такой моделью и является предлагаемое в работе уравнение (25) и его дальнейшие модификации (33).

В связи с использованием уравнения Лиувилля для прогнозирования ВФР необходимо отметить еще один теоретический аспект, возникающий при использовании этого уравнения «по прямому назначению», т.е. для динамических систем – но с хаотическим поведением. При этом хаос в динамических системах может быть не только детерминированным, т.е. порождаемым специальной схемой дискретизации некоторой непрерывной регулярной динамической системы (как в случае логистического отображения), но и появляться в результате вырожденности траектории в фазовом пространстве, т.е. в силу неединственности решений уравнений движения. Если нет указания на поведение системы в особой точке, то она может с некоторой вероятностью продолжить движение по любой из доступных ей траекторий. Статистические свойства таких систем изучались в [9, 10].

Поскольку мы строим эмпирическое уравнение Лиувилля, которому по смыслу отвечает некоторая динамическая система с хаотическим поведением, следует предусмотреть возможность того, что, получаемая численно, такая система может быть дискретной реализацией динамической системы с особыми точками в фазовом пространстве. Для такой системы уравнение Лиувилля должно быть модифицировано введением дополнительных, не содержащихся изначально в системе условий, позволяющих «пройти» особое решение.

В качестве простейшей иллюстрации проблемы построения уравнения Лиувилля для динамических систем с особенностями рассмотрим систему, порождаемую дифференциальным уравнением, не разрешенным относительно производной по времени:

                                               .                                              (36)

Достаточным условием разрешимости уравнения (36) относительно p является . В этом случае через каждую такую точку  в расширенном фазовом пространстве, называемую регулярной, проходит ровно одна интегральная кривая. В точке  определено единственное векторное поле , где  является локальным решением уравнения (36) относительно . Исследуем эволюцию ансамбля таких систем, различающихся начальными условиями, и введем плотность распределения координаты x, учитывая, что эта координата является решением уравнения (36). Тогда эволюция распределения  в силу системы (36) описывается уравнением Лиувилля

                                                   .                                             (37)

Т.к. полная вероятность для системы пребывать в фазовом пространстве  равна единице и сохраняется во времени, то должна существовать плотность динамически-инвариантной меры, такой, что , где  и есть эта плотность. Легко убедиться, что  есть требуемая плотность, т.к. тогда

                           ,

и уравнение (39) можно записать в дивергентной форме закона сохранения.

В случае же, если , такой меры в общем случае может и не быть, а уравнение вида (37) записать не удастся, т.к. скорость p может быть не определена. Если , то уравнение  можно локально разрешить относительно p, после чего исключить p и из исходного уравнения (36). Получившаяся в результате дискриминантная кривая  в расширенном фазовом пространстве  является особым решением, имеющим общую точку с каждой регулярной интегральной кривой. Движение ансамбля систем в этой точке не определено, т.к. в ней нет единственности решения уравнения (36). Следовательно, на дискриминантной кривой следует поставить граничные условия для плотности вероятности , которые могут быть сформулированы в виде некоторого правила рассеяния вероятности на особой траектории. Эти правила можно записать в том же виде, что и для движения по регулярной траектории:

                                           ,                                   (38)

где  есть скорость перехода (вероятность перейти из точки  в точку x в единицу времени). Для регулярной точки траектории скорость перехода  определяется через фазовый поток , который представляет собой символическую запись интегральных кривых векторного поля  для уравнения (36): . Именно, в этом случае . Для сингулярной точки выражение  содержит фактор  и правила продолжения траектории через эту точку. Важно подчеркнуть, что эти правила являются внешними по отношению к задаче (36), они не выводимы из этого уравнения. Например, можно потребовать, чтобы фазовый поток «не замечал» особую траекторию. Тогда, если обозначить  неявную функцию, определяемую дискриминантной кривой, получим правило

                                       .

После этого  доопределяется по непрерывности во всех точках.

Характерным примером, обобщающим случай (36), является движение лагранжевой системы с вырождающимся гессианом [10]. Рассмотрим систему, динамика которой задается лагранжианом . Если действие

имеет стационарную точку, то эта точка называется экстремалью и является решением уравнений Эйлера-Лагранжа

                                                   ,

которые представляют собой систему дифференциальных уравнений второго порядка по времени  относительно координат:

                               .                               (39)

Совокупность величин  образуют тензор массы динамической системы. Обозначим определитель соответствующей матрицы Гесса через . Согласно определениям [11], лагранжиан  называется невырожденным, если этот определитель не обращается в ноль ни в одной точке фазового пространства. Лагранжиан называется сильно невырожденным, если он невырожденный, и уравнения  имеют единственные непрерывно дифференцируемые и обратимые решения . Отображение  называется преобразованием Лежандра, если  сильно невырожден. Если в некоторой области фазового пространства  гессиан вырожден, то там ускорения могут не быть представлены как однозначные функции фазовых переменных. В работе [9] эта ситуация трактуется как хаотизация классической динамики.

Классическим примером является случай, когда лагранжиан зависит от модулей координаты и скорости: . Тогда

                                             ,                         

где  – энергия системы, а  и  – взаимно-ортогональные проекторы на направление вектора скорости в n-мерном пространстве:

                                .

Определитель матрицы  легко вычисляется, поскольку  – ее собственное значение кратности 1, а  – собственное значение кратности , так что

             .              (40)

Отсюда следует, что особыми точками являются нетривиальные (с ненулевым значением скорости) стационарные точки энергии и лагранжиана, причем последние существуют только в пространстве размерности больше 1. Стационарные точки энергии отвечают разрыву компоненты ускорения, направленной вдоль скорости и потому касательной к пространственной траектории тела, а стационарные точки лагранжиана – разрыву составляющей ускорения, ортогональной скорости. Согласно [10], сингулярным множеством  такой динамической системы называется множество точек, в которых . В данном примере  разбивается в прямую сумму двух множеств: сингулярного множества  первого типа, отвечающего разрыву компоненты ускорения, направленной вдоль движения, и сингулярного множества  второго типа, отвечающего разрыву нормальной компоненты ускорения.

Пусть система имеет набор независимых первых интегралов . В [10] показано, что траекторию в некоторых случаях можно продолжить через множество  с сохранением по непрерывности первых интегралов и зеркальным отражением направления движения, если это точка локального максимума, либо прохождением без изменения направления движения, если это точка локального минимума. В частности, при отражении

точка сохраняет координаты в физическом пространстве, но меняет на противоположные нормальную составляющую скорости. При прохождении без отражения и «конденсации» на особом множестве получаем

             .

В случае с сингулярным множеством  второго типа n – 1 компонент ускорения являются неопределенными. Это можно трактовать как рассеяние на произвольный угол с сохранением интеграла энергии.

В литературе весьма подробно исследованы вопросы возникновения хаоса в динамических системах (см., напр., монографии [12, 13]), а также изучены и классифицированы особенности динамических систем в виде дифференцируемых отображений [14-16]. В первых системах хаос возникает при определенных параметрах дискретизации дифференциальных уравнений движения, в которых не было каких-либо особенностей. Во вторых системах особенности присущи именно дифференциальному уравнению, а при дискретизации они проявляются лишь приближенно и при достаточно мелком шаге разбиения. В последнем случае уравнение Лиувилля имеет некоторые черты «эмпирического», т.к. оно дополняется некоторыми внешними условиями, обеспечивающими проход через особое решение, возможно, в виде вероятностного выбора каждой из траекторий. Тем самым статистическая механика динамических систем с особенностями образует связь с теорией случайных процессов и, в частности, временных рядов.

Заключение

Итак, в работе предложен метод анализа нестационарных временных рядов, который направлен на решение следующих задач:

- определение оптимального объема выборки для формирования квазистационарной ВФР и ее прогнозирования с заданной точностью;

- минимизация совокупного функционала ошибки аппроксимации временного ряда на заданном интервале времени, а также определение интервала времени, на котором такая ошибка минимальна;

- прогнозирование ВФР и временного ряда внутри границ квазистационарности ВФР.

Поскольку в данном случае отсутствуют классические понятия доверительных интервалов, то прогноз в нашем подходе означает совокупность вероятного значения ряда в наперед заданный момент времени и распределения его значений, выборочные квантили которого и дадут оценку доверительных интервалов. В работе показано, что для целей прогнозирования нестационарных временных рядов с использованием понятия эволюции ВФР могут быть применены кинетические уравнения классической статистической механики – уравнение Лиувилля и Фоккера-Планка.

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

Благодарности

Авторы выражают благодарность своим коллегам профессорам Г.Г. Амосову и Г.Г. Малинецкому, а также доцентам Н.А. Митину и В.Ж. Сакбаеву за плодотворные обсуждения и ценные указания.

Цитированные источники

1. Королюк В.С., Портенко Н.И., Скороход А.В., Турбин А.Ф. Справочник по теории вероятностей и математической статистике. – М.: Наука, 1985.

2. Уилкс С. Математическая статистика. – М.: Наука, 1967.

3. Орлов Ю.Н., Осминин К.П. Анализ нестационарных временных рядов. / Препринт ИПМ им. М.В. Келдыша РАН, № 36, 2007. – 24 с.

4. Орлов Ю.Н., Осминин К.П. Построение выборочной функции распределения для прогнозирования нестационарного временного ряда. // Мат. Мод., 2008, № 7 (в печати).

5. Осминин К.П. Алгоритмы построения статистик для анализа и прогнозирования нестационарных временных рядов. // ИТВС, № 3, 2008 (в печати).

6.http://moneycentral.msn.com/detail/stock_quote?Symbol=GE&pkw=PI&vendor=Paid+Inclusion&OCID=iSEMPI

7. Ахиезер Н.И. Классическая проблема моментов. – М.: Физматлит, 1961.

8. Главные компоненты временных рядов. Сб. статей / Ред. Д.Л. Данилов и А.А. Жиглявский. СПбГУ, 1997.

9. Pavlotsky I.P., Strianese M., Toscano R. Prolongation of the integral curve on the singular set via the first integral. // J. of Interdisciplinary Math. 1999. V.2. Nos.2-3. P.101-119.

10. Орлов Ю.Н. Основы квантования вырожденных динамических систем. – М.: МФТИ, 2004. – 236 с.

11. Дубровин Б.А., Новиков С.П., Фоменко А.Т. Современная геометрия. – М.: Наука, 1986. – 759 с.

12. Кузнецов С.П. Динамический хаос. – М.: Физматлит, 2001. – 296 с.

13. Малинецкий Г.Г., Потапов А.Б. Современные проблемы нелинейной динамики. – Москва-Ижевск: Регулярная и хаотическая динамика, 2000.

14. Каток А.Б., Хасселблат Б. Введение в современную теорию динамических систем. – М.: Факториал, 1999. – 767 с.

15. Арнольд В.И., Ильяшенко Ю.С. Обыкновенные дифференциальные уравнения. / Итоги науки и техники ВИНИТИ. – М., 1985.

16. Арнольд В.И., Варченко А.Н., Гусейн-Заде С.М. Особенности дифференцируемых отображений. Т.1. – М.: Наука, 1982.