Thursday, June 16, 2011

Фильтр Калмана, руководство полного идиота


Фильтру Калмана посвящено множество статей и книг, совсем недавно была опубликована статья на хабре. Однако, чтение подобных статей у меня всегда опасение: а не идиот ли я? К сожалению, разобраться в чем же собственно суть дела, за всеми этими многоэтажными индексами и закорючками довольно сложно. Да и нужно ли? Здесь я не буду в очередной раз расписывать алгоритм и его математические основы, об этом можно почитать в книжках, википедии или где угодно еще. Сегодня я постараюсь изложить практическую сторону вопроса.

Итак, что же такое фильтр Калмана. Предположим у вас есть некая величина X изменяющаяся во времени. У вас есть некий прибор, измеряющий величину X с некоторой погрешностью. Вы хотите оценить, чему на самом деле равна величина X. Например, в случае акций Xk может быть истинной ценой акции в момент времени k (она определяется доходами компании, курсом валют, ожиданиями инвесторов итп). У нас есть только один способ измерить эту величину - пойти на биржу и посмотреть котировки. Но котировки на бирже подвержены влиянию массы факторов, цены колеблются в широком диапазоне и измеренная величина Yk может отличаться от Xx. Алгоритм Калмана позволяет нам оценить истинное значение Xk на основе истории изменения измеряемой величины Yk. При этом, естественно, делается ряд предположений о связи между Y и X и характере изменения X.


Предположение номер 1


Xk = Ak * Xk-1 + Bk * Uk + Wk-1

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

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

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

Предположение номер 2


Yk = Hk * Xk + Vk

Т.е. наша измеряемая величина Y это линейная комбинация реальной величины X и погрешности измерения V. При этом, делается дополнительное предположение что погрешность измерения тоже имеет нормальное распределение со своей дисперсией dV. Очень часто на практике H = 1.

Что мы только что сделали? Мы построили модель, описывающую нашу систему, при этом, заметьте, пока никакого фильтра Калмана не упоминалось.

Что мы хотим сделать?


Мы хотим получить некое уравнение

Xk = F(Xk-1, k-2..., Yk, k-1..., dW, dV)

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

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

Xk = Kk * Yk + (1 - Kk) * Xk-1

Где Kk - это Калмановский коэффициент усиления. Именно эту величину, в процессе работы и определяет алгоритм фильтра Калмана. При этом фильтр Калмана находит оптимальную величину коэффициента усиления при заданных дисперсиях dW и dV. Чем лучше вы оцените эти дисперсии, тем лучше будет работать фильтр Калмана.

А теперь, представьте, что X и Y это не одномерные величины, а вектор. A, B и H - матрицы. Погрешность и шум V и W тоже векторы. Оказывается, что при этом, все уравнения будут работать точно так же, и алгоритм фильтра Калмана легко обобщается на многомерные величины.

А теперь, по традиции, пример на R


Фильтр Калмана можно найти в нескольких модулях R, для примера я буду использовать реализацию из модуля dlm. В данном примере мы применим фильтр Калмана к значениям индекса S&P500.

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

Итак, пример:

library(dlm)
library(quantmod)

getSymbols("^GSPC", from="2010-09-01")
y <- as.ts(Cl(GSPC))
x <- dlmFilter(y, dlmModPoly(1, dV=0.1, dW=0.01))

plot(y, type='l')
lines(x$f[-1], col = "red")
#Добавим EMA для сравнения
lines(EMA(y), col = "green")

Для сравнения, на график, зеленым цветом, было добавлено 10-дневное экспоненциальное скользящее среднее Получилась картинка:


Большое спасибо читателю по имени Kirill за найденную в коде примера ошибку.

8 comments:

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

    ReplyDelete
  2. Не совсем так, фильтра Калмана не изменяет предыдующие значения. Он адаптируется к предыдущим для вычисления последующих. Но никакой правки уже посчитанных значений нет. К его использованию, я думаю, мы еще вернемся.

    ReplyDelete
  3. Картинка вызвала подозрения- после аналогичных участков кривая вдруг ведет себя по разному.
    Проверил пример в RSudio :обрезал данные Y(в окне WorkSpace) в интересных экстремальных местах,а потом добавлял по 1 точке.
    Увы) - пересчитывает на экстремуме 7-8 предыдущих точек.Прошлое предсказывает хорошо,а будущее не лучше EMA.
    P.S.Данные брал из финама по Газпрому.

    ReplyDelete
  4. Действительно косяк :) Спасибо, что заметили. У меня там ошибка в коде, я использую калмановское сглаживание вместо калмановской фильтрации, а оно действительно пересчитывает прошлые данные. Сейчас исправлю.

    ReplyDelete
  5. Самое интересное - это "отдельными задачами являются оценки параметров дисперсий dV и dW."

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

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

    ReplyDelete
  7. У меня есть данные. Помогите по ним прогноз через фильтр калмана сделать.

    ReplyDelete
  8. Есть результат в применении данного фильтра на фондовом рынке?
    При поверхностном анализе я пришел к выводу, что он не очень подходит для этого и требует априорных данных о системе.
    На фондовом рынке наверное нужно обратное - по реальным данным оценить состояние системы.

    ReplyDelete