Калмановская фильтрация. Современные проблемы науки и образования
Фильтр Калмана
Постановка задачи
Обозначим за x k величину, которую мы будем измерять, а потом фильтровать. Это может быть координата, скорость, ускорение, влажность, температура, давление, и т.д.
Начнем с простого примера, который и приведет нас к формулировке общей задачи. Представьте себе, что у нас есть радиоуправляемый движущейся объект, который может перемещаться только вперед и назад. Мы, зная вес объекта, форму, покрытие поверхности, по которой он перемещается и т.д., рассчитали как контролирующий джойстик влияет на скорость движения v k .
Ри. 19 Движущийся объект
Тогда координата объекта будет изменяться по закону:
x k+1 = x k + v k dt (3.7)
В реальной же жизни мы не можем учесть в наших расчетах малые возмущения, действующие на объект (ветер, неровности поверхности, препятствия), поэтому настоящая скорость объекта будет отличаться от расчетной. К правой части написанного уравнения добавится случайная величина о k :
x k+1 = x k + v k dt +о k (3.8)
У нас есть установленный на объекте GPS сенсор, который пытается мерить истинную координату x k объекта, и, конечно же, не может ее померить точно, а мерит с ошибкой з л , которая является тоже случайной величиной.
В итоге с сенсора мы получаем ошибочные данные:
z k = x k + з k (3.9)
Задача состоит в том, что, зная неверные показания сенсора z k , найти хорошее приближение для истинной координаты объекта x k . Это хорошее приближение мы будем обозначать как.
В формулировке же общей задачи, за координату x k может отвечать все что угодно (температура, влажность и т.д.), а член, отвечающий за контроль системы извне, мы обозначим за u k (в примере c объектом u k = v k dt ).
Уравнения для координаты и показания сенсора будут выглядеть так:
x k+1 = x k + u k +о k
z k = x k + з k (3.10)
И так, что нам известно:
· u k - это известная величина, которая контролирует эволюцию системы. Мы ее знаем из построенной нами физической модели.
· Ошибка модели о k и ошибка сенсора з k - случайные величины. И их законы распределения не зависят от времени (от номера итерации k ).
· Средние значения ошибок равны нулю: Eо k = Eз k = 0 .
· Сам закон распределения случайных величин может быть нам и не известен, но известны их дисперсии и. Заметим, что дисперсии не зависят от k , потому что законы распределения не зависят от него.
· Предполагается, что все случайные ошибки независимы друг от друга: какая ошибка будет в момент времени k совершенно не зависит от ошибки в другой момент времени k ".
Нелишним будет отметить, что задача фильтрации - это не задача сглаживания. Мы не стремимся сглаживать данные с сенсора, мы стремимся получить наиболее близкое значение к реальной координате x k .
Алгоритм Калмана
Мы будем рассуждать по индукции. Представьте себе, что на k -ом шаге мы уже нашли отфильтрованное значение с сенсора, которое хорошо приближает истинную координату системы x k . Не забываем, что мы знаем уравнение, контролирующее изменение нам неизвестной координаты:
x k+1 = x k + u k +о k ,
поэтому, еще не получая значение с сенсора, мы можем предположить, что на шаге k+1 система эволюционирует согласно этому закону и сенсор покажет что-то близкое к. К сожалению, пока мы не можем сказать ничего более точного. С другой стороны, на шаге k+1 у нас на руках будет неточное показание сенсора z k+1 .
Идея Калмана состоит в том, что чтобы получить наилучшее приближение к истинной координате x k+1 , мы должны выбрать золотую середину между показанием z k+1 неточного сенсора и - нашим предсказанием того, что мы ожидали от него увидеть. Показанию сенсора мы дадим вес K, а на предсказанное значение останется вес (1-K):
Коэффициент K называют коэффициентом Калмана. Он зависит от шага итерации, поэтому правильнее было бы писать K k+1 , но пока, чтобы не загромождать формулы расчетах, мы будем опускать его индекс.
Мы должны выбрать коэффициент Калмана K таким, чтобы получившееся оптимальное значение координаты было бы наиболее близко к истинной координате x k+1 . К примеру, если мы знаем, что наш сенсор очень точный, то мы будем больше доверять его показанию и дадим значению z k+1 больше весу (K близко единице). Если же сенсор, наоборот, совсем не точный, тогда больше будем ориентироваться на теоретически предсказанное значение.
В общем случае, чтобы найти точное значение коэффициента Калмана K , нужно просто минимизировать ошибку:
Используем уравнения (3.10), чтобы переписать выражение для ошибки:
Доказательство:
Теперь необходимо разобраться, что такое минимизировать ошибку? Ведь ошибка, как мы видим, сама по себе является случайной величиной и каждый раз принимает разные значения.
На самом деле не существует однозначного подхода к определению того, что означает, что ошибка минимальна.
Точно как и в случае с дисперсией случайной величины, когда мы пытались оценить характерную ширину ее разброса, так и тут мы выберем самый простой для расчетов критерий. Мы будем минимизировать среднее значение от квадрата ошибки:
Распишем последнее выражение:
Доказательство:
Из того что все случайные величины, входящие в выражение (5.13) для e k+1 , независимы, следует, что все "перекрестные" члены равны нулю: E (о k з k+1 ) = E (e k о k ) = E (e k з k+1 ) =0.
Так же Eз k+1 =Eо k =0 , тогда формула для дисперсий выглядит намного проще:
Выражение (3.15) принимает минимальное значение, когда (приравниваем производную к нулю):
Подставляем в выражение (3.15) для среднеквадратичной ошибки минимизирующее ее значение коэффициента Калмана (5.17) . Получаем:
Наша задача решена. Мы получили итерационную формулу (3.18), для вычисления коэффициента Калмана. Практическая реализация в приложении А.
Рис. 20
Если проследить, как с шагом итерации k изменяется коэффициент Калмана K k , то можно показать, что он всегда стабилизируется к определенному значению K stab . К примеру, когда среднеквадратичные ошибки сенсора и модели относятся друг к другу как десять к одному, то график коэффициента Калмана в зависимости от шага итерации выглядит так:
Рис. 21
Основная идея фильтра Калмана состоит в том, что надо найти коэффициент K такой, чтобы отфильтрованное значение:
в среднем меньше всего отличалось бы от реального значения координаты x k+1 . Мы видим, что отфильтрованное значение есть линейная функция от показания сенсора z k+1 и предыдущего отфильтрованного значения. А предыдущее отфильтрованное значение является, в свою очередь, линейной функцией от показания сенсора z k и предпредыдущего отфильтрованного значения. И так далее, пока цепь полностью не развернется. То есть отфильтрованное значение зависит от всех предыдущих показаний сенсора линейно:
Поэтому фильтр Калмана называют линейным фильтром.
Фильтр Калмана - это, наверное, самый популярный алгоритм фильтрации, используемый во многих областях науки и техники. Благодаря своей простоте и эффективности его можно встретить в GPS-приемниках, обработчиках показаний датчиков, при реализации систем управления и т.д.
Про фильтр Калмана в интернете есть очень много статей и книг (в основном на английском), но у этих статей довольно большой порог вхождения, остается много туманных мест, хотя на самом деле это очень ясный и прозрачный алгоритм. Я попробую рассказать о нем простым языком, с постепенным нарастанием сложности.
Для чего он нужен?
Любой измерительный прибор обладает некоторой погрешностью, на него может оказывать влияние большое количество внешних и внутренних воздействий, что приводит к тому, что информация с него оказывается зашумленной. Чем сильнее зашумлены данные тем сложнее обрабатывать такую информацию.Фильтр - это алгоритм обработки данных, который убирает шумы и лишнюю информацию. В фильтре Калмана есть возможность задать априорную информацию о характере системе, связи переменных и на основании этого строить более точную оценку, но даже в простейшем случае (без ввода априорной информации) он дает отличные результаты.
Рассмотрим простейший пример - предположим нам необходимо контролировать уровень топлива в баке. Для этого в бак устанавливается емкостный датчик, он очень прост в обслуживании, но обладает некоторыми недостатками - например, зависимость от заправляемого топлива (диэлектрическая проницаемость топлива зависит от многих факторов, например, от температуры), большое влияние «болтанки» в баке. В итоге, информация с него представляет типичную «пилу» с приличной амплитудой. Такого рода датчики часто устанавливаются на тяжелой карьерной технике (не смущайтесь объемам бака):
Фильтр Калмана
Немного отвлечемся и познакомимся с самим алгоритмом. Фильтр Калмана использует динамическую модель системы (например, физический закон движения), известные управляющие воздействия и множество последовательных измерений для формирования оптимальной оценки состояния. Алгоритм состоит из двух повторяющихся фаз: предсказание и корректировка. На первом рассчитывается предсказание состояния в следующий момент времени (с учетом неточности их измерения). На втором, новая информация с датчика корректирует предсказанное значение (также с учетом неточности и зашумленности этой информации):
Уравнения представлены в матричной форме, если вы не знаете линейную алгебру - ничего страшного, дальше будет упрощенная версия без матриц для случая с одной переменной. В случае с одной переменной матрицы вырождаются в скалярные значения.
Разберемся сначала в обозначениях: подстрочный индекс обозначает момент времени: k - текущий, (k-1) - предыдущий, знак «минус» в верхнем индексе обозначает, что это предсказанное промежуточное значение.
Описание переменных представлены на следующих изображениях:
Можно долго и нудно описывать, что означают все эти таинственные матрицы переходов, но лучше, на мой взгляд, на реальном примере попробовать применить алгоритм - чтобы абстрактные значения обрели реальный смысл.
Опробуем в деле
Вернемся к примеру с датчиком уровня топлива, так как состояние системы представлено одной переменной (объем топлива в баке), то матрицы вырождаются в обычные уравнения:
Определение модели процесса
Для того, чтобы применить фильтр, необходимо определить матрицы/значения переменных определяющих динамику системы и измерений F, B и H:F - переменная описывающая динамику системы, в случае с топливом - это может быть коэффициент определяющий расход топлива на холостых оборотах за время дискретизации (время между шагами алгоритма). Однако помимо расхода топлива, существуют ещё и заправки… поэтому для простоты примем эту переменную равную 1 (то есть мы указываем, что предсказываемое значение будет равно предыдущему состоянию).
B - переменная определяющая применение управляющего воздействия. Если бы у нас были дополнительная информация об оборотах двигателя или степени нажатия на педаль акселератора, то этот параметр бы определял как изменится расход топлива за время дискретизации. Так как управляющих воздействий в нашей модели нет (нет информации о них), то принимаем B = 0.
H - матрица определяющая отношение между измерениями и состоянием системы, пока без объяснений примем эту переменную также равную 1.
Определение сглаживающих свойств
R - ошибка измерения может быть определена испытанием измерительных приборов и определением погрешности их измерения.Q - определение шума процесса является более сложной задачей, так как требуется определить дисперсию процесса, что не всегда возможно. В любом случае, можно подобрать этот параметр для обеспечения требуемого уровня фильтрации.
Реализуем в коде
Чтобы развеять оставшиеся непонятности реализуем упрощенный алгоритм на C# (без матриц и управляющего воздействия):class KalmanFilterSimple1D
{
public double X0 {get; private set;} // predicted state
public double P0 { get; private set; } // predicted covariance
Public double F { get; private set; } // factor of real value to previous real value
public double Q { get; private set; } // measurement noise
public double H { get; private set; } // factor of measured value to real value
public double R { get; private set; } // environment noise
Public double State { get; private set; }
public double Covariance { get; private set; }
Public KalmanFilterSimple1D(double q, double r, double f = 1, double h = 1)
{
Q = q;
R = r;
F = f;
H = h;
}
Public void SetState(double state, double covariance)
{
State = state;
Covariance = covariance;
}
Public void Correct(double data)
{
//time update - prediction
X0 = F*State;
P0 = F*Covariance*F + Q;
//measurement update - correction
var K = H*P0/(H*P0*H + R);
State = X0 + K*(data - H*X0);
Covariance = (1 - K*H)*F;
}
}
// Применение...
Var fuelData = GetData();
var filtered = new List();
Var kalman = new KalmanFilterSimple1D(f: 1, h: 1, q: 2, r: 15); // задаем F, H, Q и R
kalman.SetState(fuelData, 0.1); // Задаем начальные значение State и Covariance
foreach(var d in fuelData)
{
kalman.Correct(d); // Применяем алгоритм
Filtered.Add(kalman.State); // Сохраняем текущее состояние
}
Результат фильтрации с данными параметрами представлен на рисунке (для настройки степени сглаживания - можно изменять параметры Q и R):
За рамками статьи осталось самое интересное - применение фильтра Калмана для нескольких переменных, задание взаимосвязи между ними и автоматический вывод значений для ненаблюдаемых переменных. Постараюсь продолжить тему как только появится время.
Надеюсь описание получилось не сильно утомительным и сложным, если остались вопросы и уточнения - добро пожаловать в комментарии)
- Tutorial
В интернете, в том числе и на хабре, можно найти много информации про фильтр Калмана. Но тяжело найти легкоперевариваемый вывод самих формул. Без вывода вся эта наука воспринимается как некое шаманство, формулы выглядят как безликий набор символов, а главное, многие простые утверждения, лежащие на поверхности теории, оказываются за пределами понимания. Целью этой статьи будет рассказать об этом фильтре на как можно более доступном языке.
Фильтр Калмана - это мощнейший инструмент фильтрации данных. Основной его принцип состоит в том, что при фильтрации используется информация о физике самого явления. Скажем, если вы фильтруете данные со спидометра машины, то инерционность машины дает вам право воспринимать слишком быстрые скачки скорости как ошибку измерения. Фильтр Калмана интересен тем, что в каком-то смысле, это самый лучший фильтр. Подробнее обсудим ниже, что конкретно означают слова «самый лучший». В конце статьи я покажу, что во многих случаях формулы можно до такой степени упростить, что от них почти ничего и не останется.
Ликбез
Перед знакомством с фильтром Калмана я предлагаю вспомнить некоторые простые определения и факты из теории вероятностей.Случайная величина
Когда говорят, что дана случайная величина , то имеют ввиду, что эта величина может принимать случайные значения. Разные значения она принимает с разной вероятностью. Когда вы кидаете, скажем, кость, то выпадет дискретное множество значений: . Когда речь идет, например, о скорости блуждающей частички, то, очевидно, приходится иметь дело с непрерывным множеством значений. «Выпавшие» значения случайной величины мы будем обозначать через но иногда, будем использовать ту же букву, которой обозначаем случайную величину:В случае с непрерывным множеством значений случайную величину характеризует плотность вероятности , которая нам диктует, что вероятность того, что случайная величина «выпадет» в маленькой окрестности точки длиной равна . Как мы видим из картинки, эта вероятность равна площади заштрихованного прямоугольника под графиком:
Довольно часто в жизни случайные величины распределены по Гауссу, когда плотность вероятности равна .
Мы видим, что функция имеет форму колокола с центром в точке и с характерной шириной порядка .
Раз мы заговорили о Гауссовом распределении, то грешно будет не упомянуть, откуда оно возникло. Также как и числа и прочно обосновались в математике и встречаются в самых неожиданных местах, так и распределение Гаусса пустило глубокие корни в теорию вероятностей. Одно замечательное утверждение, частично объясняющее Гауссово всеприсутствие, состоит в следующем:
Пусть есть случайная величина имеющая произвольное распределение (на самом деле существуют некие ограничения на эту произвольность, но они совершенно не жесткие). Проведем экспериментов и посчитаем сумму «выпавших» значений случайной величины. Сделаем много таких экспериментов. Понятно, что каждый раз мы будем получать разное значение суммы. Иными словами, эта сумма является сама по себе случайной величиной со своим каким-то определенным законом распределения. Оказывается, что при достаточно больших закон распределения этой суммы стремится к распределению Гаусса (к слову, характерная ширина «колокола» растет как ). Более подробно читаем в википедии: центральная предельная теорема . В жизни очень часто встречаются величины, которые складываются из большого количества одинаково распределенных независимых случайных величин, поэтому и распределены по Гауссу.
Среднее значение
Среднее значение случайной величины - это то, что мы получим в пределе, если проведем очень много экспериментов, и посчитаем среднее арифметическое выпавших значений. Среднее значение обозначают по-разному: математики любят обозначать через (математическое ожидание или mean value), а заграничные математики через (expectation). Физики же через или . Мы будем обозначать на заграничный лад: .Например, для Гауссова распределения , среднее значение равно .
Дисперсия
В случае с распределением Гаусса мы совершенно четко видим, что случайная величина предпочитает выпадать в некоторой окрестности своего среднего значения .Еще раз полюбоваться распределением Гаусса
Как видно из графика, характерный разброс значений порядка . Как же оценить этот разброс значений для произвольной случайной величины, если мы знаем ее распределение. Можно нарисовать график ее плотности вероятности и оценить характерную ширину на глаз. Но мы предпочитаем идти алгебраическим путем. Можно найти среднюю длину (модуль) отклонения от среднего значения: . Эта величина будет хорошей оценкой характерного разброса значений . Но мы с вами очень хорошо знаем, что использовать модули в формулах - одна головная боль, поэтому эту формулу редко используют для оценок характерного разброса.
Более простой способ (простой в смысле расчетов) - найти . Эту величину называют дисперсией, и часто обозначают как . Корень из дисперсии - хорошая оценка разброса случайной величины. Корень из дисперсии еще называют среднеквадратичным отклонением.
Например, для распределение Гаусса можно посчитать, что определенная выше дисперсия в точности равна , а значит среднеквадратичное отклонение равно , что очень хорошо согласуется с нашей геометрической интуицией.
На самом деле тут скрыто маленькое мошенничество. Дело в том, что в определении распределения Гаусса под экспонентой стоит выражение . Эта двойка в знаменателе стоит именно для того, чтобы среднеквадратичное отклонение равнялось бы коэффициенту . То есть сама формула распределения Гаусса написана в виде, специально заточенном для того, что мы будем считать ее среднеквадратичное отклонение.
Независимые случайные величины
Случайные величины бывают зависимыми и нет. Представьте, что вы бросаете иголку на плоскость и записываете координаты ее обоих концов. Эти две координаты зависимы, они связаны условием, что расстояние между ними всегда равно длине иголки, хотя и являются случайными величинами.Случайные величины независимы, если результат выпадения первой из них совершенно не зависит от результата выпадения второй из них. Если случайные величины и независимы, то среднее значение их произведения равно произведению их средних значений:
Доказательство
Например, иметь голубые глаза и окончить школу с золотой медалью - независимые случайные величины. Если голубоглазых, скажем а золотых медалистов , то голубоглазых медалистов Этот пример подсказывает нам, что если случайные величины и заданы своими плотностями вероятности и , то независимость этих величин выражается в том, что плотность вероятности (первая величина выпала , а вторая ) находится по формуле:
Из этого сразу же следует, что:
Как вы видите, доказательство проведено для случайных величин, которые имеют непрерывный спектр значений и заданы своей плотностью вероятности. В других случаях идея доказательтсва аналогичная.
Фильтр Калмана
Постановка задачи
Обозначим за величину, которую мы будем измерять, а потом фильтровать. Это может быть координата, скорость, ускорение, влажность, степень вони, температура, давление, и т.д.Начнем с простого примера, который и приведет нас к формулировке общей задачи. Представьте себе, что у нас есть радиоуправляемая машинка, которая может ехать только вперед и назад. Мы, зная вес машины, форму, покрытие дороги и т.д., расcчитали как контролирующий джойстик влияет на скорость движения .
Тогда координата машины будет изменяться по закону:
В реальной же жизни мы не можем учесть в наших расчетах маленькие возмущения, действующие на машину (ветер, ухабы, камушки на дороге), поэтому настоящая скорость машины будет отличаться от расчетной. К правой части написанного уравнения добавится случайная величина :
У нас есть установленный на машинке GPS сенсор, который пытается мерить истинную координату машинки, и, конечно же, не может ее померить точно, а мерит с ошибкой , которая является тоже случайной величиной. В итоге с сенсора мы получаем ошибочные данные:
Задача состоит в том, что, зная неверные показания сенсора , найти хорошее приближение для истинной координаты машины . Это хорошее приближение мы будем обозначать как .
В формулировке же общей задачи, за координату может отвечать все что угодно (температура, влажность...), а член, отвечающий за контроль системы извне мы обозначим за (в примере c машиной ). Уравнения для координаты и показания сенсора будут выглядеть так:
(1)
Давайте подробно обсудим, что нам известно:
Нелишним будет отметить, что задача фильтрации - это не задача сглаживания. Мы не стремимся сглаживать данные с сенсора, мы стремимся получить наиболее близкое значение к реальной координате .
Алгоритм Калмана
Мы будем рассуждать по индукции. Представьте себе, что на -ом шаге мы уже нашли отфильтрованное значение с сенсора , которое хорошо приближает истинную координату системы . Не забываем, что мы знаем уравнение, контролирующее изменение нам неизвестной координаты:Поэтому, еще не получая значение с сенсора, мы можем предположить, что на шаге система эволюционирует согласно этому закону и сенсор покажет что-то близкое к . К сожалению, пока мы не можем сказать ничего более точного. С другой стороны, на шаге у нас на руках будет неточное показание сенсора .
Идея Калмана состоит в том, что чтобы получить наилучшее приближение к истинной координате , мы должны выбрать золотую середину между показанием неточного сенсора и - нашим предсказанием того, что мы ожидали от него увидеть. Показанию сенсора мы дадим вес а на предсказанное значение останется вес :
Коэффициент называют коэффициентом Калмана. Он зависит от шага итерации, поэтому правильнее было бы писать , но пока, чтобы не загромождать формулы расчетах, мы будем опускать его индекс.
Мы должны выбрать коэффициент Калмана таким, чтобы получившееся оптимальное значение координаты было бы наиболее близко к истинной координате . К примеру, если мы знаем, что наш сенсор очень точный, то мы будем больше доверять его показанию и дадим значению больше весу ( близко единице). Eсли же сенсор, наоборот, совсем не точный, тогда больше будем ориентироваться на теоретически предсказанное значение .
В общем случае, чтобы найти точное значение коэффициента Калмана , нужно просто минимизировать ошибку:
Используем уравнения (1) (те которые на голубом фоне в рамочке), чтобы переписать выражение для ошибки:
Доказательство
Теперь самое время обсудить, что означает выражение минимизировать ошибку? Ведь ошибка, как мы видим, сама по себе является случайной величиной и каждый раз принимает разные значения. На самом деле не существует однозначного подхода к определению того, что означает, что ошибка минимальна. Точно как и в случае с дисперсией случайной величины, когда мы пытались оценить характерную ширину ее разброса, так и тут мы выберем самый простой для расчетов критерий. Мы будем минимизировать среднее значение от квадрата ошибки:
Распишем последнее выражение:
ключ к доказательству
Из того что все случайные величины, входящие в выражение для , независимы и средние значения ошибок сенсора и модели равны нулю: , следует, что все «перекрестные» члены равны нулю:
.
Плюс к этому, формулы для дисперсий выглядит намного проще: и (так как )
Это выражение принимает минимальное значение, когда (приравниваем производную к нулю)
Здесь мы уже пишем выражение для коэффициента Калмана с индексом шага , тем самым мы подчеркиваем, что он зависит от шага итерации.
Подставляем в выражение для среднеквадратичной ошибки минимизирующее ее значение коэффициента Калмана . Получаем:
Наша задача решена. Мы получили итерационную формулу, для вычисления коэффициента Калмана.
Все формулы в одном месте
Пример
На рекламной картинке в начале статьи отфильтрованы данные с вымышленного GPS сенсора, установленного на вымышленной машине, которая едет равноускоренно c известным вымышленным ускорением .Еще раз посмотреть на результат фильтрования
Код на матлабе
clear all;
N=100 % number of samples
a=0.1 % acceleration
sigmaPsi=1
sigmaEta=50;
k=1:N
x=k
x(1)=0
z(1)=x(1)+normrnd(0,sigmaEta);
for t=1:(N-1)
x(t+1)=x(t)+a*t+normrnd(0,sigmaPsi);
z(t+1)=x(t+1)+normrnd(0,sigmaEta);
end;
%kalman filter
xOpt(1)=z(1);
eOpt(1)=sigmaEta; % eOpt(t) is a square root of the error dispersion (variance). It"s not a random variable.
for t=1:(N-1)
eOpt(t+1)=sqrt((sigmaEta^2)*(eOpt(t)^2+sigmaPsi^2)/(sigmaEta^2+eOpt(t)^2+sigmaPsi^2))
K(t+1)=(eOpt(t+1))^2/sigmaEta^2
xOpt(t+1)=(xOpt(t)+a*t)*(1-K(t+1))+K(t+1)*z(t+1)
end;
plot(k,xOpt,k,z,k,x)
Анализ
Если проследить, как с шагом итерации изменяется коэффициент Калмана , то можно показать, что он всегда стабилизируется к определенному значению . К примеру, когда среднеквадратичные ошибки сенсора и модели относятся друг к другу как десять к одному, то график коэффициента Калмана в зависимости от шага итерации выглядит так:
В следующем примере мы обсудим как это поможет существенно облегчить нашу жизнь.
Второй пример
На практике очень часто бывает, что нам вообще ничего не известно о физической модели того, что мы фильтруем. К примеру, вы захотели отфильтровать показания с вашего любимого акселерометра. Вам же заранее неизвестно по какому закону вы намереваетесь крутить акселерометр. Максимум информации, которую вы можете выцепить - это дисперсия ошибки сенсора . В такой непростой ситуации все незнание модели движения можно загнать в случайную величину :Но, откровенно говоря, такая система уже совершенно не удовлетворяет тем условиям, которые мы налагали на случайную величину , ведь теперь туда запрятана вся неизвестная нам физика движения, и поэтому мы не можем говорить, что в разные моменты времени ошибки модели независимы друг от друга и что их средние значения равны нулю. В этом случае, по большому счету, теория фильтра Калмана не применима. Но, мы не будем обращать внимания на этот факт, а, тупо применим все махину формул, подобрав коэффициенты и на глаз, так чтобы отфильтрованные данные миленько смотрелись.
Но можно пойти по другому, намного более простому пути. Как мы видели выше, коэффициент Калмана с увеличением номера шага всегда стабилизируется к значению . Поэтому вместо того, чтобы подбирать коэффициенты и и находить по сложным формулам коэффициент Калмана , мы можем считать этот коэффициент всегда константой, и подбирать только эту константу. Это допущение почти ничего не испортит. Во-первых, мы уже и так незаконно пользуемся теорией Калмана, а во-вторых коэффициент Калмана быстро стабилизируется к константе. В итоге все очень упростится. Нам вообще никакие формулы из теории Калмана не нужны, нам просто нужно подобрать приемлемое значение и вставить в итерационную формулу:
На следующем графике показаны отфильтрованные двумя разными способами данные с вымышленного сенсора. При условии того, что мы ничего не знаем о физике явления. Первый способ - честный, со всеми формулами из теории Калмана. А второй - упрощенный, без формул.
Как мы видим, методы почти ничем не отличаются. Маленькое отличие наблюдается только вначале, когда коэффициент Калмана еще не стабилизировался.
Обсуждение
Как мы увидели, основная идея фильтра Калмана состоит в том, что надо найти коэффициент такой, чтобы отфильтрованное значениеВ среднем меньше всего отличалось бы от реального значения координаты . Мы видим, что отфильтрованное значение есть линейная функция от показания сенсора и предыдущего отфильтрованного значения . А предыдущее отфильтрованное значение является, в свою очередь, линейной функцией от показания сенсора и предпредыдущего отфильтрованного значения . И так далее, пока цепь полностью не развернется. То есть отфильтрованное значение зависит от всех предыдущих показаний сенсора линейно:
Поэтому фильтр Калмана называют линейным фильтром.
Можно доказать, что из всех линейных фильтров Калмановский фильтр самый лучший. Самый лучший в том смысле, что средний квадрат ошибки фильтра минимален.
Многомерный случай
Всю теорию фильтра Калмана можно обобщить на многомерный случай. Формулы там выглядят чуть страшнее, но сама идея их вывода такая же, как и в одномерном случае. В этой прекрасной статье вы можете увидеть их: 1Проведено исследование использования фильтра Калмана в современных разработках комплексированных навигационных систем. Приведен и разобран пример построения математической модели, использующей расширенный фильтр Калмана для повышения точности определения координат беспилотных летательных аппаратов. Рассмотрен частичный фильтр. Сделан краткий обзор научных работ, использующих данный фильтр для повышения надежности и отказоустойчивости навигационных систем. Данная статья позволяет сделать вывод, что использование фильтра Калмана в системах определения местоположения БПЛА практикуется во многих современных разработках. Существует огромное количество вариаций и аспектов такого использования, которое дает и ощутимые результаты в повышении точности, особенно в случае отказа стандартных спутниковых навигационных систем. Это является главным фактором влияния данной технологии на различные научные области, связанные с разработкой точных и отказоустойчивых навигационных систем для различных летательных аппаратов.
фильтр Калмана
навигация
беспилотный летательный аппарат (БПЛА)
1. Макаренко Г.К., Алешечкин А.М. Исследование алгоритма фильтрации при определении координат объекта по сигналам спутниковых радионавигационных систем // Доклады ТУСУРа. – 2012. – № 2 (26). – С. 15-18.
2. Bar-Shalom Y., Li X. R., Kirubarajan T. Estimation with Applications
to Tracking and Navigation // Theory Algorithms and Software. – 2001. – Vol. 3. – P. 10-20.
3. Bassem I.S. Vision based Navigation (VBN) of Unmanned Aerial Vehicles (UAV) // UNIVERSITY OF CALGARY. – 2012. – Vol. 1. – P. 100-127.
4. Conte G., Doherty P. An Integrated UAV Navigation System Based on Aerial Image Matching // Aerospace Conference. – 2008. –Vol. 1. – P. 3142-3151.
5. Guoqiang M., Drake S., Anderson B. Design of an extended kalman filter for uav localization // In Information, Decision and Control. – 2007. – Vol. 7. – P. 224–229.
6. Ponda S.S Trajectory Optimization for Target Localization Using Small Unmanned Aerial Vehicles // Massachusetts institute of technology. – 2008. – Vol. 1. – P. 64-70.
7. Wang J., Garrat M., Lambert A. Integration of gps/ins/vision sensors to navigate unmanned aerial vehicles // IAPRS&SIS. – 2008. – Vol. 37. – P. 963-969.
Одной из актуальных задач современной навигации беспилотных летательных аппаратов (БПЛА) является задача повышения точности определения координат. Эта задача решается путём использования различных вариантов комплексирования навигационных систем. Одним из современных вариантов комплексирования является сочетание gps/глонасс-навигации с расширенным фильтром Калмана (Extended Kalmanfilter), рекурсивно оценивающего точность с помощью неполных и зашумленных измерений. В данный момент существуют и разрабатываются различные вариации расширенного фильтра Калмана, включающие разнообразное число переменных состояний . В этой работе мы покажем, насколько эффективным может быть его использование в современных разработках. Рассмотрим одно из характерных представлений подобного фильтра .
Построение математической модели
В данном примере мы будем говорить только о движении БПЛА в горизонтальной плоскости, иначе, мы рассмотрим так называемую проблему 2d локализации . В нашем случае это оправдано тем, что для многих практически встречающихся ситуаций БПЛА может оставаться примерно на одной и той же высоте. Это предположение широко используется для упрощения моделирования динамики летательных аппаратов . Динамическая модель БПЛА задается следующей системой уравнений:
где {} - координаты БПЛА в горизонтальной плоскости как функции времени, направление БПЛА, угловая скорость БПЛА, и vпутевая скорость БПЛА, функции и будем считать постоянными. Они взаимно независимы, с известными ковариациями
и , равными и соответственно и используются для моделирования изменений ускорения БПЛА, вызванных ветром, маневрами пилота и т.д. Значения и являются производными от максимальной угловой скорости БПЛА и опытных значений изменений линейной скорости БЛА, - символ Кронекера.
Данная система уравнений будет приближенной из-за нелинейности в модели и из-за присутствия шума. Самый простой способ аппроксимации в данном случае - это приближение методом Эйлера. Дискретная модель динамической системы движения БПЛА показана ниже.
дискретный вектор состояний фильтра Калмана, позволяющий аппроксимировать значение непрерывного вектора состояний. ∆ - временной интервал между k и k+1 измерениями. {} и {} - последовательности значений белого гауссовского шума с нулевым средним значением. Матрица ковариации для первой последовательности:
Аналогично, для второй последовательности:
Выполнив соответствующие замены в уравнениях системы (2), получаем:
Последовательности и взаимно независимы. Они также являются последовательностями белого гауссовского шума с нулевым средним значением и с матрицами ковариации и соответственно. Преимущество этой формы в том, что она показывает изменение дискретного шума в интервале между каждыми измерениями. В итоге получаем следующую дискретную динамическую модель:
(3)
Уравнение для :
= 
+ , (4)
где, х и y - координаты БПЛА в k-момент времени, а гауссовская последовательность случайных параметров с нулевым средним значением, которая используется для задания погрешности. Предполагается, что эта последовательность не зависит от {} и {}.
Выражения (3) и (4) служат основой для оценки местоположения БПЛА, где к-е координаты получены с помощью расширенного фильтра Калмана. Моделлирование отказа навигационных систем применительно к данному типу фильтра показывает его существенную эффективность .
Для большей наглядности приведем небольшой простой пример. Пусть некоторый БПЛА летит равноускоренно, с некоторым постоянным ускорением а.
Где, х - координата БПЛА в t-момент времени, а δ - некоторая случайная величина.
Предположим, что у нас есть gps-сенсор, который получает данные о местоположении летательного аппарата. Представим результат моделирования данного процесса в программном пакете MATLAB.
Рис. 1. Фильтрация показания сенсора с помощью фильтра Калмана
На рис. 1 видно, насколько эффективным может быть использование фильтрации по алгоритму Калмана.
Однако в реальной ситуации сигналы зачастую имеют нелинейную динамику и ненормальный шум. Именно в таких случаях и используется расширенный фильтр Калмана. В том случае, если дисперсии шумов не слишком велики (т.е. линейная аппроксимация является адекватной), применение расширенного фильтра Калмана дает решение задачи с высокой точностью. Однако в том случае, когда шумы не являются гауссовскими, расширенный фильтр Калмана применять нельзя. В этом случае обычно применяют частичный фильтр, в котором используются численные методы взятия интегралов на основе методов Монте Карло с марковскими цепями.
Частичный фильтр
Представим один из алгоритмов, развивающих идеи расширенного фильтра Калмана - частичный фильтр. Частичная фильтрация является неоптимальным способом фильтрации, который работает при выполнении объединения методом Монте-Карло на множестве частиц, которые представляют собой распределение вероятностей процесса. Здесь частица это элемент, взятый из априорного распределения оцениваемого параметра. Основная идея частичного фильтра заключается в том, что большое количество частиц может быть использовано для представления оценки распределения. Чем большее число частиц используется, тем точнее множество частиц будет представлять априорное распределение. Фильтр частиц инициализируется помещением в него N частиц из априорного распределения параметров, которые мы хотим оценить. Алгоритм фильтрации предполагает прогон этих частиц через специальную систему, а затем взвешивание с помощью информации, полученной от измерения данных частиц. Полученные частицы и связанные с ними массы представляют апостериорное распределение оценочного процесса. Цикл повторяется для каждого нового измерения, и веса частиц обновляются для представления последующего распределения. Одна из основных проблем с традиционным подходом фильтрации частиц состоит в том, что в результате такой подход обычно имеет несколько частиц, имеющих очень большой вес, в отличие от большинства остальных, вес которых очень незначителен. Это приводит к нестабильности фильтрации . Эта проблема может быть решена введением частоты дискретизации, где N новых частиц берется из распределения, составленного из старых частиц. Результат оценки получают путем получения выборки среднего значения множества частиц. Если мы имеем несколько независимых выборок, то средняя выборка будет точной оценкой среднего значения, задающей конечную дисперсию.
Даже если фильтр частиц является неоптимальным, то при стремлении количества частиц к бесконечности эффективность алгоритма приближается в байесову правилу оценивания. Поэтому желательно иметь столько частиц, сколько возможно, чтобы получить наилучший результат. К сожалению, это приводит к сильному увеличению сложности вычислений, а, следовательно, вынуждает к поиску компромисса между точностью и скоростью расчета. Итак, число частиц должно быть выбрано исходя из требований к задаче оценки точности. Еще одним важным фактором для работы фильтра частиц является ограничение на частоту дискретизации. Как упоминалось ранее, частота дискретизации является важным параметром фильтрации частиц и без него в конечном итоге алгоритм становится вырождающимся. Идея заключается в том, что если веса распределяются слишком неравномерно и порог дискретизации скоро будет достигнут, то частицы с низким весом отбрасываются, и оставшееся множество образует новую вероятностную плотность, для которой могут браться новые выборки. Выбор порога частоты дискретизации представляет собой довольно сложную задачу, ведь слишком высокая частота служит причиной чрезмерной чувствительности фильтра к шуму, а слишком низкая дает большую погрешность. Также важным фактором является плотность вероятности .
В целом, алгоритм фильтрации частиц показывает хорошую производительность расчета местоположения для стационарных целей и в случае относительно медленно движущихся целей с неизвестной динамикой ускорения. В общем случае, алгоритм фильтрации частиц является более стабильным, чем расширенный фильтр Калмана, и менее склонным к вырождению и серьезным сбоям. В случаях нелинейного, негауссового распределения данный алгоритм фильтрации показывает весьма хорошую точность определения местоположения цели, в то время как алгоритм расширенной фильтрации Калмана нельзя использовать при таких условиях. К минусам данного подхода можно отнести его более высокую сложность относительно расширенного фильтра Калмана, а также то, что не всегда очевидно, как правильно подобрать параметры для этого алгоритма.
Перспективные исследования в данной области
Использование модели фильтра Калмана, подобной той, что привели мы, можно видеть в , где он используется для улучшения характеристик комплексированной системы (GPS + модель компьютерного зрения для сопоставления с географической базой), и также моделируется ситуация отказа спутникового навигационного оборудования. С помощью фильтра Калмана результаты работы системы в случае отказа были существенно улучшены (например, погрешность в определении высоты была снижена примерно в два раза, а погрешности в определении координат по разным осям снижены практически в 9 раз). Аналогичное использование фильтра Калмана приведено также в .
Интересная с точки зрения совокупности методов задача решается в . Там также используется фильтр Калмана с 5 состояниями, с некоторыми отличиями в построении модели. Полученный результат превосходит результат приведенной нами модели за счет использования дополнительных средств комплексирования (используются фото и тепловизионные изображения). Применение фильтра Калмана в данном случае позволяет уменьшить погрешность определения пространственных координат заданной точки до значения 5,5 м.
Заключение
В заключение отметим, что использование фильтра Калмана в системах определения местоположения БПЛА практикуется во многих современных разработках. Существует огромное количество вариаций и аспектов такого использования, вплоть до одновременного применения нескольких подобных фильтров с разными факторами состояний . Одним из наиболее перспективных направлений развития Калмановских фильтров видится работа над созданием модифицированного фильтра, погрешности которого будут представлены цветным шумом, что сделает его еще более ценным для решения реальных задач. Также большой интерес в данной области представляет собой частичный фильтр, с помощью которого можно фильтровать негауссовские шумы. Названное разнообразие и ощутимые результаты в повышении точности, особенно в случае отказа стандартных спутниковых навигационных систем, являются главными факторами влияния данной технологии на различные научные области, связанные с разработкой точных и отказоустойчивых навигационных систем для различных летательных аппаратов.
Рецензенты:
Лабунец В.Г., д.т.н., профессор, профессор кафедры теоретических основ радиотехники Уральского федерального университета имени первого Президента России Б.Н. Ельцина, г. Екатеринбург;
Иванов В.Э., д.т.н., профессор, зав. кафедрой технологии и средств связи Уральского федерального университета имени первого Президента России Б.Н. Ельцина, г. Екатеринбург.
Библиографическая ссылка
Гаврилов А.В. ИСПОЛЬЗОВАНИЕ ФИЛЬТРА КАЛМАНА ДЛЯ РЕШЕНИЯ ЗАДАЧ УТОЧНЕНИЯ КООРДИНАТ БПЛА // Современные проблемы науки и образования. – 2015. – № 1-1.;URL: http://science-education.ru/ru/article/view?id=19453 (дата обращения: 06.04.2019). Предлагаем вашему вниманию журналы, издающиеся в издательстве «Академия Естествознания»
Random Forest - один из моих любимых алгоритмов data mining. Во-первых он невероятно универсален, с его помощью можно решать как задачи регрессии так и классификации. Проводить поиск аномалий и отбор предикторов. Во-вторых это тот алгоритм, который действительно сложно применить неправильно. Просто потому, что в отличии от других алгоритмов у него мало настраиваемых параметров. И еще он удивительно прост по своей сути. И в то же время он отличается удивительной точностью.
В чем же идея такого замечательного алгоритма? Идея проста: допустим у нас есть какой-то очень слабый алгоритм, скажем, . Если мы сделаем очень много разных моделей с использованием этого слабого алгоритма и усредним результат их предсказаний, то итоговый результат будет существенно лучше. Это, так называемое, обучение ансамбля в действии. Алгоритм Random Forest потому и называется "Случайный Лес", для полученных данных он создает множество деревьев приятия решений и потом усредняет результат их предсказаний. Важным моментом тут является элемент случайности в создании каждого дерева. Ведь понятно, что если мы создадим много одинаковых деревьев, то результат их усреднения будет обладать точностью одного дерева.
Как он работает? Предположим, у нас есть некие данные на входе. Каждая колонка соответствует некоторому параметру, каждая строка соответствует некоторому элементу данных.
Мы можем выбрать, случайным образом, из всего набора данных некоторое количество столбцов и строк и построить по ним дерево принятия решений.
Thursday, May 10, 2012
Thursday, January 12, 2012
Вот собственно и всё. 17-ти часовой перелет позади, Россия осталась за океаном. А в окно уютной 2-ух спальной квартиры на нас смотрит Сан-Франциско, знаменитая Кремниевая долина, Калифорния, США. Да, это и есть та самая причина, по которой я практически не писал последнее время. Мы переехали.
Всё это началось еще в апреле 2011 года, когда я проходил телефонное интервью в компании Zynga. Тогда это все казалось какой-то игрой не имеющей отношения к реальности и я и представить себе не мог, во что это выльется. В июне 2011 года Zynga приехали в Москву и провели серию собеседований, рассматривалось около 60 кандидатов прошедших телефонное интервью и из них было отобрано около 15 человек (точное число не знаю, кто-то потом передумал, кто-то сразу отказался). Интервью оказалось неожиданно простым. Ни тебе задачек на программирование, ни заковыристых вопросов про форму люков, в основном проверялись способности болтать. А знания, на мой взгляд, оценивались лишь поверхностно.
А дальше началась канитель. Сначала мы ждали результатов, потом офера, потом одобрение LCA, потом одобрения петиции на визу, потом документы из США, потом очередь в посольстве, потом дополнительную проверку, потом визу. Временами мне казалось, что я готов все бросить и забить. Временами я сомневался, а нужна ли нам эта Америка ведь и в России не плохо. Весь процесс занял где-то около полугода, в итоге, в середине декабря мы получили визы и начали готовиться к отъезду.
В понедельник был мой первый рабочий день на новом месте. В офисе созданы все условия для того чтобы не только работать, но и жить. Завтраки, обеды и ужины от собственных поваров, куча разнообразнейшей еды распиханной по всем уголкам, спортзал, массаж и даже парикмахер. Все это совершенно бесплатно для сотрудников. Многие добираются на работу на велосипеде и для хранения транспорта оборудовано несколько комнат. В общем, ничего подобного в России мне встречать не доводилось. Всему, однако, есть своя цена, нас сразу предупредили, что работать придется много. Что такое "много", по их меркам, мне не очень понятно.
Надеюсь, однако, что несмотря на количество работы, в обозримом будущем смогу возобновить ведение блога и, может быть, расскажу что-нибудь о американской жизни и работе программистом в Америке. Поживем - увидим. А пока, поздравляю всех с наступившим новым годом и рождеством и до новых встреч!
Для примера использования, распечатаем дивидендную доходность российских компаний. В качестве базовой цены, берем цену закрытия акции в день закрытия реестра. Почему-то на сайте тройки этой информации нет, а она ведь гораздо интересней чем абсолютные величины дивидендов.
Внимание! Код выполняется довольно долго, т.к. для каждой акции требуется сделать запрос на сервера finam и получить её стоимость.
Result <- NULL for(i in (1:length(divs[,1]))){ d <- divs if (d$Divs>0){ try({ quotes <- getSymbols(d$Symbol, src="Finam", from="2010-01-01", auto.assign=FALSE) if (!is.nan(quotes)){ price <- Cl(quotes) if (length(price)>0){ dd <- d$Divs result <- rbind(result, data.frame(d$Symbol, d$Name, d$RegistryDate, as.numeric(dd)/as.numeric(price), stringsAsFactors=FALSE)) } } }, silent=TRUE) } } colnames(result) <- c("Symbol", "Name", "RegistryDate", "Divs") result
Аналогично можно построить статистику для прошлых лет.
