Математическое моделирование
морских систем
УДК 551.46.02
И.Е Тимченко, Е.М. Игумнова
Ассимиляция данных наблюдений и адаптивный прогноз природных процессов
Предложена динамическая модель для прогнозирования случайных составляющих природных процессов, основанная на системной концепции адаптивного баланса влияний - Adaptive Balance of Causes (ABC-модель). Модель содержит динамические уравнения для коэффициентов влияний, которые адаптируются к корреляционным связям прогнозируемого процесса. С целью уточнения прогнозов рассмотрены две возможные схемы ассимиляции данных наблюдений в уравнениях ABC-модели: Колмогорова и Калмана. Обе схемы ориентированы на использование выборочных коэффициентов корреляции при прогнозировании временных рядов наблюдений и, следовательно, учитывают нестационарность реальных природных процессов. Приводятся примеры прогнозирования имитированных временных рядов, поясняющие алгоритмы ассимиляции данных наблюдений. Делается вывод о перспективности системного моделирования и адаптивного прогноза случайных процессов ABC-методом.
Введение. Создание динамических моделей природных процессов, которые позволяют получать оценки этих процессов в будущие моменты времени, является одной из основных прикладных проблем геофизических исследований. Природные процессы - сложные объекты моделирования, поскольку наряду с детерминированными составляющими в них, как правило, присутствуют случайные компоненты. Применение общих методов теории случайных функций при построении прогностических моделей дает возможность оценить статистику ошибок прогнозов по отношению к данным наблюдений и использовать эту информацию для улучшения качества прогнозов. На этом основана ассимиляция данных наблюдений [1,2] в моделях природных процессов, которая играет важную роль при построении сценариев развития процессов в интегральных эколого-экономических системах морской среды [3, 4] и в природно-хозяйственных системах прибрежной зоны моря [5].
При интегральном описании морских систем используются наиболее распространенные прогностические модели, предложенные А.Н. Колмогоровым [6] и Р. Калманом [7]. В этих моделях будущие значения процесса формируются его прошлыми значениями под влиянием корреляционных связей, которые предполагаются известными из наблюдений либо могут быть найдены из динамических уравнений, специально создаваемых для этой цели. Оптимальной по точности прогностической оценкой процесса является его среднее значение в будущий момент времени, условное по отношению к на© И.Е Тимченко, Е.М. Игумнова, 2009
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
47
блюдениям процесса в прошлом. Достаточно полный перечень задач, связанных с вероятностными алгоритмами ассимиляции данных океанологических наблюдений, можно найти в работах [8 - 13].
Наиболее реалистичной представляется модель случайной составляющей процесса с переменными статистическими связями между сечениями ансамбля реализаций. Для построения подобной модели необходимо перейти к выборочным корреляциям, зависящим от времени, что становится возможным при накоплении архивных данных наблюдений о реальной динамике природного процесса и оперативной оценке коэффициентов корреляции. В этих условиях происходит непрерывная адаптация прогностических оценок к выборочным или прогнозируемым коэффициентам корреляции.
Для практической реализации этой идеи оказываются полезными системные принципы моделирования природных процессов [14]. В основе системной методологии моделирования лежит принцип адаптивного баланса влияний, согласно которому любая сложная система стремится к установлению динамического баланса с переменными внешними воздействиями, приложенными к этой системе. Будущие значения процесса формируются под воздействием внешних сил, уже оставивших следы своего влияния в данных наблюдений изменчивости процесса в прошлом. Поэтому будущие значения процесса можно рассматривать как результат адаптации к ожидаемым внешним воздействиям на систему и к известным наблюдениям процесса в прошлом.
Заметим, что концепция влияния известных из прошлого данных наблюдений случайного процесса на прогнозируемые будущие значения заключает в себе несколько иной смысл, чем концепция статистической зависимости между двумя соответствующими сечениями ансамбля реализаций, разнесенными по времени. Различие между влиянием и корреляцией такое же, как между физикой и математикой: влияние не может существовать вне времени, так как причина всегда предшествует следствию, тогда как корреляция -мгновенная фиксация того факта, что имеется (в среднем) сходство тенденций в изменчивости двух случайных величин.
На этом соображении основан системный принцип адаптивного баланса влияний, который послужил отправным пунктом при создании одноименного метода моделирования причинно-следственных связей в сложных системах: АВС-метода [14]. Используя этот принцип в задаче прогнозирования природных процессов, предполагают, что будущие значения процесса представляют собой результат адаптивной подстройки под систему влияний, исходящих из прошлого. Этим объясняются применяемые в данной работе термины «системное моделирование» и «адаптивный прогноз».
Принцип адаптивного баланса влияний позволяет построить общую прогностическую модель для любых процессов. Поэтому с методической точки зрения целесообразно рассмотреть возможности системного моделирования природных процессов на имитированных временных рядах наблюдений, имеющих такие же корреляционные связи, какие отмечаются, например, в случайных флуктуациях возвышений уровня морской поверхности или в случайных отклонениях температуры поверхности моря от известного суточного хода.
48 ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
Системное моделирование природных процессов. Модель для прогноза природных процессов должна учитывать в явном виде скорость их изменения. Универсальное правило построения формальных моделей динамических процессов было найдено И. Ньютоном: мгновенная скорость изменения любого процесса есть некоторая функция от самого этого процесса и от приложенных к нему внешних влияний. Поэтому динамическое уравнение для любого природного процесса х(1), находящегося под действием внешнего влияния у(), должно иметь следующий вид:
I=У) ■ ®
Все искусство моделирования состоит в том, чтобы, используя дополнительную информацию о рассматриваемом процессе, найти явный вид функции Е(х, у). Если допустить, что внешние влияния отсутствуют и выбрать наиболее простую зависимость Е(х) = х, то решением уравнения будет экспоненциально растущая функция х = Се1. Такой сценарий развития процесса объясняется наличием положительной обратной связи между скоростью изменения процесса и самим процессом. Поскольку в данном случае они равны друг другу, с увеличением скорости растет значение процесса, что, в свою очередь, еще больше увеличивает его скорость и т.д. Если в качестве правой части уравнения (1) выбрать зависимость Е (х)= — х, это уравнение дает пример отрицательной обратной связи: с ростом скорости процесса его значения будут убывать, и его сценарий с течением времени будет представлен неограниченно убывающей функцией х = Се— ■
Системный принцип адаптивного баланса влияний [14] постулирует существование динамического равновесия в системе, которое должно быть обеспечено балансом положительных и отрицательных обратных связей. Применяя этот принцип, следует учесть в правой части уравнения (1) одновременно обе зависимости х и — х с некоторыми коэффициентами Е^\х) и
Е^+)(х), которые могут быть названы базовыми функциями влияния. Тогда уравнение процесса принимает вид
Л
Потребуем далее, чтобы с ростом х базовая функция монотонно
убывала, а базовая функция монотонно росла. Для обеспечения обще-
го баланса достаточно поставить дополнительное условие симметрии базовых функций влияния:
Е Н(х)+Е(+)(х) = 1. (3)
С учетом этого условия универсальное уравнение системной модели процесса х принимает следующую форму:
— = ЕН(х)х — ^+)(х)х . (2)
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
49
% = х[1 - 2*<+>(х)], (4)
где ^+)(х) - любая монотонно растущая функция. Выбор базовой функции влияния ^+)(х) в уравнении (4) определяется исключительно удобством математического решения этого уравнения. Универсальный характер уравнения объясняется системным принципом баланса влияний, что подтверждает целесообразность использованного выше названия «системное моделирование». Таким образом, системное моделирование дает возможность получить некоторую общую прогностическую модель природного процесса.
Представляет интерес применить наиболее простую базовую функцию влияния, выбрав ее в виде ^+)(х) = х . Тогда уравнение (4) принимает форму нелинейного уравнения Бернулли
йх
— = х(1 - 2х) . (5)
Подобные уравнения (с учетом коэффициентов, уравнивающих размерности) широко используются для моделирования динамических процессов в естествознании. Например, в математической биологии одно их них известно как уравнение динамики популяции или уравнение логистического типа [15]. Общее решение этого уравнения получается в виде
х =-Т^-ТГ , (6)
2 х0 +(1 - 2х0 )е
где через х0 обозначено начальное условие при t = 0. При любых положительных х0 решение стремится к постоянному равновесному значению х = 0,5. Таким образом, процесс, представляемый уравнением (5), приходит в состояние равновесия и находится в нем до тех пор, пока на него не будет оказано внешнее влияние.
В АВС-методе было предложено учитывать внешние влияния путем добавления соответствующих функций в аргументы базовых функций влияния
^+)(х) [12]. Рассмотрим в качестве примера формирование будущего значения процесса хг под влиянием суммы п — 1 известных прошлых значений у ■.
Полагая, что каждое известное значение дает вклад, пропорциональный его величине, т.е. у = а^] , получим
йх,
- = хг
Ж
( п Л
1 — 2^ (+) х — ^
I ]=1
х — У аг]х}
г * ] . (7)
С учетом условия ^+)(х) = х система уравнений (7) принимает оконча-
'(+)(х) = -
тельный вид прогностической АВС-модели для нахождения оценок будущих значений процесса хг:
50
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
Жх,
— = х, Ж
1 - 2
--I
ч ;=!
х - I а1]х]
- * ] . (8)
Как показывают исследования [3 - 5, 14, 16], АВС-модели вида (8) при правильном выборе коэффициентов влияний а^ обладают быстрой сходимостью к стационарным решениям. В условиях изменяющихся внешних воздействий а^х^ происходит непрерывная подстройка решений под сумму этих
воздействий, благодаря заложенной в модели системной концепции адаптивного баланса влияний. Поэтому принципиальное значение имеет объективный выбор коэффициентов влияний а .
Роль внешних воздействий могут выполнять другие процессы, объединенные с процессом х в общую систему причинно-следственных связей и описывающие динамику природной среды. В этом случае уравнения (8) представляют собой динамическую модель природной среды, интегрально описывающую эволюцию ее состояний, представленных множеством процессов {хг}.
Адаптивный прогноз с использованием информации о статистических свойствах процесса в уравнениях АВС-модели. Коэффициенты влияний а в уравнениях модели (8) могут быть идентифицированы при наличии
временных рядов наблюдений за прогнозируемыми природными процессами, выполненных на протяжении некоторого времени и накопленных, например, в архивах. Для статистической оценки коэффициентов уравнений будем предполагать, что временные ряды наблюдений представляют собой случайные отклонения х'- от известного среднего значения. Они вызывают отклонения от этого среднего значения будущих значений процесса, которые являются следствием суммы влияний следующего вида:
п
х;=! ах , - * ]. (9)
]=1
Введем в рассмотрение коэффициенты взаимной корреляции отклонений
Щ = Е{х'кх';}. (10)
Умножая равенство (9) поочередно на х'к, выполняя осреднение полученных выражений и используя обозначения (10), получим
аИ = Кл ^
(
Щ-I а1кЯ,
к к,]
к=1
к * -. (11)
Выражение (11) представляет собой систему линейных алгебраических уравнений для нахождения неизвестных коэффициентов влияний в уравнениях АВС-модели (8). Таким образом, при наличии архивных данных наблюдений для нахождения элементов корреляционной матрицы (10) оказывается
п
п
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
51
возможной объективная оценка коэффициентов влияний а . Заметим, что
уравнения для оценки коэффициентов влияний через коэффициенты корреляции совпадают по форме с известными уравнениями оптимальной интерполяции стационарных случайных функций, полученными Колмогоровым [6].
В качестве примера использования модели (8) - (11) для прогнозирования случайной составляющей природного процесса рассмотрим задачу оценки будущих значений случайного временного ряда наблюдений. В проведенном вычислительном эксперименте были использованы отклонения от д е-терминированной составляющей временного ряда, корреляционная функция которых показана на рис. 1, а. Условимся называть «отсчетами» значения отклонений процесса от среднего значения, «точками отсчета» - моменты времени выполнения наблюдений, а «шагами» по времени - временные интервалы между наблюдениями. Поставим задачу прогноза временного ряда отклонений на 20 шагов вперед по 2, 3 и 4 известным отклонениям, измеренным через 20 шагов. Обоснованность использования подобного интервала между измерениями будет рассмотрена ниже.
Для оценки коэффициентов влияний был выбран случайный временной ряд наблюдений, содержавший 1000 отсчетов, и рассчитана его корреляционная матрица по значениям в 200 начальных точках ряда. Прогнозирование ряда выполнялось с 220-й точки отсчета, как показано на рис. 1, б. Точность прогнозирования характеризуют графики дисперсий ошибок прогнозов, вычисленные путем сравнения оценок, полученных по модели (8) - (11), с истинными значениями ряда. На рис. 1, в эти графики показывают изменчивость выборочных дисперсий ошибок прогнозов по мере перемещения операции прогнозирования вдоль ряда. Дисперсии рассчитывались путем осреднения по 200 текущим точкам ряда и нормировались на выборочные дисперсии истинных значений этого ряда. Рис. 2, а, б позволяют судить о качестве прогнозов, полученных в этом эксперименте.
Точность прогнозирования оказалась различной по мере перемещения группы из четырех измерений вдоль имитированного ряда наблюдений. Если судить по средним вдоль всего ряда величинам относительных дисперсий ошибок, качество прогнозов по трем отсчетам заметно выше, чем по двум. Вместе с тем эффект от добавления данных еще одного (четвертого) измерения оказался не столь существенным.
Постоянные коэффициенты влияний, идентифицированные с помощью системы уравнений Колмогорова (11), имели следующие значения: а12 = 0,15, а13 = 0,21, а14 = 0,5, а15 = 0,5. Судя по их величинам, основную роль в формировании будущих значений процесса играли третье и четвертое измерения, в то время как наибольшие корреляционные связи с точками прогнозов имели первое и второе измерения. Объяснение этого эффекта, возможно, связано с тем, что в этом эксперименте были использованы постоянные выборочные коэффициенты корреляции, построенные по первым 200 точкам ряда. Поэтому представляло интерес провести дальнейшие эксперименты с переменной вдоль всего ряда наблюдений матрицей выборочных коэффициентов корреляции.
52
ISSN 0233-7584. Мор. гидрофиз. журн, 2009, № 6
Х5 Х4 Х3 Х2 X 1
140 160 180 200 220
б
в
Р и с. 1. Постановка задачи прогнозирования: а - корреляционная функция прогнозируемого ряда наблюдений; б - положение точек измерений (черные) по отношению к точке прогноза (белая); в - дисперсии ошибок прогнозов, нормированные на дисперсию прогнозируемого ряда (осреднение по 200 точкам), при использовании данных измерений х1 и х2 (1), измерений хь х2 и х3 (2), измерений хь х2, х3 и х4 (3)
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
53
0.6л
Р и с. 2. Сопоставление прогнозов (1) с истинным процессом (2): а - прогноз по данным двух измерений на 20 шагов вперед с интервалом 20 шагов между измерениями; б - тот же прогноз по данным четырех измерений (см. рис. 1, б)
Построение адаптивной прогностической модели с переменными коэффициентами влияний. Из общего алгоритма оценки коэффициентов влияний (11) следует, что они не зависят от данных наблюдений, используемых при прогнозе, так как определяются только значениями коэффициентов корреляции между точкой прогноза и точками измерений [17]. Если изменить схему расположения измерений, изменятся коэффициенты корреляции и, как следствие, изменятся весовые коэффициенты экстраполяции ряда наблюдений. Эту причинно-следственную зависимость можно выразить с помощью
54 ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
метода адаптивного баланса влияний следующим образом. В соответствии с системой уравнений Колмогорова (11) коэффициенты влияний в уравнениях АВС-модели (8) подстраиваются под схему корреляционных связей, заданную выбором точек расположения измерений. Поэтому по отношению к ним может быть вновь применен системный принцип адаптивной подстройки под сумму внешних воздействий, в качестве которых теперь выступают коэффициенты корреляции, заданные выбором расположения точек измерений. Это дает основание для построения динамической АВС-модели коэффициентов влияний следующего вида [14, 16]:
Уравнения (12) позволяют прогнозировать коэффициенты влияний с использованием для этой цели корреляционной матрицы системы измерений x ■ и введением предположения о том, что коэффициенты корреляции сохраняют свои значения на интервале времени прогноза будущих значений процесса x . Здесь и в дальнейшем изложении мы опускаем штрихи в обозначениях отклонений x' от среднего значения.
Таким образом, уравнения (12) вводят переменные коэффициенты влияний в АВС-модель (8) и, следовательно, открывают возможность ассимилировать в модели информацию об изменчивости статистических связей между точкой прогноза и точками измерений. Иными словами, прогностическая модель (8), (12) становится адаптивной, поскольку в ней используются выборочные коэффициенты корреляции, учитывающие нестационарность случайных составляющих природных процессов, которая чаще всего имеет место на практике. Для того чтобы оценить роль переменных коэффициентов влияний при адаптивном прогнозе вероятностных процессов, были проведены вычислительные эксперименты, результаты которых приведены на рис. 3 и 4.
Случайный ряд отклонений, изображенный на рис. 2, был использован дважды для прогноза его значений сначала на 40, а затем на 60 шагов вперед. Прогнозирование выполнялось на каждом шаге вычислений, начиная с 200-го отсчета. Для получения прогностических оценок были использованы четыре последовательных отсчета ряда, расположенных через 20 шагов. На каждом временном шаге рассчитывалась корреляционная матрица по 200 точкам ряда, которые предшествовали группе из четырех измерений, используемых для прогнозов. Выборка из 200 отсчетов перемещалась вдоль ряда вместе с 4 отсчетами, по которым выполнялись прогнозы. Начальные 200 точек имитировали имеющиеся архивные данные. Последующие значения имитировали поступление новых данных наблюдений, которые пополняли архив для пересчета матрицы выборочных корреляций на каждом шаге. Поэтому в первом эксперименте находились прогностические оценки в 240, 241, ... , 900 точках ряда, а во втором - в 260, 261, ... , 900 точках.
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6 55
к Фi . (12)
0.6-1
в
Р и с. 3. Адаптивный прогноз на 40 шагов вперед: а - исходный ряд (1), прогноз (2); б - динамика коэффициентов влияний в процессе адаптации; в - дисперсия ошибки прогноза на 40 шагов вперед, нормированная на дисперсию прогнозируемого ряда (осреднение по 200 точкам)
56
ISSN 0233-7584. Мор. гидрофиз. журн, 2009, № 6
0.6 -|
в
Р и с. 4. Адаптивный прогноз на 60 шагов вперед: а - исходный ряд (1), прогноз (2); б - динамика коэффициентов влияний в процессе адаптации; в - дисперсия ошибки прогноза на 60 шагов вперед, нормированная на дисперсию прогнозируемого ряда (осреднение по 200 точкам)
ISSN 0233-7584. Мор. гидрофиз. журн, 2009, № 6
57
На каждом шаге решалась следующая система уравнений с переменными коэффициентами влияний, которая вытекает из общей модели (8), (12):
йх1 йг
= х
1 - 2
4
х1
V J=1
а1]+1х]+1
йа
йг
12 = а^1 - 2
(
а12 Я22
Я12 а1,к +1Як
12
V к=2
•к+1,2
йа13 йг
1 - 2
а - Я -
а13 Я33
Я
13
а1
к=1 к * 2
а к+1Як+1,3
(13)
йа14 йг
1 - 2
а - Я -
а14 Я44
Я
14
а1,
к=1 к * 3
а к+1Як+1,4
йа15
йг
1 - 2
а15 - Я55
Я15 а1
к+1Як+1,5
к=1
Система уравнений АВС-модели (13) обеспечивает адаптацию коэффициентов влияний а друг к другу и к переменным элементам корреляционной
матрицы Яи на каждом шаге вычислений. В свою очередь, переменные коэффициенты влияли на формирование прогностических оценок, получаемых по уравнению для х1, тем самым обеспечивая адаптивный прогноз процесса.
Динамика адаптации коэффициентов влияний показана на рис. 3, б и 4, б. Из рисунков следует, что эти коэффициенты испытывали существенные изменения, обусловленные нестационарностью выборочных коэффициентов корреляции. На каждом шаге прогнозирования происходило перераспределение ролей каждого из четырех использованных измерений в формировании прогностической оценки процесса (см. рис. 1, а). Увеличение заблаговремен-ности прогноза на 20 шагов во втором эксперименте не привело к значительным изменениям в графиках каждого из коэффициентов влияний. И в первом и во втором экспериментах наибольшее значение для формирования прогнозов имел отсчет х2, так как его коэффициент влияния принимал максимальные значения. Эти результаты позволяют предположить, что определяющую роль в формировании прогностической оценки играет не близость расположения точки измерения к точке прогноза, а совместный (эмерджентный [16]) эффект системы статистических связей, представляемых матрицей выборочных корреляций.
= а13<
\
/
= а14<
\
/
= а15<
V
/
58
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
Качественное сравнение прогнозов на 40 и 60 шагов с истинными значениями процесса можно выполнить по их временным графикам, приведенным на рис.3, а и 4, а. Количественные оценки - выборочные дисперсии ошибок прогнозов, нормированные на соответствующие дисперсии исходного ряда, показывают заметное отличие в точности прогнозирования (см. рис. 3, в и 4, в). Если при прогнозе на 40 шагов эти оценки колеблются в интервале значений 0,2 - 0,6, то при увеличении заблаговременности до 60 шагов они перемещаются в интервал 0,5 - 0,9.
Ассимиляция данных наблюдений по схеме оптимального фильтра Колмогорова. Последовательное прогнозирование случайной составляющей природного процесса в точки будущих наблюдений дает возможность с течением времени накопить временные ряды ошибок или невязок прогнозов по отношению к данным наблюдений. Изучение статистики невязок позволяет использовать заключенную в них информацию для уточнения прогностических оценок процесса. Рассматривая невязки прогноза в точках наблюдений как «измерения» компоненты ошибок прогноза, можно выполнить прогноз их значений и добавить прогностические оценки невязок к прогностическим оценкам самого процесса. Эта операция ассимиляции данных наблюдений должна увеличивать точность прогностических оценок процесса за счет новой информации, содержащейся в измерениях невязок.
Оптимальный прогноз невязок должен быть обеспечен корреляционной матрицей невязок. Естественно использовать при этом метод оптимальной интерполяции и назвать эту процедуру ассимиляцией данных наблюдений по схеме оптимального фильтра Колмогорова. В рассмотренном выше случае прогнозирования временного ряда на 40 шагов вперед по данным четырех измерений корреляционная матрица невязок прогнозов могла быть рассчитана непосредственно по ряду наблюдений ошибок прогнозов, накопленному за некоторый период времени, например, за 200 шагов. Поэтому применение схемы Колмогорова, основанной на АВС-модели (13), не представляло труда.
Поставим своей целью уточнить результаты прогнозирования путем усвоения четырех невязок прогнозов, которые появлялись в точках отсчетов (см. рис. 1, б) всякий раз, когда наступал момент времени прогноза. Корреляционная функция ряда невязок прогнозов, рассчитанная по первым 200 точкам ряда, показана на рис. 5, а. Как видно из этого рисунка, интервал корреляции ряда невязок оказался существенно короче, чем у исходного ряда (см. рис. 1, а), так как первое пересечение нулевой отметки сместилось с 90-го на 22-й шаг по времени. С учетом установленного интервала прогнозирования на 40 шагов вперед было необходимо увеличить интервал корреляции невязок с тем, чтобы его величина превышала интервал прогнозирования. С этой целью к ряду невязок прогнозов была применена операция скользящего осреднения [17] по 250 точкам. Эта операция увеличила интервал корреляции ряда невязок, так как первое пересечение нуля графиком корреляционной функции сместилось на 57-й шаг по времени (см. рис. 5, а).
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
59
б
Р и с. 5. Статистические характеристики временных рядов наблюдений: а - корреляционные функции исходного ряда (1), невязок прогнозов до предварительной фильтрации (2), невязок прогнозов после предварительной фильтрации (3); б - спектральные плотности исходного ряда (1), невязок прогнозов до предварительной фильтрации (2), невязок прогнозов после предварительной фильтрации (3), (масштаб по горизонтальной оси 5 ■ 10-3 безразмерных единиц)
Как известно, точность прогнозирования случайных процессов зависит не только от того, насколько интервал корреляции ряда превышает заблаго-временность прогноза, но и от правильного выбора интервалов между измерениями, используемыми в прогнозах. В соответствии с теоремой Котельни-кова - Шеннона [18] для каждого случайного процесса существует некоторая
60
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
верхняя граничная частота ас соответствующей ему функции спектральной плотности £ (а), с которой связан максимальный допустимый интервал между отсчетами процесса А/ г Увеличение интервала между измерениями сверх этого значения приводит к появлению среднеквадратичной ошибки прогнозов е(ас) . Ее среднее значение может быть оценено по формуле
ад
ё(ас) = |£(а)с1а . (14)
а
Заметим, что для уточнения прогнозов путем ассимиляции данных наблюдений в этом случае необходимо выполнить операцию предварительной фильтрации (префильтрации [8 - 10, 18]) невязок прогнозов, смысл которой заключается в сглаживании ряда невязок и приведении верхней граничной частоты спектральной плотности этого ряда в соответствие с выбранным интервалом между отсчетами.
Для того чтобы убедиться в правильности выбранного выше интервала в 20 шагов по времени между отсчетами исходного ряда, а также ряда невязок, достаточно построить кривые спектральных плотностей для соответствующих корреляционных функций. Так как корреляционные функции, изображенные на рис. 5, а, относятся к типу экспоненциально-косинусных, все они могут быть аппроксимированы зависимостями вида:
Я(т) = е-аг оовДт. (15)
Оценка параметров в формуле (15) методом наименьших квадратов дала следующие результаты: для корреляционной функции исходного ряда а1 = = 0,01, Д = 0,03; для ряда невязок прогнозов а2 = 0,15, Д = 0,08; для сглаженного ряда невязок а3 = 0,05, Д = 0,018.
Выполняя преобразование Фурье от корреляционной функции (15), можно получить выражение для спектральных плотностей рассматриваемых временных рядов:
£ (а) = — 2л
а а
■ + -
а + (а+Д)2 а + (®-ДТ
(16)
Построенные по этой формуле графики спектральных плотностей приведены на рис. 5, б. По этому рисунку несложно найти безразмерные оценки аы 5 ■ 10-3 - верхних граничных частот спектров: для исходного ряда ас1 =
= 15, для невязок прогнозов ас2 = 40 и для сглаженных невязок прогнозов
ас3 = 20.
По теореме Котельникова - Шеннона интервалы между измерениями процесса, которые используются для восстановления значений этого процесса в произвольной точке, не должны превышать величину А/ которая может быть определена по формуле:
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
61
Поэтому для исходного ряда Д^ = 42 шага по времени, для ряда невязок
Д^2= 16 шагов, а для сглаженного ряда невязок Дор^ = 31 шаг. Выбранный
для прогноза интервал между отсчетами исходного ряда составлял 20 шагов и, следовательно, был меньше предельно допустимого. Однако для ряда невязок он превысил предельно допустимую величину 16 шагов. Поэтому была выполнена префильтрация ряда невязок, которая увеличила предельно допустимую величину интервала между невязками до значения в 31 шаг.
Для реализации алгоритма фильтра Колмогорова на каждом шаге прогнозирования по системе уравнений (13) рассчитывалась корреляционная матрица невязок прогнозов р Ее элементы существенным образом менялись во времени, о чем можно судить по рис. 6, а, на котором показаны временные графики весовых коэффициентов ассимиляции данных наблюдений.
100 200 300 400 500 600 700 б
Р и с. 6. Общие характеристики прогнозов: а - динамика весовых коэффициентов ассимиляции данных наблюдений; б - динамика среднеквадратичных ошибок прогнозов до ассимиляции (1) и после ассимиляции (2)
62
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
Все весовые коэффициенты усвоения невязок прогнозов, за исключением коэффициента £13, демонстрируют значительную изменчивость, обусловленную изменениями выборочных коэффициентов корреляции. Взвешивание с этими коэффициентами четырех усваиваемых невязок прогнозов позволило выполнить оптимальный прогноз невязок на 40 шагов вперед. На рис. 6, б среднеквадратичные ошибки прогнозов с усвоением данных по схеме фильтра Колмогорова сопоставлены с ошибками прогнозов без усвоения данных. Из этого рисунка следует, что весовые коэффициенты ассимиляции, несмотря на внешне неупорядоченный временной ход, в каждый отдельный момент времени были достаточно хорошо согласованы между собой, что и привело к уточнению прогнозов.
Ассимиляция данных наблюдений по схеме оптимального фильтра Калмана. Рассмотрим теперь общий случай, когда информации о фактических ошибках прогнозов недостаточно для непосредственной оценки элементов корреляционной матрицы невязок по ряду наблюдений. Выход из подобного положения был найден Р. Калманом [7], который предложил строить прогностические модели для оценки элементов корреляционной матрицы невязок прогнозов, используя для этой цели динамическую модель исходного процесса. В этом заключается смысл известного фильтра Калмана, который широко применяется в приложениях [2, 8 - 13].
Обозначим через z ■ невязки прогнозов в точках *. , а весовые коэффициенты, с которыми они должны быть добавлены к прогнозу процесса в момент времени , - через . Для коэффициентов корреляции ряда невязок в точках
и ^ будем снова использовать обозначение р . Предположим далее, что для прогноза в точку ^ используется щ отсчетов поля, а для прогноза в точку *. используется щ отсчетов. Чтобы уточнить прогноз в точке ^ производится
ассимиляция данных наблюдений в т точках, предшествующих этой точке.
Тогда на основе прогностической модели процесса (8), (12) может быть построена следующая динамическая АВС-модель ассимиляции данных наблюдений в т точках, реализующая фильтр Калмана через прогнозирование элементов корреляционной матрицы невязок прогнозов в этих точках:
—1 = Р Ж и
1 - 2
Р- -- + 2 а'кК-к + 2 а-1К'1 -22а>ка-1ЕЬ
к=1
I=1
к=1 I=1
- = 1,2,...,п , к Ф г, - ФI,
(17)
Ж
= ё- ^ - 2
£ - Р
ОН -1 11
Р- 2 £гкРк-
к=1
г,- = 1,2,...,т, - Ф г, к Ф г
(18)
opt
Ж
■ = х,
opt
1 - 2
'°Р*- (X +2 ё«2')
1=1
г = 1, 2,...,п, I Ф г, (19)
п
п
п
п
1
т
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
63
ССР'
сргк _ р'
С = Рк
1 - 2
рк рк + ^ § ]кР]к
7=1
г, к = 1,2,..,т, к ф 7 . (20)
Поясним последовательность операций при ассимиляция данных наблюдений по схеме фильтра Калмана в АВС-модели (13), (17) - (20).
1. По известной из анализа архивных данных корреляционной матрице }
процесса х, которая связывает точку прогноза 1г с точками наблюдений 1. , из системы уравнений (12) находятся коэффициенты влияний а^- , а из уравнений (8) вычисляется адаптивная прогностическая оценка х в будущий момент времени.
2. АВС-модель коэффициентов влияний (12) и корреляционная матрица {^ | используются для прогноза элементов корреляционной матрицы {р }
отклонений прогностических оценок процесса от данных наблюдений - невязок прогнозов процесса. Этой цели служат адаптационные уравнения (17), которые обеспечивают подстройку будущих значений матрицы р } под значения элементов матрицы {^ }.
3. С помощью полученной корреляционной матрицы {р } для намеченной к ассимиляции группы данных наблюдений (невязок прогнозов) по системе уравнений (18) вычисляются коэффициенты ассимиляции ^ .
4. Из уравнения (19) получается оптимальная прогностическая оценка процесса х°рР в будущий момент времени 1г, адаптированная к п наблюдениям процесса в прошлом и улучшенная путем ассимиляции т данных наблюдений ошибок (невязок) прогнозов в предшествующие моменты времени.
5. С помощью уравнений (20) корректируются прогностические оценки корреляционной матрицы {р }, которые затем могут быть использованы на
следующем этапе ассимиляции данных. При г = к эти же уравнения позволяют оценивать погрешности ассимиляции данных, связанные в выбором коэффициентов ассимиляции ^ .
Рассмотрим пример ассимиляции данных наблюдений по приведенной выше схеме фильтра Калмана в варианте прогноза процесса на 40 шагов вперед по данным четырех точек измерений, отстоящих друг от друга на 20 шагов. В целях упрощения примем, что для уточнения прогнозов должны быть использованы только две текущих невязки прогнозов: во второй и третьей точках (см. рис. 1, б). Адаптивный прогноз с ассимиляцией данных двух предшествующих прогнозу наблюдений дает уравнение (19):
Сх
С1
— = Х°Р111 - 2 х1^1 - (х1 + §12^2 + §13*3 )]}.
(21)
Весовые коэффициенты ассимиляции §12 и §13 должны удовлетворять системе уравнений (18), которая в данном случае сводится к виду:
т
\
У
64
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
Жё2 = ё12{1 - 2[ё12 - РЛР12 - £13Р32)]},
^13
Л
= ё13{1 -2[ё13 -Р33_1(Р13 - ё12Р23)]}.
Коэффициенты взаимной корреляции невязок прогнозов, необходимые для решения этой системы уравнений, находим из адаптационных АВС-уравне-ний (17) следующего вида:
ЖР
Ж
12 = Р12{1 - 2[Р12 - Р12 + а23Р13 + а12(Р22 - а23Р23) + а13(Р32 - а23Р33)]} ,
ЖР
Ж
13 = Р13{1 - 2[Р13 - р13 + а32Р12 + а12(Р23 - а32Р22) + а13(Р33 - а32Р32)]} ,
йР
Ж
32 = Р32{1 - 2[Р32 - Р32 + а23Р33 + а32(Р22 - а23Р23)]} , (23)
ЖР
Ж
22 = Р22{1 - 2[Р22 - Р22 + а23(2Р23 - а23Р33)]} ,
ЖР
33 = Р33{1 - 2[Р33 - Р33 + а32(2Р32 - а32Р22)]} .
Ж
Для нахождения прогностических оценок коэффициентов корреляции из системы уравнений (23) должна быть построена еще одна система уравнений, позволяющая адаптировать коэффициенты влияний а^ к выборочным коэффициентам корреляции исходного временного ряда. С учетом схемы ассимиляции данных двух текущих наблюдений, которые предшествуют прогнозам, система АВС-уравнений для коэффициентов влияний принимает вид:
Жа12
= а12{1-2[(а12 - Р22 1(Р12 = а13{1 - 2[(а13 -Р33-1(Р13 - al2R23)}}, ^ = а23{1-2[(а23- Р33-1(Р23-а21Р13)]} (24)
23 _ _ _
= а23 I1 2[(а23 Ц33 (Ц23 а21Ц13)
Ш
= а21{1 -2[(а21 - РпЛР21 а23Р31)] }
Жа21
= а2Ц1 2[(а21
М
32
= а32{1-2[(а32- Р22 1(Р32- а3^Р21)]}.
Таким образом, системы уравнений (21) - (24) образуют адаптивную прогностическую модель фильтра Калмана, построенную АВС-методом.
ISSN 0233-7584. Мор. гидрофиз. журн, 2009, № 6
65
Результаты применения ний показаны на рис. 7 и 8.
этого алгоритма ассимиляции данных наблюде-
б
Р и с. 7. Изменчивость коэффициентов влияний (а, б) и результаты прогнозирования элементов корреляционной матрицы ряда невязок прогнозов (в)
66
ISSN 0233-7584. Мор. гидрофиз. журн, 2009, № 6
в
Р и с. 8. Результаты ассимиляции данных наблюдений: а - изменчивость весовых коэффициентов ассимиляции g12, g13; б - исходный ряд наблюдений (1), прогноз на 40 шагов вперед по четырем отсчетам ряда с усвоением ближайших двух отсчетов (2); в - среднеквадратичная ошибка прогноза исходного ряда на 40 шагов вперед по четырем отсчетам ряда без ассимиляции данных наблюдений (1), того же прогноза с ассимиляцией данных четырех наблюдений по схеме фильтра Колмогорова (2), того же прогноза с ассимиляцией данных двух наблюдений по схеме фильтра Калмана (3)
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6 67
Рис. 7, а и б дают представление о заметной изменчивости всех коэффициентов влияний а12, а13, а23, а21, а32 в процессе перемещения вдоль ряда наблюдений группы отсчетов, используемых при прогнозировании значений этого ряда. На рис. 7, в приведены результаты прогнозирования элементов корреляционной матрицы ряда невязок прогнозов Р12, Р13, Р32, Р22, Р33. Как следует из этого рисунка, коэффициенты корреляции невязок прогнозов имели общую тенденцию к убыванию по величине по мере увеличения количества ассимилированных данных наблюдений.
Результаты ассимиляции иллюстрируют графики процессов, изображенные на рис. 8. Весовые коэффициенты ассимиляции данных наблюдений £2 и £ 3, показанные на рис. 8, а, имели существенную изменчивость вдоль временного ряда невязок прогнозов. Добавление к прогнозу процесса двух невязок с этими коэффициентами привело к некоторому улучшению прогностических оценок. Достигнутую итоговую точность прогнозов с ассимиляцией данных характеризуют рис.8, б и в. Первый их них дает качественное представление о результатах прогнозирования процесса на 40 шагов вперед с ассимиляцией данных двух наблюдений, второй - показывает, насколько уменьшилась выборочная среднеквадратичная ошибка прогнозов при усвоении данных двух наблюдений по схеме фильтра Калмана.
Наилучшим из выполненных вариантов прогнозирования оказался прогноз с усвоением данных четырех измерений по схеме Колмогорова с использованием переменных (выборочных) коэффициентов корреляции ряда наблюдений. Несколько меньшей оказалась точность прогнозирования по схеме фильтра Калмана с ассимиляцией данных двух измерений. Заметим, однако, что в последнем случае корреляционная функция ряда невязок прогнозов рассчитывалась по модели и усваивались данные только двух измерений процесса вместо четырех.
Заключение. Прогнозирование природных процессов состоит в оценке их будущих значений с помощью модели, в которой используется конкретная информация о прошлом поведении процесса, заключенная в данных его измерений, а также информация о динамике процесса и об ожидаемых условиях протекания процесса в будущем. Системное моделирование природных процессов дает возможность строить их динамические модели, используя общий метод адаптивного баланса влияний (АВС-метод). Для процессов, имеющих случайный характер изменчивости, уравнения АВС-метода позволяют применить метод оптимальной фильтрации Колмогорова для получения наиболее точных прогностических оценок. Повышение точности достигается за счет адаптации будущих значений процесса к прогнозируемым статистическим связям между прошлыми и будущими значениями процесса.
Особенностью применения оптимальной фильтрации Колмогорова в АВС-моделях является динамическая модель (12), которая вводит переменные коэффициенты влияний в уравнения АВС-метода и таким путем обеспечивает адаптацию прогнозов к динамике выборочных статистических связей процессов. В проведенных вычислительных экспериментах переменные коэффициенты влияний в АВС-модели (8) позволили ассимилировать в этой модели информацию об изменчивости статистических связей между точкой
68
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6
прогноза и точками измерений, что привело к повышению качества прогнозирования. Результаты экспериментов позволяют сделать следующие общие выводы.
1. Системное моделирование природных процессов, имеющих случайный характер изменчивости, позволяет создавать практические схемы ассимиляции данных наблюдений, основанные на оптимальных методах фильтрации Колмогорова и Калмана. Теоретически эти методы обеспечивают наилучшую точность прогнозирования случайных процессов при условии известных из наблюдений или из модельных расчетов корреляционных матриц процессов.
2. АВС-модель (8), (12) может быть успешно применена для прогноза процессов с ассимиляцией данных наблюдений по схеме оптимального фильтра Колмогорова, когда имеется возможность вычислять выборочные коэффициенты корреляции самого процесса и ряда невязок прогнозов непосредственно по наблюдениям.
3. АВС-модель (8), (12), (17) - (20) позволяет прогнозировать процессы с усвоением одиночных измерений, когда отсутствует возможность непосредственной оценки коэффициентов корреляции невязок прогнозов и статистика невязок вычисляется по динамической модели (12), (17).
Полученные результаты применимы к системам взаимосвязанных процессов, развивающихся в интегральных моделях природной среды. Они могут быть распространены также на случаи ассимиляции данных наблюдений в узлах сеточной области, которая обычно производится при четырехмерном анализе полей океана и атмосферы.
СПИСОК ЛИТЕРАТУРЫ
1. Беляев В.И., Тимченко И.Е. О применении объективного и четырехмерного анализа в океанографии // Морские гидрофизические исследования. - Севастополь: МГИ АН УССР, 1972. - № 2(58). - С. 80 - 92.
2. Тимченко И.Е. Прогнозирование гидрофизических процессов на основе фильтра Калмана // Там же. - Севастополь: МГИ АН УССР, 1973. - № 1(60). - С. 99 - 109.
3. Игумнова Е.М., Тимченко И.Е. Моделирование процессов адаптации в экосистемах// Морской гидрофизический журнал. - 2003. - № 1. - С. 46 - 57.
4. Еремеев В.Н., Игумнова Е.М., Тимченко И.Е. Моделирование эколого-экономических систем. - Севастополь: НПЦ «ЭКОСИ-Гидрофизика», 2004. - 320 с.
5. Иванов В.А., Игумнова Е.М., Латун В.С., Тимченко И.Е. Модели управления ресурсами прибрежной зоны моря. - Севастополь: НПЦ «ЭКОСИ-Гидрофизика», 2007. - 258 с.
6. Колмогоров А.Н. Интерполирование и экстраполирование стационарных случайных последовательностей // Изв. АН СССР. Серия матем. - 1941. - 5. - С. 3 -13.
7. Kalman R.E. A new approach to linear filtering and prediction problems // J. Basic Engen. -Trans. of ASME. - 1960. - March. - P. 35 - 45.
8. Timchenko I.E. Stochastic Modeling of Ocean Dynamics. - Chur - London -Paris -New-York: Harwood Acad. Publ., 1984. - 320 p.
9. Кочергин В.П., Тимченко И.Е. Мониторинг гидрофизических полей океана. - Л.: Гидро-метеоиздат, 1987. - 280 с.
ISSN 0233-7584. Мор. гидрофиз. журн., 2009, № 6 69
10. Тимченко И.Е. Системные методы в гидрофизике океана. - Киев: Наук. думка, 1988. -180 с.
11. Modern Approaches to Data Assimilation in Ocean Modeling / Ed. P. Malanotte-Rizzoli // Elsevier Oceanography Series. - 1996. - 455 p.
12. Robinson A.R., Lermusiaux P.F.J. Overview of data assimilation // Harvard Rep. Phys. Interdiscip. Ocean Sci. - Cambridge, Massachusetts: Harvard University. - 2000. - № б2. - 28 р.
13. Коротаев Г.К., Еремеев В.Н. Введение в оперативную океанографию Черного моря. -Севастополь: НПЦ «ЭКОСИ-Гидрофизика», 2006. - 382 с.
14. Тимченко И.Е., Игумнова Е.М., Тимченко И.И. Системный менеджмент и АВС-техноло-гии устойчивого развития. - Севастополь: НПЦ «ЭКОСИ-Гидрофизика», 2QQQ. - 225 с.
15. Математические модели в биологической океанографии / Под ред. Т. Платта, К.Х. Манна, Р.Е. Улановича. - Париж: Изд-во ЮНЕСКО, 1984. - 195 с.
16. Тимченко И.И., Игумнова Е.М., Тимченко И.Е. Образование и устойчивое развитие. Системная методология. - Севастополь: НПЦ «ЭКОСИ-Гидрофизика», 2QQ4. - 527 с.
17. Яглом А.М. Корреляционная теория стационарных случайных функций. - Л.: Гидроме-теоиздат, 1981. - 280 с.
18. Petersen D., Middlton D. Sampling and reconstruction of wave-number-limited functions in N-dimensional Euclidean spaces // Inform. and Contr. - 1962. - 5. - P. 81 - 104.
Морской гидрофизический институт НАН Украины, Материал поступил
Севастополь в редакцию 10.07.08
После доработки 07.08.08
ABSTRACT Dynamic model for predicting random components of natural processes based on the systems conception of Adaptive Balance of Causes (ABC-model) is proposed. The model contains dynamic equations for cause coefficients which adapt themselves to correlation dependences of the forecasted process. To improve prediction accuracy, two possible schemes of data assimilation in the ABC-model equations, namely Kolmogorov and Kalman ones, are considered. Both of them are supposed to use sample correlation coefficients in prediction of measurement time series and, consequently, take into account non-stationarity of real natural processes. Examples of prediction of the simulated time series are presented to clarify the algorithms of data assimilation. The conclusion is drawn that systems modeling and adaptive prediction of random processes by the ABC-method show considerable promise.
70
ISSN 0233-7584. Мор. гидрофиз. журн, 2009, № 6