Энтропия и прогноз временных рядов в теории динамических систем

ЭНТРОПИЯ И ПРОГНОЗ ВРЕМЕННЫХ РЯДОВ В ТЕОРИИ ДИНАМИЧЕСКИХ СИСТЕМ

AM. Лоскутов, А.А. Козлов, Ю.М. Хаханов

В работе дается современный взгляд на такие понятия, как размерность и энтропия динамических систем. Описание данных характеристик включает в рассмотрение другие представления и свойства, относящиеся к сложному поведению нелинейных систем -размерность вложения, горизонт предсказуемости и др., которые также используются в работе. Изучается вопрос о возможности применения этих концепций к реальным наблюдаемым экономического происхождения - ценам акций компаний Schlumberger, Deutsche Bank, Honda, Toyota, Starbucks, BP. С помощью метода сингулярного спектрального анализа дается прогноз цен акций в различных фазах экономического цикла - предкризисный и кризисный периоды. Определены основные ограничения использования этого метода.

Ключевые слова: Энтропия, временные ряды, хаос, размерность.

Введение

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

В экономической теории многие финансовые показатели ввиду сложности описания рыночных систем часто рассматриваются как случайные величины. Поэтому финансовые временные ряды традиционно анализируются с помощью статистических методов и моделей. Среди них наиболее известны авторегрессия и скользящее среднее. Такой подход очень распространен, но он не дает ключ к пониманию природы изменчивости рынков капитала и вопроса о предсказании колебания цен. Для этой задачи более успешно, по-видимому, может быть использована теория динамических систем [1], также именуемая теорией хаоса [2-4] или нелинейной

динамикой [5]. Эти исследования насчитывают историю значительно более короткую, нежели статистический анализ, особенно в применении к экономической теории [6,7]. Однако сейчас можно надеяться, что методы нелинейной динамики в сочетании с уже выполненными экономическими исследованиями приведут к созданию новых моделей рынков капитала.

Необходимо также отметить, что в настоящее время активно развивается такой раздел математической физики, как анализ временных рядов [8,9], вобравший в себя методы теории динамических систем и фрактального анализа [10]. Подобный анализ все чаще применяется в экономике для обработки информации в виде временных рядов котировок акций компаний. Это связано с тем, что рынок ценных бумаг нестабилен, так что предсказать динамику котировок чрезвычайно сложно. Однако для хаотических динамических систем это типичная ситуация. Экспоненциальное разбегание близких траекторий для класса хаотических аттракторов определяется старшим ляпуновским показателем, а скорость «расползания» элемента фазового объема - энтропией. Вычислив любую из этих величин, можно получить значение горизонта прогноза для данного временного ряда и/или поведения системы в целом. Кроме того, методы теории динамических систем, применяемые к финансовым временным рядам, позволяют найти размерность вложения, то есть минимальное число динамических переменных, однозначно описывающих наблюдаемый процесс. Если размерность вложения конечна, это указывает на конечномерность рассматриваемого процесса и дает возможность определить число главных факторов на исследуемом сегменте рынка.

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

1. Нелинейная динамика и теория временных рядов

Скалярным временным рядом называется массив из N чисел, представляющих собой значения некоторой измеренной (наблюдаемой) динамической переменной x(t) с некоторым постоянным шагом At по времени, то есть взятые в моменты ti = to + (i — 1)At : Xi = x(ti), i = 1,..., N. В теории временных рядов возможны два различных случая, связанных с наличием информации о порождающей их системе. Если математические модели, описывающие систему, известны (например, в астрономии), то обработка временных рядов представляет собой нахождение параметров, входящих в известные нам уравнения. Во втором случае возникают принципиальные трудности, так как математический объект, который следует поставить в соответствие полученным данным (например, финансовым временным рядам), неизвестен. Таким образом, выделяются две основные задачи анализа временных рядов: задача идентификации и задача прогноза.

Задача идентификации предполагает ответ на вопрос, каковы параметры системы, породившей данный временной ряд - размерность вложения, корреляционная размерность, энтропия и др. Рассмотрим смысл каждой искомой величины.

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

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

tp — K1"«?

где е - относительная погрешность данных о ее состоянии. Вторая трактовка связана с понятием показателей Ляпунова -j. В этом случае также можно ввести время предсказуемости

T- — --log-.

—max е

