Навигация по расстояниям до точечных ориентиров адаптивным методом существенной выборки Текст научной статьи по специальности «Физика»

-►

Математическое моделирование: методы, алгоритмы, технологии

УДК 519.6

Д.Г. Арсеньев, В.М. Иванов, Н.А. Берковский

НАВИГАЦИЯ ПО РАССТОЯНИЯМ ДО ТОЧЕЧНЫХ ОРИЕНТИРОВ АДАПТИВНЫМ МЕТОДОМ СУЩЕСТВЕННОй ВЫБОРКИ

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

В статье описано применение адаптивного метода существенной выборки [1-3] к вычислению интегралов, входящих в формулу для оптимальной байесовской оценки, и проведено сравнение результатов с вычислениями по схеме, применяемой для расчета апостериорных средних в методе important sampling [6, 9], являющемся одним

из основных последовательных методов Монте-Карло (Sequential Monte Carlo Methods). Численные примеры свидетельствуют о преимуществах предлагаемого адаптивного метода.

Постановка задачи

Наша цель - определить координаты объекта (xp x2) по известным результатам измерений расстояний до двух точечных неподвижных ориентиров, координаты которых (xj, xj) и (x2, x2) предполагаются известными. Предположим, что сделано m пар измерений из одного положения (по одному измерению для каждого ориентира в каждой паре). Вводя измерительные шумы, обозначим результаты этих измерений как

У-=}/(х1-х^ + (х2-х^ +У,2; (г = 1т).

Ошибки измерений предполагаем одинаково распределенными, гауссовскими, независимыми и центрированными, т. е. V* ~ N(0,г), (г = 1,т; к = 1,2). Здесь 1Г(0,г) обозначает нормальное распределение на прямой с центром в нуле и дисперсией г2 . В рамках байесовского подхода мы должны задать некоторое априорное распределение положения объекта на плоскости, в данной статье оно задается функцией р(хр х2) = N (х1, х2, х°°, х°2, о1, о2, 0), где N (х1, х2, х0, х2, о1, о2, 0) обозначает нормальное распределение на плоскости с независимыми компонентами, распределениями которых являются N (х^, о1) и N (х2, о2). Апостериорная плотность распределения вектора (х1, х2) при условии известных результатов серии измерений у1 и у2 (г = 1, т ) задается формулой

р(х1;х2 I y};yf)=

— Х1 Г + (х2 ~ х2 Г '

к=\ ^

~х\1+(х2 - >Д

значение функции плотности распределения нормального закона N \ д/^ - х\ ^ + (х:2 - х\) ,г | в

где N

точке уьк , (к = 1, т s = 1,2), а С - нормирующий множитель. Отсюда имеем формулы для опти-« К Л2 , ( 1 с мальных байесовских оценок х1 и х2координат

У к. V 1*1 - ] +{Х2-Х2),Г\ обозначает объекта. 1 2 Р

_к=\ V_¿_\_)

Ш^'-^П1* у\, ^(х, -х})2 +(х2-к1)2, ^ , ^(х, , г^ёъ

к=\ \ ) V )

(1)

(2)

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

Схема применения АСВ к вычислению оптимальной оценки положения объекта

Рассмотрим применение АСВ к вычислению оценки первой координаты объекта, определенной формулой (1). Для краткости обозначим х = (х1, х2), х е Я2, и рассмотрим функции этого векторного аргумента:

к=1

Здесь и далее символ йх означает элементарный элемент площади на плоскости. Предположим, что статистическая оценка величины вычисляется с помощью метода Монте-Карло в виде

N

а I

х = , где /1=—> —.—^

1 ~2

