УДК 629.7.085.2
АЛГОРИТМ РЕКОНСТРУКЦИИ ФАЗОВОГО ПРОСТРАНСТВА ДИНАМИЧЕСКОЙ СИСТЕМЫ И ЕГО ПРИМЕНЕНИЕ ДЛЯ РАЗРАБОТКИ ПРОГНОЗНЫХ МОДЕЛЕЙ
Г. Н. Мальцева, доктор техн. наук, профессор
A. В. Назарова, канд. техн. наук, доцент
B. Л. Якимова, канд. техн. наук, заместитель начальника кафедры аВоенно-космическая академия им. А. Ф. Можайского, Санкт-Петербург, РФ
Цель: разработка прогнозной модели состояния нелинейной динамической системы на основе реконструкции ее фазового пространства. Результаты: сформулирована общая постановка задачи прогнозирования состояния динамической системы на основе реконструкции аттрактора (предельного множества точек) ее фазового пространства с использованием системы дискретных отображений наблюдаемых параметров. Разработан алгоритм реконструкции фазового пространства динамической системы по временному ряду одного из ее параметров, включающий разбиение исходной временной реализации на обучающую и тестовую выборки и их предварительную обработку; оценку размерности фазового пространства; определение необходимого для реконструкции аттрактора количества точек временного ряда; анализ сечений аттрактора; определение размерности и построение прогнозной модели в виде системы дискретных отображений; оптимизацию модели и ее тестирование. Приводится пример прогнозной модели телеметрического параметра бортовой динамической системы космического аппарата. Для рассмотренного примера определено время упреждения прогноза и продемонстрированы возможности предварительной обработки анализируемых временных рядов с использованием алгоритмов исключения аномальных отсчетов и сглаживания. Практическая значимость: реконструкция фазового пространства динамической системы позволяет аналитическим путем оценить параметры синтезируемой прогнозной модели и ускорить решение задачи оптимизации ее характеристик. Полученные алгоритмы и модель использованы в программном комплексе обработки телеметрической информации.
Ключевые слова — динамическая система, прогнозная модель, фазовое пространство, аттрактор, телеметрическая информация.
Введение
Одной из задач технической диагностики сложных технических систем является задача прогнозирования изменений параметров их технического состояния (ТС) в процессе функционирования. Данная задача решается, например, при испытаниях и летной эксплуатации объектов ракетно-космической техники в процессе обработки телеметрической информации [1, 2]. При этом в динамике функционирования сложных технических систем некоторые параметры их ТС могут иметь сложный нелинейный характер как следствие шумов измерений, а также отказов и неисправностей отдельных функциональных узлов, что требует использования нелинейных прогнозных моделей [3-6]. В настоящей работе рассматривается один из подходов к построению прогнозных моделей параметров ТС динамических систем, относящихся к классу самоорганизующихся систем, основанный на реконструкции фазового пространства, характеризующего изменение их ТС.
Общая постановка задачи разработки алгоритма реконструкции фазового пространства динамической системы
В нелинейной динамике свойства сложных технических систем рассматривают в фазовом пространстве состояний [5-7]. В соответствии
с теоремой Такенса, формулирующей требования к реконструкции фазового пространства динамических систем, их поведение можно описать по временной реализации одного из параметров системы с использованием метода задержек [8-10]. При этом параметр динамической системы должен быть представлен дискретными отсчетами: Хо, Хр Х2, ..., хм-1, где М — количество отсчетов экспериментальной реализации параметра, полученных с интервалом Д£, который обычно выбирается равным времени первого пересечения нуля автокорреляционной функцией временного ряда. Для цифровых телеметрических систем значение Дí кратно периоду дискретизации теле-метрируемого параметра [9, 10].
Под состояниями динамической системы будем понимать точки Sо = [Хо, х^, Х2, ..., х^^],
S 1 = [х1, х2, х3, "•, ХЫ], "•, SМ-N = [хМ-ЛГ, ХМ-Ы+1' Хм-ы+2, "•, Хм—1] в ^мерном фазовом пространстве, где число состояний составляет М - N + 1, а координаты хъ, х^+2, ..., х^+^1 представ-
ляют собой значения параметров й-го состояния и образуют траекторию. Порядок формирования векторов состояния из элементов временного ряда параметра динамической системы должен быть одинаковым для всех й-х состояний. Если анализируемая динамическая система обладает свойством самоорганизации, то с течением времени траектории ее состояний, выходящие из различных начальных точек фазового пространства,
стремятся собраться в некоторых сравнительно небольших областях фазового пространства и в дальнейшем их не покидают. Такие области получили названия аттракторов [9, 10]. Аттрактор характеризует процесс, происходящий в самоорганизующейся динамической системе, который может рассматриваться как выделение некоторого числа параметров (существенных переменных) системы, определяющих ее динамику.
Аттрактор представляет собой множество точек фазового пространства, на которое выходит самоорганизующаяся система при устремлении времени протекания процесса ее функционирования в бесконечность, и может быть использован при прогнозировании ее состояния. Алгоритм реконструкции позволяет найти размерность фазового пространства N, в котором необходимо рассматривать аттрактор (размерность вложения) и получить математическую модель динамической системы. Если модель найдена, то она должна с заданной точностью воспроизводить имеющиеся экспериментальные реализации параметров системы и давать прогноз последующего состояния динамической системы по предыдущему состоянию на заданное время упреждения прогноза Т. Такая прогнозная модель параметра ТС динамической системы имеет вид системы дискретных отображений [10, 11]
Х1+1 = А (Х1, Х1-1 , •••, х1-N+1);
Х1+2 = А2 (х' Х;-1, •••, Хi-N+1 );
xi+N - fN (xi, xi-1' •••' xi-N+1)'
(1)
где х1 — дискретные отсчеты параметра ТС, рассмотренные в моменты времени Ш; I — номер дискретного отсчета параметра ТС, I = N - 1, N, ..., М - N - 1.
Функции в правой части выражения (1) связывают параметры предыдущего и последующего состояний в фазовом пространстве динамической системы, разнесенных во времени на интервал Т = NAt. Данные функции можно аппроксимировать полиномами степени п:
А (х , Х1-1 ,•••, Х1-N+1 ) =
п N-1 N-1
= X А-^-г П Х^, X 4 * п (2)
l0'h.'—'lN-1-0
z-0
В выражении (2) Аи,^^ — коэффициенты и-го полинома, которые могут быть найдены с помощью метода наименьших квадратов; Ьг — составляющие степени компонентов полиномов. В модели (1) размерность фазового пространства N определяет время упреждения прогноза Т = NAt. Параметрами модели (1) являются размерность фазового пространства N и степень полиномов п.
Для повышения точности прогноза зашум-ленных параметров ТС необходимо использовать алгоритмы исключения аномальных отсчетов и сглаживания. В настоящей работе в качестве таких алгоритмов рассматриваются алгоритм исключения аномальных отсчетов, основанный на вычислении соседних разностей в пределах окна заданной размерности, и алгоритм сглаживания на основе преобразования Фурье. Данные алгоритмы обработки были синтезированы как итерационные и имеют следующие настраиваемые параметры: K — количество итераций алгоритма исключения аномальных отсчетов; Fгр — значение граничной частоты (частоты среза) фильтра нижних частот.
Для оценки точности прогноза используются следующие статистические показатели: Ei — среднеквадратическая относительная погрешность между прогнозными и сглаженными значениями параметра; E2 — среднеквадратическая относительная погрешность между прогнозными и исходными несглаженными значениями параметра. Первый показатель характеризует возможности используемой модели (1), а второй — вклад в точность прогноза алгоритмов исключения аномальных отсчетов и сглаживания. Кроме того, для оценки прогноза используются вероятностные показатели достоверности прогноза Gi и G2, представляющие собой вероятности того, что значения Ei и E2 не превысят заданного значения Езадан = 0,05 с доверительной вероятностью
P = о 99
задан
В алгоритме реконструкции фазового пространства должна быть предусмотрена возможность оптимизировать параметры модели (1) и используемых алгоритмов обработки в соответствии с критериями Ei, E2(K, F^, n, N) ^ min и Gi, G2(K, F , n, N) ^ max при ограничении Ei, E2 - -^адан.
Описание алгоритма и результаты моделирования
В качестве примера рассмотрим реконструкцию фазового пространства одной из бортовых динамических систем космического аппарата, обладающей свойством самоорганизации, по значениям временного ряда одного из ее параметров. Данный параметр нелинейный, характеризуется сложным поведением во времени, содержит аномальные отсчеты и шумы. Исходная выборка значений параметра была получена по реальным данным телеметрирования бортовой аппаратуры космического аппарата с интервалом дискретизации At = 3 мин.
Структурная схема алгоритма реконструкции фазового пространства динамической системы по временному ряду анализируемого параметра
и построения прогнозной модели вида (1) представлена на рис. 1.
На шаге 1 исходная временная реализация значений анализируемого параметра разбивается на две — обучающую и тестовую. Первая необходима для разработки модели вида (1), а вторая — для проверки ее работоспособности. Затем выполняется первичная обработка параметра с использованием алгоритмов исключения аномальных
отсчетов и сглаживания случайной составляющей и по исходной временной реализации параметра с использованием ее автокорреляционной функции уточняется значение временного интервала Д^ При необходимости осуществляется прореживание временного ряда значений телеметри-руемого параметра.
На шаге 2 по имеющейся временной реализации анализируемого параметра оценивается
с
КОНЕЦ
3
■ Рис. 1. Схема алгоритма реконструкции фазового пространства динамической системы по временному ряду параметра и синтеза прогнозной модели
размерность аттрактора. При этом используется корреляционная размерность [9]
где С(е) = lim -1— X е(е-|Sk - Sy —показатель, m k,j=o
характеризующий вероятность того, что две точки на фазовой траектории динамической системы будут находиться друг от друга на расстоянии, не превышающем заданного; е — заданное расстояние между парой точек в ^-мерном фазовом пространстве; Sfe, Sj — векторы состояния динамической системы из начала коорди-
D = lim
m-1
ln (C(e)) In (е) '
(3)
— функция Хевисайда;
|<0
m — число точек фазовой траектории. Для определения расстояния е между парами точек в ^мерном фазовом пространстве при вычислении показателя С часто используют евклидову метрику и скорректированную евклидову метрику, при которой евклидово расстояние между парами точек нормируется на размерность фазового пространства N.
Корреляционная размерность (3) получает широкое практическое применение в нелинейной динамике и входит в семейство размерностей Реньи (в это же семейство входят размерность Хаусдорфа — Безиковича, информационная раз-
1, (е- Sk -j
0, (е- Sk -S/|
мерность и др.). При значительном времени наблюдения за объектом можно полагать, что
D
ln (C (е)) ln (е) '
С(е) «-1- Y 0(е-|\Sk -Sy||).
m k,j=o
m-1
На практике для приближенного расчета показателя С(е) используется формула [9, 10]
С(ф
2
m-2 m-1
m(m -1) k=0 y=k+1
X X е(е-||Sk -Sy||).
(4)
Расчет значений С(е) по формуле (4) основывается на определении числа точек фазовой траектории, расстояние между которыми меньше заданной величины е или равно ей. Значения С(е) рассчитываются для различных значений е. Полученная зависимость в логарифмическом масштабе отображается на координатной сетке, причем по оси абсцисс откладываются значения 1п(е), а по оси ординат — значения 1п(С(е)). При малых значениях е данная зависимость имеет линейный участок. Указанная процедура повторяется для различных значений размерности фазового пространства N. Предел полученной величины при устремлении размерности фазового пространства N в бесконечность (на практике N изменяют в диапазоне 1 < N < 10) является корреляционной размерностью D.
Вид в трехмерном пространстве аттрактора исследуемой динамической системы получен по отсчетам анализируемого телеметрического параметра в отсутствие (рис. 2, а) и при наличии
xi-1 0,2
xi-2 0,2
4-2
i-i
■ Рис. 2. Аттрактор исследуемой динамической системы, полученный по исходным (а) и сглаженным (б) значениям параметра
(рис. 2, б) сглаживания и исключения аномальных отсчетов.
Логарифмические зависимости показателя С от заданного расстояния е при различных размерностях фазового пространства N (рис. 3) позволяют определить по сглаженным значениям временного ряда максимальное значение корреляционной размерности аттрактора: Б = 2,7. Кривые соответствуют значениям N = 1, 2, ..., 14; с увеличением размерности N кривые С(е) опускаются ниже.
На шаге 3 определяется необходимое для реконструкции аттрактора динамической системы число точек временного ряда М. Существуют приблизительные оценки данной величины, например М > Мт1п = 10а+рБ, где Б — размерность аттрактора; а = 2, р = 0,4 — эмпирические параметры [12]. При М < Мт1п малое число точек не позволит правильно определить размерность аттрактора и описать его аналитически. Более того, при малом числе точек траектория динамической системы может не иметь аттрактора и восприниматься как случайная. В рассматриваемом случае, применительно к анализируемому параметру, условие выполняется: М и 10 000,
Мшп = 1203.
На шаге 4 проводится визуальный анализ сечений аттрактора (сечений Пуанкаре), который может помочь при поиске размерности пространства вложения для построения модели вида (1). В результате такого анализа можно получить дополнительную априорную информацию о сложности аттрактора динамической системы [9]. Это бывает необходимо в случае различного рода наложений и складок фазового пространства, которые имеют место в реальных системах и влияют на величину корреляционной размерности, определяемой выражением (3).
1п С(е)
N =1^'
N =
щтм = 14
Г 1 |
5,5
-5
4,5
-4
1п £
Рис. 3. Логарифмические зависимости показателя С от заданного расстояния е при различных размерностях фазового пространства
На шаге 5 осуществляется определение размерности пространства вложения N и построение модели вида (1) с полиномами степени п. В соответствии с теоремой Такенса для определения размерности пространства вложения N необходимо использовать следующее соотношение [9]:
N > 2Б +1.
(5)
Однако для реальных динамических систем, параметры которых содержат шум и искажения, оценка (5) является завышенной, и для сильно зашумленных рядов при выборе N целесообразно ограничиться ближайшим к Б сверху целым числом. Исходя из этого для рассматриваемого примера получены две оценки размерности пространства вложения: нижняя оценка Nmin = 3 и верхняя оценка Nmax = 7. В результате экспериментальных исследований установлено, что при N = Nmax существует риск получить плохой прогноз на тестовой выборке, а при N = Nmin имеет место малое время упреждения прогноза. Поэтому целесообразно выбирать среднее значение размерности пространства вложения N = ^^^ + N^^/2. Для рассматриваемого примера соответствующее значение размерности пространства составляет N= 5.
Для найденного значения размерности фазового пространства N строится модель вида (1) с использованием метода наименьших квадратов. На обучающей реализации оцениваются ее характеристики Е^ Е2, О-^, О2. Полученные на очередной ^-й итерации алгоритма значения Е1([, Е25, О^, О2(1 показателей Е^, Е2, О^ О2 сравни-
О2
ваются со значениями Е1д-1, Е2д_1, О1д-1, О2q-\, полученными на предыдущей итерации, и путем изменения степени полиномов п в правой части уравнений модели (1), а также параметров алгоритмов исключения аномальных отсчетов К и сглаживания Егр осуществляется поиск условий оптимизации этих характеристик в соответствии с выбранными критериями.
На шаге 6 полученная на обучающей реализации прогнозная модель проверяется на тестовой реализации. Если характеристики прогноза удовлетворяют заданным условиям (Е1, Е2 < Езадан), работа алгоритма заканчивается. В противном случае необходимо уменьшить размерность пространства вложения N и повторить шаги 5 и 6 заново. Для рассматриваемого примера при N = 5 были найдены оптимальная степень полиномов в правой части уравнений модели (1) п = 2; оптимальные значения параметров алгоритмов исключения аномальных отсчетов и сглаживания К = 1, Егр = 0,25Ед, где — частота дискретизации параметра.
Прогнозные модели получены на тестовой реализации телеметрического параметра исследуемой динамической системы в дискретном
■
a)
0,75
0,7
0,65
0,6
0,55
0,5
0,45
0,4
0,35
0,3
0,25
&
2
1 л
\
* 4 1) 1
) Щ
\
( \
50
100
150
б)
0,65
0,6
0,55
0,5
0,45
0,4
0,35
0,3
0,25
и
1 1 1
1
3
| 1
120 140 160 180 200 220 240
■ Рис. 4. Зависимость значений прогнозных 1 и исходных 2 (а), прогнозных 1 и сглаженных 3 (б) отсчетов параметра ТС от номера отсчета
времени I (интервал дискретизации Дt = 3 мин) без сглаживания (рис. 4, а) и со сглаживанием (рис. 4, б). Приведенные результаты показывают, что разработанная прогнозная модель позволяет получить удовлетворительный краткосрочный прогноз анализируемого телеметрического параметра (Е = 0,05; = 0,92; Е2 = 0,005; в2 = 0,99) на время упреждения Т = 15 мин и может быть использована для многошагового прогноза на время упреждения Т = 30 мин. Дальнейшее увеличение времени упреждения в рассмотренном случае ведет к значительному ухудшению качества прогноза.
Заключение
Рассмотренный алгоритм создания прогнозных моделей динамических систем по результатам анализа их параметров ТС является одним из
Литература
1. Сухорученков Б. И., Меньшиков В. А. Методы анализа характеристик летательных аппаратов. — М.: Машиностроение, 1995. — 368 с.
2. Назаров А. В., Козырев Г. И., Шитов И. В. Современная телеметрия в теории и на практике. — СПб.: Наука и Техника, 2007. — 672 с.
3. Якимов В. Л., Назаров А. В. Прогнозирование технического состояния малых космических аппаратов с использованием многослойных нейронных сетей // Изв. вузов. Приборостроение. 2006. № 1. С. 7-12.
4. Якимов В. Л. Прогнозирование параметров технического состояния стартового комплекса с использованием нейронных сетей // Изв. вузов. Приборостроение. 2006. № 7. С. 7-10.
аналитических подходов к синтезу нелинейных прогнозных моделей. Данный алгоритм может быть использован при обработке измерительной информации, характеризующей ТС сложных технических систем, а также на этапе создания систем со сложной динамикой функционирования. В последнем случае построение прогнозной модели позволяет оценить динамические свойства создаваемой системы и определить критические значения некоторых параметров. Эта информация дает возможность избежать необходимости полного перебора всех возможных значений и уменьшить тем самым время на создание системы. Особенно это важно при больших объемах измерительной информации, что имеет место, например, при испытаниях и в ходе летной эксплуатации объектов ракетно-космической техники.
5. Бутырский Е. Ю. Метод аппроксимации нелинейной динамической системы линейной с переменными параметрами // Научное приборостроение. 2013. Т. 23. № 2. С. 89-98.
6. Борисов Ю. Ю. Построение прогнозирующих моделей динамических систем на основе исследования окрестностей реконструированных аттракторов // Автоматизация и современные технологии. 2007. № 2. С. 32-37.
7. Букреев В. Г. Выявление закономерностей во временных рядах в задачах распознавания состояний динамических объектов. — Томск: Изд-во Томского политехнического университета, 2010. — 254 с.
8. Каток А. Б., Хасселблат Б. Введение в современную теорию динамических систем. — М.: Факториал, 1999. — 768 с.
9. Малинецкий Г. Г., Потапов А. Б. Современные проблемы нелинейной динамики. — М.: Едиториал УРСС, 2002. — 360 с.
10. Анищенко В. С., Астахов В. В. Нелинейные эффекты в хаотических и стохастических системах. — Москва-Ижевск: Институт компьютерных исследований, 2003. — 529 с.
11. Павлов А. Н., Янсон Н. Б., Анищенко В. С. Применение статистических методов при решении задачи глобальной реконструкции // Письма в ЖТФ. 1997. Т. 23. № 8. С. 7-13.
12. Nerenberg M. A. H., Essex C. Correlation dimension and systematic geometric effects // Physical Review A. 1990. Vol. 42. N 12. P. 7065-7074.
UDC 629.7.085.2
A Reconstruction Algorithm for a Dynamic System Phase Space and its Application for Development of Predictive Models
Maltsev G. N.a, Dr. Sc., Tech., Professor, [email protected]
Nazarov A. V.a, PhD, Tech., Associate Professor, [email protected]
Yakimov V. L.a, PhD, Tech., Deputy Head of a Sub-faculty, [email protected]
aA. F. Mozhayskii Military Space Academy, 13, Zhdanovskaia St., 197082, Saint-Petersburg, Russian Federation
Purpose: Development of a predictive model of state of a non-linear dynamic system on the basis of reconstruction of its phase space. Results: There has been formulated general problem definition of predicting a state of a dynamic system on the basis of attractor reconstruction (a limit set of points) of its phase space using a system of discrete mapping of observed parameters. There has been developed a reconstruction algorithm of phase space of a dynamic system according to time series of one of its parameters which includes splitting initial time sampling in training and test samples as well as their preliminary processing; assessment of phase space dimension; definition of an amount of time series points required for attractor reconstruction; the analysis of attractor sections; determination of dimension and creation of a predictive model in a form of a system of discrete displays; optimization of model and its testing. There has been given an example of a predictive model of a telemetric parameter of the spacecraft onboard dynamic system. For the given example time of prediction anticipation has been defined and possibilities of preliminary processing of analyzed time series have been shown using algorithms of exception of abnormal values and smoothing. Practical significance: Reconstruction of phase space of a dynamic system allows to evaluate parameters of a synthesized predictive model and to accelerate the solution of the problem of optimization of its characteristics. The obtained algorithms and the model have been used in a software complex for processing of telemetric information. Keywords — Dynamic System, Predictive Model, Phase Space, Attractor, Telemetric Information.
References
1. Sukhoruchenkov B. I., Men'shikov V. A. Metody analiza kharakteristik letatel'nykh apparatov [Methods of the Analysis of Characteristics of Aircraft]. Moscow, Mashi-nostroenie Publ., 1995. 368 p. (In Russian).
2. Nazarov A. V., Kozyrev G. I., Shitov I. V. Sovremennaia telemetriia v teorii i na praktike [Modern Telemetry in the Theory and in Practice]. Saint-Petersburg, Nauka i Tekhni-ka Publ., 2007. 672 p. (In Russian).
3. Yakimov V. L., Nazarov A. V. Prediction of Technical State of Small Spacecrafts with Usage Multilayer Neural Networks. Izvestiia vuzov. Priborostroenie, 2006, no. 5, pp. 7-12 (In Russian).
4. Yakimov V. L. Neural Networks Usage for Parameters Prediction of Technical State of the Launch Complex. Izvestiia vuzov. Priborostroenie, 2006, no. 7, pp. 7-10 (In Russian).
5. Butyrskii E. Iu. Approximations Method of the Nonlinear Dynamic System Linear with Variable Parameter. Nauchnoe priborostroenie, 2013, vol. 23, no. 2, pp. 89-98 (In Russian).
6. Borisov Iu. Iu. Prediction Models Designing of the Dynamic Systems on the Basis of the Investigation of the Reconstruction Environs for Attract. Avtomatizatsiia i sovremennye tekhnologii, 2007, no. 2, pp. 32-37 (In Russian).
7. Bukreev V. G. Vyiavlenie zakonomernostei vo vremennykh riadakh v zadachakh raspoznavaniia sostoianii dinami-
cheskikh ob"ektov [Detection of Regularities in Temporary Ranks in Problems of Recognition of Conditions of Dynamic Objects]. Tomsk: Tomskii politekhnicheskii universitet Publ., 2010. 254 p. (In Russian).
8. Katok A. B., Khasselblat B. Vvedenie v sovremennuiu teo-riiu dinamicheskikh sistem [Introduction in the Modern Theory of Dynamic Systems]. Moscow, Faktorial Publ., 1999. 768 p. (In Russian).
9. Malinetskii G. G., Potapov A. B. Sovremennye problemy nelineinoi dinamiki [Modern Problems of Nonlinear Dynamics]. Moscow, Editorial URSS Publ., 2002. 360 p. (In Russian).
10. Anishchenko V. S., Astakhov V. V. Nelineinye effekty v kha-oticheskikh i stokhasticheskikh sistemakh [Nonlinear Effects in Chaotic and Stochastic Systems]. Moscow-Izhevsk, Institut komp'iuternykh issledovanii Publ., 2003. 529 p. (In Russian).
11. Pavlov A. N., Ianson N. B., Anishchenko V. S. Application of Statistical Methods at the Solution of a Problem of Global Reconstruction. Pis 'ma v ZhTF, 1997, vol. 23, no. 8, pp. 7-13 (In Russian).
12. Nerenberg M. A. H., Essex C. Correlation Dimension and Systematic Geometric Effects. Physical Review A, 1990, vol. 42, no. 12, pp. 7065-7074.