По порядку величины обе эти оценки совпадают. На практике, однако, значение времени предсказуемости считают по несколько упрощенным формулам: Tp — K-1 и T- — -max-1. В данной работе мы будем использовать первую из них.

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

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

• сингулярный спектральный анализ - глобальный метод [11];

• локальная аппроксимация - локальный метод [12];

• локальная аппроксимация с использованием сингулярного спектрального анализа - глобально-локальный подход [13].

Исследование временных рядов базируется на идее реконструкции аттрактора динамической системы, породившей данный временной ряд. Дело в том, что удовлетворительную геометрическую картину странного аттрактора можно получить, если вместо переменных Xj, входящих в уравнения динамической системы x = f (x), использовать m-мерные векторы следующего вида (то есть векторы задержек) [14]:

zi — {xij xi+1j ...j xi+m-1}.

Такой подход к анализу временных рядов впервые был математически обоснован в работах Ф. Такенса [15]. Чтобы более наглядно представить содержание этой теории, дадим вначале некоторые пояснения.

Пусть динамическая система задана как ф*(х) с фазовым пространством Р. Величины, образующие временной ряд, являются значениями некоторой скалярной функции состояния динамической системы х(£) на многообразии М* из Р: хг = Л,(х(£г)) = (х0)). Рассмотрим построение векторов задержки z, если временной шаг равен т. Очевидно, для векторов состояния х(^) = хг верны равенства

Хг+1 = фТ(х»), хг+2 = ф2Т(хг), . . . , хг+т_1 = ф(т_1)т(хг). Поэтому для компонент z-векторов задержки получим

Хг = Н(х) = Фо(хг),

Хг+1 = й(хг+1) = ^(фт(хг)) = Ф1(хг),

Хг+2 = ^(хг+2) = Л,(ф2т(хг)) = Ф2(хг),

Хг+т_1 = Л,(хг+т_1) = Н(ф(т 1)т(хг)) = Ф^-^).

Таким образом, поскольку все компоненты вектора zг = {хг, хг-1,..., хг+т-1} можно связать с одним и тем же состоянием динамической системы хг, то существует векторная функция Л, отображающая векторы хг £ М* в векторы т-мерного евклидова пространства Мт

Zг = Л(хг), хг £ Мл, zг £ Мт

Приведенные рассуждения и составляют основное содержание теоремы Та-кенса [15], утверждающей, что при т ^ 2! + 1 типичным свойством отображения Л является вложение МЛ в Мт. Другими словами, образ МЛ в Мт (Б* = Л(Ма)) не имеет самопересечений. При этом функция Л дифференцируема и имеет обратную дифференцируемую функцию Л-1, определенную на Бл : Мл = Л-1(Б*). Тогда на Б* можно определить динамическую систему, так как хг = Л-1^г), а хг+1 = фт(хг),

^+1 = Л(хг+1) = Л(фТ(Л-1Ы)) = Ф^г), Zг £ Б*,

где функция Ф отображает Б* в Б*. Легко видеть, что, если в этом соотношении оставить последнюю компоненту, то получим другой вариант - покомпонентный -записи временного ряда: хг = Г(хг-1, ...,хг-т). Данное соотношение можно использовать для решения задачи прогнозирования, так как для этого требуется знать только величину т - размерность вложения.

Таким образом, мы имеем два отображения: хг+1 = фт(хг) = Ф(хг), х £ Ми Zг+l = Ф^г), z £ Бл. Так как z = Л(х), то это всего лишь различные представления одного и того же отображения. Поэтому характеристики, инвариантные относительно невырожденной замены, у обеих систем должны совпадать. Отсюда следует, что фрактальную размерность аттрактора, энтропию динамической системы и некоторые другие инвариантные характеристики можно определить по измеряемым величинам, то есть по соответствующим временным рядам.

2. Размерность вложения и метод корреляционного интеграла

Первыми работами, где была рассмотрена задача идентификации для дис-сипативных динамических систем, явились статьи П. Грассбергера и И. Прокаччиа [16,17]. Используя понятие корреляционного интеграла, они разработали метод оценки размерности вложения и фрактальной размерности аттрактора динамической системы по временному ряду некоторой наблюдаемой.

Будем регистрировать некоторую скалярную величину х через равные промежутки времени т. Тогда получим соответствующую последовательность длины N: хг=х(ь + (г — 1)т), г = 1,..., N. Выделим из этой последовательности все подпоследовательности длины т и построим набор из п = N — (т — 1) векторов задержек ъг = {хг, хг+1,..., хг+^т-1)}. Теперь определим расстояние рт(г, ]) для любой пары Zj, таких векторов

