CC BY
102
УДК 519.876.5 ББК 22.161.61 В 58
Власенко А.В.
Кандидат технических наук, доцент кафедры компьютерных технологий и информационной безопасности института компьютерных систем и информационной безопасности, начальник управления организации научных исследований Кубанского государственного технологического университета, Краснодар, e-mail: [email protected]
Жданов А.А.
Аспирант кафедры компьютерных технологий и информационной безопасности института компьютерных систем и информационной безопасности Кубанского государственного технологического университета, Краснодар, e-mail: [email protected] Сизоненко А.Б.
Доктор технических наук, доцент, начальник кафедры информационной безопасности Краснодарского университета МВД России, подполковник полиции, Краснодар, e-mail: [email protected]
Разработка алгоритма прогнозирования временных рядов нелинейных динамических систем
(Рецензирована)
Аннотация. При анализе временных рядов главной задачей является реконструкция породившей этот ряд динамической системы. Задача построения модели динамической системы по одномерному временному ряду решается методами усвоения данных: вариационным или с помощью Калмановской фильтрации. Работу каждого шага фильтра Калмана можно разделить на два этапа: прогноз и корректировка. Этап прогноза вычисляет вектор состояния системы по его же значению на предыдущем шаге работы фильтра. На этапе корректировки в алгоритм поступают данные текущих измерений (наблюдений), которые используются для уточнения прогнозного значения вектора состояния и вычисления собственно оценки вектора состояния динамической системы. Алгоритм последовательно обрабатывает вновь поступающие векторы измерений, учитывая при этом значения, вычисленные на предшествующем цикле. На следующем шаге с помощью обрабатываемых на данном цикле измерений уточняются начальные условия. Для этого алгоритм вычисляет вес поправок к ним на основе ковариационных матриц оценки состояния и измерений. Исправленные таким образом начальные условия и являются выходными данными фильтра Калмана на каждом цикле.
Ключевые слова: временной ряд, нелинейная динамическая система, прогнозирование, алгоритм, фильтр Калмана, бурение скважин.
Vlasenko A.V.
Candidate of Technical Sciences, Associate Professor of Computer Technologies and Information Security Department, Institute of Computer Systems and Information Security, Head of Research Organization Department, Kuban State University of Technology, Krasnodar, e-mail: [email protected]
Zhdanov A.A.
Post-graduate student of Computer Technologies and Information Security Department, Institute of Computer Systems and Information Security, Kuban State University of Technology, Krasnodar, e-mail: Zhdanovan-drey234@mail. ru Sizonenko A.B.
Doctor of Technical Sciences, Associate Professor, Head of Information Security Department of the Krasnodar University of the MIA of Russia, Police Lieutenant Colonel, Krasnodar, e-mail: [email protected]
Development of the algorithm of forecasting time serie of nonlinear dynamic systems
Abstract. The main task in the analysis of time series is the reconstruction of the dynamical system that generated this series. The problem of constructing a model of a dynamical system from a one-dimensional time series can be solved by data assimilation methods: variational or by Kalman filtering. The work of each step of the Kalman filter can be divided into two stages: forecast and correction. The forecast phase calculates the state vector of the system by its value in the previous step of the filter operation. At the stage of correction, the algorithm receives data from current measurements (observations), which are used to refine the predicted value of the state vector, and calculate the actual estimate of the state vector of the dynamic system. The algorithm sequentially processes newly arriving measurement vectors, taking into account the values calculated on the previous cycle. In the next step, the initial conditions are refined by means of the measurements being processed on this cycle. For this, the algorithm calculates the weight of corrections to them based on the covariance matrices of the state and measurement estimates. The initial conditions thus corrected are the output of the Kalman filter on each cycle.
Keywords: time series, nonlinear dynamic system, forecasting, algorithm, Kalman filter, drilling of the wells.
При анализе временных рядов главной задачей является реконструкция породившей этот ряд динамической системы. В соответствии с теорией Такенса-Мане приемлемое описание фазового пространства динамической системы можно получить, если взять вместо реальных переменных системы (которые могут быть неизвестны вообще) ^-мерные векторы задержек, составленные из значений ряда в последовательные моменты времени [1].
Решение обратной задачи нелинейной динамики, то есть восстановление в заданном классе системы дифференциальных или разностных уравнений по скалярному временному ряду наблюдаемого процесса можно разделить на три основных этапа:
1. Построение базиса независимых переменных, временная эволюция которых удовлетворяла бы некоторой системе дифференциальных уравнений.
2. Построение общей модели исследуемой системы в выбранном классе дифференциальных уравнений и при заданном типе нелинейности.
3. Редукция и уточнение модельных уравнений на основе требования глобальной устойчивости исследуемой системы.
В работах по восстановлению уравнений динамики излагается метод построения базиса производных динамического процесса, которые в качестве независимых переменных используются при восстановлении динамической системы по одномерному временному ряду.
Из результатов (2-4) [2] вытекает следующий алгоритм построения временной последовательности данных с задержками:
1. Исходя из данной временной последовательности, построить корреляционную функцию, рассматривая последовательно возрастающие значения размерности фазового пространства п.
2. Получить наклон с1 вблизи начала координат и оценить характер изменения результата при возрастании п.
3. Если величина с1 в зависимости от п выходит на плато выше некоторого относительно небольшого п, то представленная данной временной последовательностью система должна иметь аттрактор. Вышедшую на насыщение величину с1 следует рассматривать как размерность аттрактора, представленного временной последовательностью. Значение п, выше которого наблюдается насыщение, представляет собой минимальное число переменных, необходимых для моделирования поведения, соответствующего данному аттрактору [3].
Согласно [2], система линейных уравнений
ИИ] 2 + Иг гИц + Иг 2И ]3 + ИИ 2 = 0
ИИ +И3И1 = 0
ИпИц +И2И 2 +И3И 3 = 0, (1)
+и2+И3 =1
■ц
ß2n 2 +ß2ß = 1
позволяет построить по одномерному временному ряду модель динамической системы, где задача сводится к поиску минимума функционала вида:
® = Z
(% - F (Хх> У> z>а dt
(2)
где N0 - число точек исследуемого временного ряда.
Задача (2) решается методами усвоения данных: вариационным или с помощью Калма-новской фильтрации.
Фильтр Калмана использует известную математическую модель динамики объекта, которая описывает, какие вообще изменения состояния объекта возможны, чтобы устранить погрешности измерения и представить с хорошей точностью положение объекта в данный момент (фильтрация), в будущие моменты (предсказание) или в какие-то из прошедших моментов (интерполяция или сглаживание). Теоретическим основам и истории развития фильтров Калмана посвящены работы [4-6].
2
1=1
Работу каждого шага фильтра Калмана можно разделить на два этапа: прогноз и корректировка.
Этап прогноза вычисляет вектор состояния системы по его же значению на предыдущем шаге работы фильтра.
На этапе корректировки в алгоритм поступают данные текущих измерений (наблюдений), которые используются для уточнения прогнозного значения вектора состояния и вычисления собственно оценки вектора состояния динамической системы.
Условием оптимальности построенной оценки состояния является минимум ее средней квадратической ошибки. Алгоритм последовательно обрабатывает вновь поступающие векторы измерений, учитывая при этом значения, вычисленные на предшествующем цикле.
На следующем шаге с помощью обрабатываемых на данном цикле измерений уточняются начальные условия. Для этого алгоритм вычисляет вес поправок к ним на основе ковариационных матриц оценки состояния и измерений. Чем меньшей погрешностью характеризуются измерения по сравнению с оценкой состояния системы, тем больший вес они получают. Относительные веса неизвестных, определяющих вектор состояния системы, зависят от степени их влияния на вектор измерений: больший вес получат те переменные, вклад которых в измерения больше.
Уточнение начальных условий на основе поступивших на данном цикле измерений в общем случае приводит к уменьшению неопределенности в оценке состояния системы. Исправленные таким образом начальные условия и являются выходными данными фильтра Калмана на каждом цикле.
На заключительном этапе работы алгоритма происходит подготовка к поступлению нового вектора измерений. На основе заданного линейного преобразования, связывающего последующий вектор состояния с предыдущим, прогнозируется оценка состояния системы, отнесенная к моменту следующего измерения. При построении ковариационной матрицы прогнозируемого вектора состояния фильтром Калмана учитывается возможность искажения модели, описывающей поведение системы, некоторым случайным процессом с известными статистическими параметрами. Поскольку конкретные значения возмущающего эффекта не могут быть известны, данное обстоятельство способствует повышению неопределенности прогноза. По мере последовательной обработки новых измерений происходит накопление фильтром полезной информации, поэтому если элементы вектора состояния уверенно выражаются через измеренные величины, то суммарная погрешность оценок, как правило, должна снижаться. Однако поскольку вместе с улучшением точности оценок на этапе их уточнения имеет место ее снижение при построении прогноза, то эти тенденции, компенсируя друг друга, впоследствии приведут к стабилизации неопределенности, характеризующей оценку состояния системы. В случае отсутствия фактора, вносящего возмущения в процесс перехода системы из одного состояния в другое, погрешность оценок в итоге достигнет нуля. Изменяющаяся в процессе работы алгоритма степень неопределенности оценки состояния системы влечет за собой и изменение весов, вычисляемых на втором шаге; данное обстоятельство выделяет фильтр Калмана как алгоритм с переменными весами.
Если состояние рассматриваемой системы неизменно, то алгоритм фильтра Калмана сводится к последовательной форме классического метода наименьших квадратов, в котором матрица, обратная ковариационной, выступает в качестве весовой. Другими словами, фильтр Калмана является, по существу, рекуррентным способом решения задачи уравнивания по методу наименьших квадратов [7].
Рассмотрим модифицированный ЬБТКБ-алгоритм [8].
Вводим обозначение уп = х^), ^етп, для траектории в пространстве состояний моде-
ли в окне т =
At At
h - у; к +—
где tn - время анализа, х(^) - вектор конечномерного образа состояния системы на сетке модели.
Имеем два входа ЬБТКБ-алгоритма - вектор наблюдений у°п и ансамбль К анализов
х^, к = 1,2,...,К , полученных для момента предыдущего анализа = 1п - Аt.
На шаге прогноза каждый участник ансамбля интегрирован вперед на интервал време-
3
ни —Аt, с использованием К членов ансамбля анализов ха(!), к = 1,2,...,К , как началь-
2 п-1
ных условий для получения ансамбля траекторий первого приближения к\ к = 1,2,...,К .
В конкретной реализации алгоритма членами ансамбля первого приближения являются траектории 6-циклового прогноза, которые начинаются во время 3-циклового прогноза и заканчивается в момент 9-циклового прогноза относительно tn-1.
Разработка шага обновления состояния основывается на предположении, что мы знаем оператор наблюдения к(у п), который удовлетворяет соотношению
У0 =
к(/П)+sn, t . (3)
Здесь уп является модельным представлением (неизвестной) истинной системы траекторий и еп - вектор гауссовских случайных наблюдений шума с нулевой средней и матрицей ковариаций ошибок Яп.
В описываемой реализации компонент интерполяции по времени И(/п) для всех наблюдений формируется с помощью сохранения траекторий первого приближения ГпЬ(к), к = 1,2,...,К, с разрешением 1 цикл; для получения ансамбля модельных состояний к сохраненным полям модели применяется линейная интерполяция по времени наблюдения с точностью до одного цикла.
Определим локальный вектор состояний х1, формируемый из переменных модели в точке модельной сетки I. Фильтр генерирует К членов ансамбля локальных анализов ха(к), к = 1,2,...,К, при помощи вычисления ансамбля весовых векторов wa(l), к = 1,2,...,К, таким образом, чтобы
ха(к)=х+х>а( к). (4)
Здесь хь означает среднее по ансамблю локальных векторов состояния первого приближения хЬ(к), к = 1,2,...,К, и Хь - матрица ансамбля возмущений первого приближения, к-я колонка которого - это к-й член ансамбля возмущений первого приближения: (хЬ(к) - хь).
Наилучшей оценкой локального состояния является среднее по ансамблю анализов:
ха = хЬ + х№ . (5)
Члены ансамбля глобальных анализов ха(к) и глобальный анализ ха получаются путем объединения локальных анализов ха(к) и ха для всех локализаций I.
Вычисляем вес векторов wa(k), к = 1,2,...,К , и среднее по ансамблю ^ в следующем порядке.
Предварительный (нулевой) шаг. Определяется оператор наблюдений к(уп) .
Шаг 1. Оператор наблюдений к(уп) применяется к каждому члену ^к), к = 1,2,...,К, ансамбля траекторий первого приближения для получения ансамбля уЬ(к), к = 1,2,...,К , прогнозных значений в точках наблюдения.
Вычисляем среднее уь по ансамблю уЬ(к), к = 1,2,...,К , и строим матрицу Уь, принимая за ее столбцы вектора, полученные путем вычитания уь из каждого члена ансамбля уЬ(к), к = 1,2,...,К .
Шаг 2. Выполняется локализация. Для каждой точки сетки I выбираются для усвоения наблюдения, которые, как мы думаем, несут полезную информацию о состоянии системы в точке сетки I. Выбранные наблюдения образуют вектор локальных наблюдений у°. Век- 146 -
тор и матрицы У^ и Щ формируются путем выбора тех компонент вектора и элементов матрицы, которые связаны с выбранным набором наблюдений на I. Шаг 3. Весовой вектор ^^ вычисляется следующим образом:
w = Pа Yb К1 (у0 - yb)
где
ра — 11 ~
(к -1))
Р
■(у? ) R-X
(6)
(7)
р> 1 - фактор множественной ковариации и I - единичная матрица. Шаг 4. Вычисляется матрица Ща = [(( -1)а ^ .
Шаг 5. Весовой вектор
w,
добавляется к каждой строке W"
рующей матрицы являются члены ансамбля весовых векторов w
а ( к )
. Столбцами результи-к = 1,2,...,K .
Эффективной иллюстрацией сложного проявления частной модели поведения коммерческой скорости бурения скважин является бифуркационная диаграмма (дерево), на которой представляется зависимость установившихся значений переменной х от параметра а. Проведенное аналитическое исследование позволяет нарисовать начальный участок дерева. Дальнейшее прямое аналитическое исследование циклов и их устойчивости становится невозможным, что приводит к необходимости численного эксперимента с моделью, представленного на рисунках 1-4.
Рис. 1. Бифуркационная диаграмма при а=1,5
Рис. 2. Бифуркационная диаграмма при а=2
-1
Рис. 3. Бифуркационная диаграмма при а=3,24 Рис. 4. Бифуркационная диаграмма при а=3,97
Для нелинейной системы с диссипацией практически невозможно предсказать конкретный ход ее развития, так как реальные начальные условия никогда не могут быть заданы с абсолютной точностью, а наличие точек бифуркации (ветвления) приводит к тому, что даже малые возмущения могут сильно повернуть направление эволюции.
В работах [2, 9, 10] показано, что даже при четко определенных начальных условиях и однозначно заданных параметрах система может иметь несколько устойчивых состояний, которые она принимает последовательно с периодом, зависящим от параметров. При этом даже незначительное изменение параметров системы приводит к изменению динамического
режима. Переход к новому периоду происходит через область хаотического поведения. Таким образом, любая неопределенность в начальных условиях или параметрах системы не позволяет точно прогнозировать ее состояние в последующие моменты времени.
Одним из способов перехода от детерминистической модели к вероятностной является ансамблевое моделирование [11]. При использовании методов ансамблевого моделирования делается попытка количественно оценить чувствительность ситуации к начальным условиям и таким образом определить степень неопределенности модели, возникающей за счет этой причины.
Примечания:
1. Полунин Ю.А., Тимофеев И.Н. Нелинейные политические процессы. М.: МГИМО-Университет, 2009. 204 с.
2. Власенко А.В. Восстановление модели динамической нелинейной системы по порождаемому ей временному ряду // Политематический сетевой электронный научный журнал Кубанского государственного аграрного университета. 2017. № 129. С. 79-92.
3. Ситникова О.В. Математическая модель динамики фьючерсных контрактов на основе методов теории детерминированного хаоса: дис. ... канд. техн. наук. Томск, 2004. 138 с.
4. Cohn S. An introduction to estimation theory // Journal of the Meteorological Society of Japan. 1997. Vol. 75. P. 257-288.
5. A local ensemble Kalman filter for atmospheric data assimilation / E. Ott, B.R. Hunt, I. Szunyogh [et al.] // Tellus. 2004. Vol. 56A. P. 415-498.
6. Evensen G. Data assimilation. The ensemble Kalman filter. Berlin; Heidelberg: Springer-Verlag, 2007. 279 pp.
7. Гончарова Е.Н. Численное дифференцирование временных рядов с использованием фильтра Кал-мана-Бьюси: дис. ... канд. физ.-мат. наук. Ставрополь, 2004. 124 c.
8. Красюк Т.В. Усвоение данных: конкуренция методов и проблема усвоения спутниковых наблюдений. URL:
http://method.meteorf.ru/publ/tr/tr348/krasuk.pdf (дата обращения: 21.08.2017)
9. Власенко А.В., Жданов А.А. Некоторые аспекты прогнозирования скорости бурения нефтяных и газовых скважин // Актуальные научные исследования в современном мире. 2017. № 5-4 (25). С. 3640.
10. Власенко А.В., Жданов А.А. Дискретные отображения как способ описания нелинейных динамических систем // Актуальные научные исследования в современном мире. 2017. № 5-4 (25). С. 4145.
11. Власенко А.В., Жданов А.А. Методика вероятностного прогноза нестационарных временных рядов // Наука нового времени: от идеи к результату: сб. науч. ст. по итогам междунар. науч.-практ. конф., г. Санкт-Петербург. СПб.: КультИнформ-Пресс, 2017. 192 с.
References:
1. Polunin Yu.A., Timofeev I.N. Nonlinear political processes. M.: MGIMO-University Publishing House, 2009. 204 pp.
2. Vlasenko A.V. Restoration of the model of dynamic nonlinear system from the time series generated by it // Polythematic Network Electronic Scientific Journal of the Kuban State Agrarian University. 2017. No. 129. P. 79-92.
3. Sitnikova O.V. Mathematical model of dynamics of futures contracts on the basis of methods of the theory of deterministic chaos: diss. for the cand. of tech. degree. Tomsk, 2004. 138 pp.
4. Cohn S. An introduction to estimation theory // Journal of the Meteorological Society of Japan. 1997. Vol. 75. P. 257-288.
5. A local ensemble Kalman filter for atmospheric data assimilation / E. Ott, B.R. Hunt, I. Szunyogh [et al.] // Tellus. 2004. Vol. 56A. P. 415-498.
6. Evensen G. Data assimilation. The ensemble Kalman filter. Berlin; Heidelberg: Springer-Verlag, 2007. 279 pp.
7. Goncharova E.N. Numerical differentiation of time series using Kalman-Bucy filter: diss. for the cand. of physics and mathematics degree. Stavropol, 2004. 124 pp.
8. Krasyuk ^V. Data assimilation: the competition of methods and problem of assimilation of satellite observations. URL:
http://method.meteorf.ru/publ/tr/tr348/krasuk.pdf (date of access: August 21, 2017)
9. Vlasenko A.V., Zhdanov A.A. Some aspects of predicting the drilling speed of oil and gas wells // Current Scientific Research in the Modern World. 2017. No. 5-4 (25). P. 36-40.
10. Vlasenko A.V., Zhdanov A.A. Discrete mappings as a way of describing nonlinear dynamic systems // Actual Scientific Research in the Modern World. 2017. No. 5-4 (25). P. 41-45.
11. Vlasenko A.V., Zhdanov A.A. The technique of probabilistic forecast of non-stationary time series // Science of the new time: from idea to result: coll. of scient. articles according to the results of the international scient. and pract. conf., St. Petersburg. SPb.: KultlnformPress, 2017. 192 pp.