хн[ ~Х1? +{х2 ~х12)2,г |х

V т

хт у1 ,-\1(х1 -х\ )2 + (х2 - х\^, ^ (1x^x2;

т

(*) = ' х2 )ПХ к=1

--«{Г +[Х2 ~х1)2 >г

(3)

(**)'

последовательность независимых случайных векторов в Я2, распределенных согласно некоторой функции плотности распределения вероятности р(х). В соответствии с обобщенной центральной предельной теоремой (дельта-методом) [7, 8 ]~можно считать, что распределение частного X при достаточно больших N близко к нормальному распределению

вида N

X],

2 2 а? ст7

Ч 2

[1

2СОУ(71,72)

\

и От

где

(4)

хы

У к > — )2 + [х2 ~х2^' Н

В этих обозначениях формула (1) перепишется следующим образом:

1 среднеквадратические отклонения

случайных величин 11 и 12. Вводя обозначения

с = 14 ' и т) = . . , где вектор х распределен с

р(х) р(х)

ф.п.р.в. р(х), можно доказать, что

соу(/~,72)= МЫ~Мт = СОУ(^Л) , (6)

откуда имеем приближенное распределение для

С г^-—-;-г\

X вида N

Л Л А т —0,5

х1, Х1N <

I* 1г ЬЬ

_ А _ й2

Таким образом, для среднеквадратической ошибки с!^} случайной величин^1 а1 можем написать приближенное равенство:

X

ïjx,}- JCj

N

-0,5

V

+ 2соу(Ы _ (7)

v

hh

Чем меньше эта дисперсия, тем точнее величина Xj аппроксимирует искомую байесовскую оценку первой координаты объекта xj. Легко заметить, что при фиксированном N дисперсия ojjcj} зависит от плотности p(x) не самым простым образом. Часто в задачах теории фильтрации, в частности, в методе important sampling в качестве p(x) берут априорную плотность распределения координат объекта, в данном случае это p(x) = p(xj, x2) = N (xj, x2, xj, xj, o1, o2, 0). Такой выбор удобен для проведения дальнейшей процедуры оценивания координат объекта (как правило, движущегося) в режиме «on-line» в случае марковской модели [6]. В случае же измерений, проводимых из одного положения, такой выбор плотности p(x) не является оптимальным.

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

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

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

В-третьих, поскольку подынтегральные функции в (1) и (2) весьма близки к нулю на большинстве точек области определения и имеют резкие максимумы, то алгоритм бисекции применяется не к самим интеграндам, а к их сдвигам вверх на некоторую положительную константу, так же как это делалось в при решении одномерной тестовой задачи навигации в [3]. Без этой «меры предосторожности», к сожалению, АСВ в случае близких к нулю функций с ненулевой вероятностью может генерировать непригодные для вычислений плотности. Видимо, сдвиг интеграндов на константу не самый эффективный, хотя и простой способ избежать этих нежелательных явлений, однако альтернативный алгоритм пока находится в стадии разработки.

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

Результаты численного эксперимента

Рассмотрим следующую навигационную задачу. Требуется найти оптимальную байесовскую оценку положения судна по результатам пяти пар измерений до двух маяков с известными координатами (3000; 0) и (0; 3000), единица измерения — метр. Результаты измерений представлены в виде векторов (2981, 2978, 2985, 3017, 2993) и (2969, 3002, 2977, 3021, 2995) для первого и второго ориентиров. Ошибки измерений предполагаем одинаково распределенными, гауссов-скими, независимыми и центрированными, т. е. V к ~ #(0,30), г = 1,5 , к = 1,2 . Априорная плотность распределения р(х) вектора координат судна х = (х1, х2) предполагается гауссовской, р(х) = Ы(х, 0, 0, 100, 100, 0). На рис.1 изображена рассматриваемая ситуация: априори судно с вероятностью, большей 0,999 находится внутри окружности радиуса 400 м с центром в начале координат. По результатам измерений оно может находиться либо в районе точки М, либо в районе точки N. Однако вариант с точкой М практически невозможен, т. к. эта точка лежит вне указанной выше окружности. Следовательно, апостериорная плотность должна иметь острый максимум в районе точки N при практически нулевых значениях в остальных точках плоскости. На рис. 2 приведен график под-

Рис. 1. Взаимное расположение судна и точечных ориентиров у1 и у2-результаты измерений расстояний до первого и второго маяков

Рис. 2. График подынтегральной функции £,(х)

ынтегральной функции f2(x) из формулы (4) применительно к данной задаче, эта функция пропорциональна апостериорной плотности. На графике хорошо виден острый максимум над точкой наиболее вероятного положения судна, определяемого по результатам измерений.

Для вычислений использовались два метода: АСВ и вариант метода Монте-Карло, в котором в качестве плотности распределения используется априорная плотность распределения координат судна. Этот метод выбран для сравнения, поскольку он используется для расчета искомой оценки в классической схеме important sampling, часто применяемой при решении навигационных задач. Для краткости назовем этот метод «МК_Априор». Сравнение проводилось по величине

Сравнение методов оценивания

2соуМ(8)

"72" Т^---77—,(8)

, 71 1г 12

которая при фиксированном числе итераций прямо пропорциональна относительной погрешности вычислений. Чем меньше ЯЕ, тем метод эффективнее. Входящие в подкоренное выражение (8) точные величины заменены их статистическими оценками. Вычисления проводились на доверительном уровне 95 % для каждой координаты, точность вычислений — 0,1 м для каждой координаты. Калькуляция останавливалась при условии, что оценки обеих координат находятся в преде-

Метод оценивания МК_Априор АСВ

Оценки координат судна 9,06 9,06

Ч 7,03 7,07

Кол-во итераций МО6 1,M05

Показатель дисперсии RE к 5,53 1,76

х2 7,21 2,29

Корреляция числителя и знаменателя к 0,7 0,25

*г 0,6 0,18

лах заданной точности.

Отсюда легко вывести, что обе координаты находятся в пределах заданной точности с вероятностью, не меньшей 90 %. Для АСВ, ввиду его вероятностной природы, используются средние значения всех статистических показателей, входящих в (8), взятые по итогам 300 просчетов. Количество итераций для АСВ указано с учетом итераций, требуемых для генерации случайной кусочно-постоянной в квадрате К = [-400, 400] х [-400, 400] плотности, количество шагов бисекции — 10000. В таблице приведены результаты численных экспериментов.

Анализируя результаты, видим что АСВ тратит приблизительно в 10 раз меньше ите-

Распределение точек в методе МК_й>приор Распределение вброшенных точек в методе АСВ

Рис. 3. Распределения 5000 точек, сгенерированных методами АСВ и МК_Априор

раций, чем МК_Априор, что подтверждается также сравнением соответствующих величин ЯЕ: (5,53/1,76)2=9,87 и (7,21/2,29)2=9,91. Вместе с тем, видно, что метод МК_Априор обеспечивает большую корреляцию между оценками числителя и знаменателя, чем АСВ, что благоприятно отражается на вычислениях, учитывая соотношение (8). Однако, вследствие уменьшения дисперсий числителя и знаменателя за счет хорошо подобранной плотности распределения генерируемых точек, метод АСВ все же имеет значительно меньшие показатели ЯЕ, чем метод МК_Априор. На рис. 3. изображено распределение генерируемых точек для обоих методов.

На рис. 3 хорошо видно, что метод АСВ сгущает точки сильнее, причем именно в той относительно небольшой подобласти квадрата, где находится максимум апостериорной плотности.

Это обстоятельство и обеспечивает большую эффективность АСВ по сравнению с методом МК_Априор.

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

СПИСОК ЛИТЕРАТУРЫ

1. Арсеньев, Д.Г. Адаптивное управление в стохастических методах вычислительной математики и механики [Текст]/Д.Г. Арсеньев, В.М. Иванов, М.Л. Кореневский.-СПб.:Наука, 2008.-2-е изд., испр. и доп.-423 с.

2. Арсеньев, Д.Г. Анализ эффективности адаптивного метода существенной выборки [Текст]/Д.Г. Арсеньев, В.М. Иванов, Н.А. Берковский//Научно-технические ведомости СПбГПУ -2009. -№ 4 (88). -С. 43-50.

3. Арсеньев, Д.Г. Адаптивный метод существенной выборки при ограниченном числе шагов бисек-ционного процесса [Текст]/Д.Г. Арсеньев, В.М. Иванов, Н.А. Берковский//Научно-технические ведомости СПбГПУ-2010.-№ 2-С. 59-68.

4. Чепанов, О.А. Основы теории оценивания с приложениями к задачам обработки навигационной информации [Текст]/О.А. Степанов//Ч. 1. Введение в теорию оценивания. -СПб.: ГНЦ РФ ЦНИИ «Электроприбор», 2009.-496 с.

5. Stepanov, O.A. Investigation of Linear Opti-

mal Estimator [TeKCT]/O.A. Stepanov, A.B. Toropov// Proc. of XVII World Congress. -Seoul, 6-11 July 2008. -P. 2750-2755.

6. Doucet, A. Sequentiual Monte-Carlo methods in practice [TeKCT]/A. Doucet., N. Freitas, N. Gordon. -N.Y. Shpringer-Verlag, 2001.-581p.

7. Kong, A. A note on Importance Sampling using Standardized Weights [TeKCT]/A. Kong//Technical report № 348.-Department of Statistics-Chicago, Illinois: The University of Chicago. 60637.-1992.

8. Geweke, J. Bayesian inference in econometric models using Monte-Carlo integration [TeKCT]/J. Geweke// Econometrica.-1989.-Vol. 57.-№. 6.-P. 1317-1339.

9. Bergman, N. Recursive Bayesian estimation. Navigation and TrackingApplications [TeKCT]/N. Bergman// Linkoping Studies in Science and Technology. Dissertations, 579. -Linkoping, Sweden, 1999.

10. Borenstein, J. Mobile Robot Positioning & Sensors and Techniques [TeKCT]/J. Borenstein, H.R. Everett, L. Feng [et al.]//J. of Robotic Systems, Special Issue on Mobile Robots.-1997.-Vol. 14.-№ 4.-P. 231-249.