рт(г,]') — \\ъг || —

т— 1

(хг+1 — хз+г)2

I 1=0

1/2

Далее выберем некоторое значение г и рассчитаем корреляционный интеграл

С(т)(г) = 11ш п-2 V 0[г — рт(г,])]

г,3=1

где 0(х) - ступенчатая функция Хевисайда. Важное свойство корреляционного интеграла, отмеченное в работах [16,17] - его скейлинговое поведение при малых г: С (г) — г4, где V - корреляционная размерность. Эта величина тесно связана с фрактальной размерностью аттрактора й, которая может быть введена следующим образом.

Предположим, что нам нужно покрыть аттрактор (гипер)кубами с длиной стороны I. Так как в типичном случае аттрактор является фракталом, то число всех кубов, содержащих точки аттрактора, будет порядка М(I) — 1-а. Пусть (г = 1, 2,..., М(1)) - число точек из набора (] = 1, 2,..., п), которые лежат в г-м кубе. Тогда можно записать

м

С (о - ^Е у

М

г=1

п2

(У2).

Используя неравенство Шварца (у2) ^ (у)2, получаем

с«) > М21 «2=

п2 • М

'м

г=1

М

Отсюда, поскольку С(I) — Г и I ^ 1, получим: V ^ й. Таким образом, корреляционная размерность всегда не превышает значение фрактальной размерности аттрактора. Однако на практике было замечено, что значения V и й часто очень близки, особенно когда точки ъг распределены по аттрактору практически равномерно. В этом случае метод корреляционного интеграла дает хорошую оценку для фрактальной

2

1

1

а

I

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

Вернемся к описанию метода. Согласно определению, корреляционный интеграл С(т)(г) учитывает только те пары векторов Zг и Zj, расстояние между которыми меньше или равно г. Поэтому при малых г имеется тенденция С(т)(г) ^ 0 по закону С (т)(г) ~ гУт. Таким образом, величины ут можно найти, построив графики зависимостей 1п С(т)(г) от 1п г для каждого т, начиная с т = 1. Тогда с ростом т величина ут, начиная с некоторого номера, перестает возрастать. Это значение т дает размерность вложения !е, а предельный тангенс угла наклона ут - корреляционную размерность !с аттрактора исследуемой динамической системы.

3. Общие подходы к определению энтропии и размерностей динамических систем

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

Н = рг ^ рг .

г

Похожим образом можно дать определение метрической энтропии динамической системы или энтропии Колмогорова-Синая. Благодаря этому понятию удалось строго сформулировать абсолютный критерий хаотической динамики системы как неустойчивого по Ляпунову движения с положительной метрической энтропией.

Пусть динамическая система представлена отображением хп+1 = £(хп) на компактном множестве А. Пусть также задано разбиение А на конечное число измеримых множеств Аг. Будем обозначать множество точек, отображаемых в Аг преобразованием £к, как £_к(Аг). Не выписывая все разбиения, порождаемые таким обратным отображением, опишем общий случай: на т-м шаге разбиение

на A™; = An П f-l(Ai2) П f-2(Агз) n ... П f-(m-1)(Aim) -

Hli2...im %1 n П * \Aii'3

это те точки из

множества A;1, которые на следующих (m — 1) шагах попадут последовательно в

Ai2 > Ai3 > • • • > Aim •

Вычислим энтропию каждого разбиения по аналогии с формулой Шеннона

H(m) = — £ P(Aiii2...irn )l0g p(Aiii2...im )•

ili2...im

Если обозначить e = maxi diamAi, то энтропия динамической системы определится как

K = lim lim

s—>0 m—ж s—>0 m—ж m

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

аппроксимирующая ее величина K2 через корреляционный интеграл. Эту величину можно ввести следующим образом.

Рассмотрим так называемую обобщенную энтропию (или энтропию Реньи) Hq:

н = I—q log Pi^j. о)

Энтропия Шеннона (1) является ее частным случаем при q — 1. Используя эту величину, можно также ввести обобщенную размерность Dq и величину Kq, которая называется обобщенной энтропией динамической системы, наиболее общим образом характеризующей меру неупорядоченности системы. Именно, если подставить в определение энтропии K (2) величину Hq вместо H, то

Kq = lim lim (H(m+1) - H(m)) = lim lim — H(m).

£—>0 т—ж V q q / £—>0 т—ж m q

Одно из важнейших свойств этого инварианта можно записать, как

K = Kl < K2 < ... < Kq-1 < Kq < ...,

где K - энтропия Колмогорова-Синая. Таким образом, энтропию Колмогорова-Синая можно оценить сверху, если найти Ki, i = 2, 3,....

Обобщенная размерность Dq определяется следующим образом:

D„ = Um (4 £—о\ ln е

Поясним смысл обобщенной размерности Dq для различных q. Если q = 0, то ppq = 1 и, следовательно, получим определение фрактальной размерности D0 = dp. Если q = 1, то числитель и знаменатель в соотношении (3) обращаются в нуль, поскольку Pi = 1. Однако нетрудно показать, что в этом случае

Di = lim £ piln pi

£—>0 1п £

Поскольку соотношение ^ рг 1п рг определяет энтропию, то величина показывает, как при £ ^ 0 возрастает количество информации о рассматриваемом множестве. Поэтому Б1 называется информационной размерностью. Если д = 2, то

^ = 1^^.

£—0 1п £

Выражение рг2 - это вероятность того, что по крайней мере две точки находятся в г-й ячейке. Следовательно, значение ^2 определяется через бинарные корреляции в распределении вероятностей. Поэтому величину Б2 называют корреляционной размерностью.

Для расчетов обобщенных величин вводится понятие обобщенного корреляционного интеграла

С. (£) = (£ р! Г"4

Так как -1п Ст)(е) = Н(т)(е), то Сд - это всего лишь переобозначение Нд, и связи с другими обобщенными величинами сохраняются. При этом обобщенный корреляционный интеграл удобен тем, что для случая ц = 2 его вид существенно упрощается

С2(е) = С(е) = ^ £ кг(е);

С2 называется корреляционным интегралом, введенным выше.

Вернемся к методу определения обобщенной размерности Од и энтропии Кд. Отмеченная связь —1п Сдт\е) = Ндт\е) позволяет выразить обобщенные величины через обобщенный корреляционный интеграл:

D = lim

e^o \ ln e J ln e

( ln C(m)(e)\ Kq = lim lim m-lH(m) = lim lim--.

Исходя из полученных соотношений можно ожидать, что при малых e и достаточно больших m выполняется

ln cim)(e) ^ Dg ln e — Kgm + const. Отсюда, кроме оценки Dg, фиксируя масштаб e и исследуя зависимость от окна реконструкции m, можно получить оценку обобщенной энтропии динамической системы Kg.

В работах [17, 18] была предложена иная схема расчета K2-энтропии, c ис-

ту гг(т+1) тт(т) п

пользованием представления Kg в виде предела разности Hg — Hg при e — 0

и m оо:

lim lim Km (e), K(m) = 1 ln C.

e^0 т^ж 2 W' 2 kT C

Здесь С (т)(е) - это корреляционный интеграл С2т)(е). Зависимость к2т) от размерности вложения т можно аппроксимировать по методу наименьших квадратов следующей функцией:

К2т) (е) = К2(е) + тс,

где Ь и с - некоторые положительные параметры, т - текущая размерность вложения.

Итак, в случае малого фиксированного е энтропия К2 определяется как предел К2т) при стремлении размерности вложения к бесконечности (т ^ ж).

4. Метод прогноза сингулярный спектральный анализ-«Гусеница»

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

которых позволяет судить об его особенностях и продолжить его на прогнозируемый отрезок времени. Метод ССА—«Гусеница» [11] относится к глобальным методам прогноза и используется для выделения периодических и квазипериодических составляющих из временного ряда. Кроме того, метод ССА можно применять для улучшения соотношения сигнал/шум.

Преобразуем исходный временной ряд {xi} x2,..., xn} во множество M-мерных векторов задержек zi = {x^ xi+i;...; xi+M-i}, которые запишем в виде матрицы задержек

XMх (N-M+i) = (zb z2,..., zN-M+i) =

( xi x2 ... xn-M+i \ x2 x3 ... xn-M+2

\ XM XM+i ... XN /

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

Этап нахождения главных компонент сводится к поиску собственных векторов и собственных значений матрицы X • Хт. Для этого находятся матрица собственных векторов Vм хм, в столбцах которой стоят координаты собственных векторов, и матрица собственных значений 2мхм = diag аг. При этом предполагаем, что собственные векторы упорядочены по убыванию соответствующих собственных значений, то есть а1 > а2 > ... > ам.

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

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

Аппроксимацию матрицы X по главным компонентам можно провести как

Х = Vмхг • ^хг • X,

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

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

Выпишем искомый вектор-столбец х^+1 = (хм-м+2;...; хм; Хм+1}. Условие минимума квадратичного отклонения искомого вектора от своей проекции на выбранную гиперплоскость главных компонент имеет вид

Оптимизация проводится по первой координате этого вектора хм+1. Таким образом можно построить прогноз следующего неизвестного значения ряда. Для дальнейшего прогнозирования достаточно добавить координаты найденного вектора хм+1 в матрицу X и снова рассчитать прогноз на следующий шаг. Кроме этого можно повторять весь алгоритм ССА заново, однако для больших рядов это нецелесообразно.

С помощью описанных методов был проведен анализ котировок акций следующих компаний: Toyota Motor Corp. (TM), Deutsche Bank, Honda Motor Co., Ltd., Schlumberger Limited, Starbucks, British Petroleum plc. Выбор этих компаний обусловлен тем, что их акции торгуются на фондовых биржах в течение длительного периода времени (не менее 10 лет), что предоставляет достаточное количество данных для их обработки. Ниже представлен подробный анализ только компании Toyota Motor Corp. Данные для остальных компаний сведены в таблице.

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

• 2007 год - так называемый предкризисный период;

• 2008 год - кризисный период.

Это соответствует информации о том, что очередной пик экономической активности в США был достигнут в декабре 2007 года [19].

Анализ достоверности уже построенного прогноза очень важен для оценки эффективности предложенного метода (ССА-«Гусеница»). Поэтому он проводился как по отклонениям предсказанных значений от реальных, так и по разнице соответствующих величин доходностей. Первая оценка предусматривает нахождение максимума модуля отклонения предсказанного значения от реального на горизонте прогноза. Вторая оценка, которая называется еще и логарифмической прибылью, вводится следующим образом [6,20]. Если x(ti) - цена акции в момент ti, x(ti-1) -цена акции в предыдущий момент ti_i, то доходность акции суть

2

^ mm .

XN+1

5. Результаты анализа реальных данных

Легко видеть, что для небольших изменений курса x(ti) — x(ti-l) = Ax(ti) доходность Si - это просто относительное изменение Ax(ti)/x(ti-l). Кроме того, знак Si, так же как и знак относительного изменения, показывает динамику котировок. Сравнивая знаки предсказанной и реальной доходностей, можно анализировать прогноз по динамике: чем больше совпадений доходностей по знаку, тем лучше прогноз.

Рассмотрим подробно результаты анализа временного ряда котировок акций компании Toyota Motor Corporation за период с 13.04.1993 по 31.12.2008 (3956 значений, см. рис. 1). Для выбранного временного ряда были проделаны все этапы метода корреляционного интеграла. При каждом расчете корреляционный интеграл считался для нескольких значений r при заранее введенном значении предполагаемой размерности вложения m (m = 1, 2,...). При этом для каждой кривой зависимости ln C (m)(l n r) необходимо было выделить линейный участок, на котором предстояло затем построить линейную аппроксимацию и найти корреляционную размерность vm. Окончательные результаты получались только после построения графика зависимости vm(m). Значение m, на котором достигалось насыщение, давало искомую размерность вложения de, а само значение насыщения vm - корреляционную размерность dc.

Полученные по методу Грассбергера-Прокаччиа результаты представлены на рис. 2. Из этих данных нетрудно определить размерность вложения и корреляционную размерность, соответственно: de = 6 и dc ~ 2.64.

1993-04-13 1997-01-21

2001-01-16

2005-01-11 2008-12-31

In с{тЛ

Рис. 1. Динамика цен акций Toyota Motor Corp. v„

а

-2.0 -1.0 0.0 1.0 ln г

т

ill vm

1 2 3 4 5 6 7 8 9 10 0.96785 1.63486 2.08106 2.32634 2.48358 2.6046 2.64294 2.65783 2.66748 2.63761

1 2 3 4 5 6 7

9 10 да

Рис. 2. Графики зависимостей 1п С(т\г) от 1пг при различных т (а) и ут от т (б) для временного ряда котировок ТМ

0 2 4 Рис. 3. График зависимости

(m) 2

от m для TM

Опишем теперь метод расчета энтропии K2. Для этого применяем методику Грассбергера-Прокаччиа, связанную с фиксированием величины r и рассмотрением зависимости от т. При увеличении т значения K2 уменьшаются до некоторой постоянной величины, равной энтропии K2. Анализ полученных результатов показан на рис. 3.

Отсюда следует, что энтропия K2 к к 0.1153. Следовательно, время предсказуемости котировок временного ряда, образованного ценами акций TM, составляет Tp ^ K2_1 к 8.67. Таким образом, рассматриваемый временной ряд теоретически можно прогнозировать на 8 точек (значений). Результаты такого прогноза представлены ниже.

Для прогноза мы использовали метод ССА с количеством главных компонент, равным размерности вложения de = 6. При этом общее количество компонент положим равным 20, то есть в 4 раза превосходящим ожидаемый период значимых колебаний на рынке (рабочая неделя - 5 дней).

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

В качестве примеров прогнозов в течение 2007 года были выбраны следующие отрезки временного ряда (после 24.04.2007 и после 12.09.2007), так как на них проявляются различные тенденции поведения котировок акций, что позволяет продемонстрировать возможности метода прогноза ССА-«Гусеница».

На рис. 4, а представлен прогноз, отличительной чертой которого является идентификация продолжения нисходящего тренда после небольшой коррекции вверх. При этом он имеет максимальное отклонение от реального ряда 0.85, а совпадение знаков доходностей происходит в 6 случаях из 8.

Второй прогноз в течение 2007 года представлен на рис. 4, б. На нем довольно четко прослеживается тенденция продолжающихся колебаний примерно на одном уровне после первоначального падения. При этом предсказанные значения максимально отклоняются от реальных на 1.92, а доходности совпадают по знаку в 6 из 8 случаев.

Прогнозы на кризисный 2008 год. В течение 2008 года для прогноза использовались значения временного ряда после 19.03.2008 и после 27.06.2008, так как на этих участках котировок акций проявляются возможности метода прогноза ССА-«Гусеница» по определению локальных экстремумов.

Рис. 4. Реальные данные и прогноз временного ряда котировок ТМ после 24.04.2007

Рис. 5. Прогноз для временного ряда ТМ после 19.03.2008

На рис. 5, а представлен прогноз, на котором точно определен момент достижения локального максимума значений временного ряда на горизонте прогноза. При этом наблюдаются серьезные отклонения реальных значений от прогноза, достигающие 1.91. Совпадение знаков доходностей наблюдается в 6 из 8 случаев. Второй прогноз на кризисный период представлен на рис. 5, б. На нем был зафиксирован конец длительного падения котировок и обозначен локальный минимум на горизонте прогноза. При этом спрогнозированные значения максимально отклоняются от реальных на 1.72, а совпадение знаков доходностей происходит в 5 случаях из 8.

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

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

Таблица

Компания (е (с К2 Т

Schlumberger Limited 7 3.12 0.1656 6.04

Deutsche Bank 6 3.26 0.1208 8.28

Honda Motor Co. Ltd. 7 2.80 0.1495 6.69

Toyota Motor Corp. 6 2.64 0.1153 8.67

Starbucks 7 2.34 0.1389 7.20

BP plc. 6 2.45 0.1403 7.13

Заключительные замечания

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

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

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

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

А.Ю. Лоскутов признателен профессору В.С. Анищенко за предложение написать данную работу и полезные дискуссии в течение недавно проходившей в СГУ международной конференции StatInfo-2009.

Библиографический список

1. Лоскутов А.Ю., Михайлов А.С. Основы теории сложных систем. РХД, 2007.

2. Кроновер Р.М. Фракталы и хаос в динамических системах. Основы теории. М.:

Постмаркет, 2000.

3. Кузнецов С.П. Динамический хаос, 2001.

4. Шустер Г.Г. Детерминированный хаос. М.: Мир, 1988.

5. МалинецкийГ.Г., Потапов А.Б. Современные проблемы нелинейной динамики.

М.: УРСС, 2002.

6. Петерс Э. Хаос и порядок на рынках капитала. Новый аналитический взгляд на циклы, цены и изменчивость рынка. М.: Мир, 2000.

7. Пу Т. Нелинейная экономическая динамика. УРСС, 2002.

8. Kantz H., Schreiber T. Nonlinear time series analysis. Cambridge University Press, 1997. 304 p.

9. Abarbanel H.D.I., Brown R., Sidorowich J.J., Tsimring L.S. The analysis of observed chaotic data in physical systems // Rev. Mod. Phys. 1993. № 65. P. 1331.

10. Mandelbrot B. The Fractal Geometry of Nature. W.H. Freeman and Company. New York, 1983.

11. Главные компоненты временных рядов: метод «Гусеница». Сб. статей / Ред. Д.Л. Данилов, А.А. Жиглявский. СПб: СПб университет, 1997. 308 с.

12. Лоскутов А.Ю., Журавлев Д.И., Котляров О.Л. Применение метода локальной аппроксимации для прогноза экономических показателей // Вопросы анализа и управления риском. 2003. Т. 1, № 1.

13. Истомин И.А., Котляров О.Л., Лоскутов А.Ю. К проблеме обработки временных рядов: расширение возможностей метода локальной аппроксимации посредством сингулярного спектрального анализа // Теоретическая и математическая физика. 2005. Т. 142, № 1.

14. Packard N.H., Crutchfield J.P., Farmer J.D., Shaw R.S. Geometry from a time series // Phys. Rev. Lett. 1980. Vol. 45. P. 712.

15. Takens F. Dynamical systems and turbulence // Lect. Notes in Math, Berlin, Springer. 1981. № 898. P. 336.

16. Grassberger P., Procaccia I. Characterization of strange attractors // Phys. Rev. Lett. 1983. Vol. 50. P. 346.

17. Grassberger P., Procaccia I. Measuring the strangeness of strange attractors // Physica D. 1983. Vol. 9. P. 189.

18. Grassberger P., Procaccia I. Estimation of the Kolmogorov Entropy from a chaotic signal // Phys. Rev. A. 1983. Vol. 28, 4. P. 2591.

19. Илларионов А. Досрочная рецессия. Smart Money. 2008, № 46 (136).

20. Романовский М.Ю., Романовский Ю.М. Введение в эконофизику. Статистические и динамические модели. М.; Ижевск: Институт компьютерных исследований, 2007.

Поступила в редакцию 8.07.2009

ENTROPY AND FORECASTING OF TIME SERIES IN THE THEORY OF DYNAMICAL SYSTEMS

A.Yu. Loskutov, A.A. Kozlov and Yu.M. Khakhanov

A contemporary consideration of such concepts as dimension and entropy of dynamical systems is given. Description of these characteristics includes into the analysis the other notions and properties related to complicated behavior of nonlinear systems as

embedding dimension, prediction horizon etc., which are used in the paper. A question concerning the application of these ideas to real observables of the economical origin, i.e. market prices of the companies Schlumberger, Deutsche Bank, Honda, Toyota, Starbucks, BP is studied. By means of the method of singular spectrum analysis the forecasting of the market prices of these companies in different phases of the economical cycle - just before crisis and during the crisis - is given. Main advantages and demerits of the method used are found.

Keywords: Entropy, time series, chaos, dimension.

Лоскутов Александр Юрьевич - доктор физико-математических наук, профессор физического факультета МГУ. Область профессиональных интересов -теория динамических систем и ее приложения. Имеет более 130 публикаций в научных журналах, 2 учебника, 4 монографии (2 из которых с соавторами). Лауреат премии им. И.И. Шувалова по математике и физике (1998). Член редколлегий журналов «Проблемы анализа риска», «Вестник МГУ, серия физика-астрономия», «Nonlinear Biomedical Physics». Приглашенный профессор многих ведущих университетов мира. 119992 Москва, Воробьевы горы

Московский государственный университет, физический факультет E-mail: loskutov@polly.phys.msu.ru

Козлов Александр Александрович - родился в 1986 году в г. Нововоронеж Воронежской области. В 2009 году окончил физический факультет Московского государственного университета имени М.В. Ломоносова, получив диплом магистра физики с отличием. 119992 Москва, Воробьевы горы

Московский государственный университет, физический факультет E-mail: kozlov@polly.phys.msu.ru

Хаханов Юрий Михайлович - родился в 1986 году в Москве. В 2009 году окончил физический факультет Московского государственного университета имени М.В. Ломоносова, получив диплом магистра физики с отличием. 119992 Москва, Воробьевы горы

Московский государственный университет, физический факультет E-mail: khakhanov@polly.phys.msu.ru