ВЕСТНИК ЮГОРСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА
2013 г. Выпуск 3 (30). С. 3-10
УДК 517.9;004; 551.58
АНАЛИЗ И ПРОГНОЗИРОВАНИЕ ЦИКЛИЧЕСКИХ ВРЕМЕННЫХ РЯДОВ С ИСПОЛЬЗОВАНИЕМ ВЕЙВЛЕТОВ И НЕЙРОСЕТЕВЫХ НЕЧЕТКИХ ПРАВИЛ ВЫВОДА
В. И. Алексеев
1. Введение
Известно, что многие природные явления, климатические переменные: солнечная активность, температура, концентрация пыли в атмосфере, метана CH4, углекислого газа C02, атмосферной пыли, изотопов кислорода 180 и водорода (дейтерия D), а также уровень мирового океана, толщина льда и многое другое изменяются циклически [1-4]. Изменяется циклически и ширина древесных колец [5]. На регулярность изменчивости этих переменных накладываются случайные факторы. Весьма актуальна необходимость эффективного прогнозирования изменений всех циклических, в том числе и климатических переменных. Обычно прогнозирование таких процессов основано на использовании экстраполяции трендов временных рядов динамики [2]. Однако такой традиционный подход не позволяет в полной мере использовать предысторию длинных и плохо структурированных временных рядов, характерных для природных явлений. В связи с этим перспективным является использование нового подхода на основе сочетания методов вейвлет-анализа и нечеткого нейро-сетевого моделирования. Однако вопросы такого подхода в настоящее время не разработаны.
Целью выполненной работы является разработка методических вопросов использования для анализа и прогнозирования сложных циклических временных процессов вейвлетов и нечетких нейросетевых технологий.
Рисунок 1 - График отклонений приземной атмосферной температуры от средней в интервале времени от 1850 до 2009 года по инструментальным измерениям
2. Кратномасштабное вейвлет-разложение
В качестве примера для исследования выбран временной ряд динамики изменчивости атмосферной температуры в интервале времени от 1850 до 2009 года, полученный по инструментальным измерениям [2]. Ее график приведен на рис. 1. Для анализа временного ряда используется кратномасштабное вейвлет-разложение [6, 7], а его прогнозирование производится с использованием алгоритма Малла [8] и адаптивной системы нейро-нечеткого вывода [9]. В работе [1] есть примеры вейвлет-разложения климатических переменных.
В вейвлет-анализе сложные функции f (t) представляются в виде совокупностей масштабирующих функций: 3
3
В. И. Алексеев
W(f(а,4)) = {f,/,„) = } f (t)■//'—t|■ dt.
(1)
В формуле /аЪ =—■ щ -—4 | - вейвлетные функции или «волновые пакеты», в кото-
\а\
а
рой а - называется масштабирующим параметром, at - параметром сдвига, множитель
—1— - параметр нормировки вейвлет-функции, обеспечивающий II/JI = 1.
|а| ’
Прямое вейвлет - преобразование (1) можно рассматривать как разложение сигнала по всем возможным сдвигам и растяжениям/сжатиям сигнала f (t), при этом параметры а и t могут принимать любые значения в пределах областей их определения а, t е R ; задание значения 0 < И << 1 соответствуют очень узким вейвлет-окнам и служат для точной локализованной регистрации высокочастотных составляющих сигнала f (t), при задании |а| >> 1 реги-
стрируются медленные или длинноволновые составляющие сигнала f (t). Таким образом, преобразование (1) при задании вектора а в некотором диапазоне (атт ^ атах ) позволяет получить из исходной функции f (t) множество функций W(f (а, b)) (в дискретной области матрицу), изменяющихся с разными частотами (от высокочастотных до низкочастотных). Этот набор функций далее может быть обработан как в отдельности, так и в совокупности и является основой так называемого кратномасштабного анализа (КМА), предложенного Мал-ла [8]. При таком подходе обработки сигнала f (t) характерные детали, которые могут оставаться незаметными при одном разрешении, легко могут быть обнаружены при другом. Использование вейвлетов для анализа функций связано с тем, что они позволяют обнаруживать такие локальные особенности сигнала (функции), которые не могут быть обнаружены любыми другими методами.
3. Вычисление структуры временного ряда на примере приземной глобальной температуры
Из графика исходного ряда динамики температуры (рис. 1) следует, что в изменении температуры присутствует цикличность, но установить закономерность этой цикличности по этим кривым не представляется возможным. Такая возможность появилась в связи с интенсивным развитием теории и техники вейвлет-анализа и применением этой техники в самых различных областях науки и техники [1]. В работе [1] для анализа климатических данных, полученных в современную инструментальную эпоху и палеоданным, получены картины вейвлет-коэффициентов WT (амплитудно-частотные характеристики (АЧХ)), которые показывают, что долговременная динамика рядов хорошо описываются квазигармонической функцией.
В настоящей работе для анализа рядов динамики температур используются комплекснозначные вейвлеты щ -—4 |, в частности вейвлеты Гаусса (“cgau5”), позволяющие вычис-
I а )
лить фазо-частотные характеристики (ФЧХ) временных рядов.
4
Анализ и прогнозирование циклических временных рядов с использованием вейвлетов...
Рисунок 2 - Структура фазо-частотной характеристики приземной атмосферной температуры в интервале времени от 1850 до 2009 года
□ .4
□
-□.4
-1 .2 -1 .6
График на рис. 2 получен при разложении исходной функции f (t) на 40 уровней детальности (совокупное использование вейвлет-деталей) с использованием комлекснозначной вейвлет-функции “cgau5”.
В таблице 2 приведены результаты анализа фазочастотной характеристики временного ряда температуры (рис. 1), соответствующей современной инструментальной эпохе (с 1850 по 2009 годы) измерений температуры, где х - уровни детальности вейвлет-разложения, средние значения периодов цикличности mT, стандартные ошибки аТ значений периодов цикличностей
Таблица 2
х (уровнидетальности вейвлета) mT, лет аТ, лет
40 54,50 ±3,53
20 30,75 ±6,20
8 13,27 ±4,31
3 4,97 ±1,27
1 3,46 ±1,22
Из полученных изображений АЧХ и ФЧХ и анализа приведенных выше таблиц следуют следующие выводы:
1) ряд динамики изменений приземной глобальной температуры атмосферы является поликвазигармонической функцией времени, в ней содержатся множества цикличностей с периодами от нескольких десятков до нескольких лет, обусловленное динамикой взаимодействия атмосферы Земли с сушей и океанами Земли, с планетами солнечной системы и галактиками Вселенной;
2) изменение глобальной температуры атмосферы является нелинейной функцией, показателем такого вывода являются бифуркации, редко и утроений периодов (частот) вейвлет фазо-частотных характеристик ряда температурных изменений. Из приведенной выше таблицы цикличностей температуры атмосферы Земли следует, что среднее значение отношений предыдущих периодов цикличностей к последующим составляет 2,22 с точностью некоторого значения среднеквадратического отклонения; эти же результаты подтверждают, что атмосфера Земли является самоорганизующейся неравновесной термодинамической системой [11];
3) характерным для частотных характеристик временных рядов динамики температур является то, что в силу влияния различных случайных факторов, удвоение частот, а иногда и утроение частот на разных временных интервалах происходят в разное время, которое проявляется на разных уровнях масштабирующих вейвлет-коэффициентов. 5
5
В. И. Алексеев
4. Прогнозирование временной динамики на примере приземной глобальной температуры атмосферы
Решение задачи прогнозирования динамики температуры атмосферы, которую можно представить как совокупность квазигармонических функций, производится на основе двух пар дискретных вейвлет-преобразований [10]:
W j, k) = vM Z ■f (t)' ^ k(t^ (2)
Ww{ j, k) = j= Zf(t) Jjk(tX (3)
где системы функций y/j k(t) = 2J /2щ(2Jt — k) и pj0,k(t) = 2]°12 p(2Jo t — k) образуют ортонормированные базисы, в которых индекс k определяет положение функций (Pj0jk (t) и Щ j k (t) на оси t, индексы j0 и j - ширины функций PJo,k (t) и Щ j ,k (t) , соответственно, вдоль оси t; В формулах (2) и (3) f (t), Pj k (t) и Щ j,k (t) - функции дискретной переменной t = 0,1,2,...,M — 1, при этом коэффициенты Wp( j0, k) называются коэффициентами приближения или аппроксимации, а WJ(J,k) - коэффициентами деталей; 1 л/M - нормировочный множитель.
В теории вейвлетов [6, 7] реконструкция или восстановление исходного сигнала f (t) по вычисленным масштабирующим и детализирующим коэффициентам производится по формуле
f (t)
ZWpPJo,k) 'P0,k(t) +
1
VM
ZZWj( j , k) 'Jj ,k(t)
j=J0 k
(4)
Если ввести обозначения: a1,: a2,..., am - аппроксимации сигнала на уровнях (1 + m),
d1,d2,...,dm - детали сигнала на тех же уровнях, то исходный сигнал f(t) (4) на нулевом
уровне с некоторой точностью (с точностью заданного числа уровней детализации) может быть представлен в виде равенств [7,10]:
f = а + d = а + d + d = а + d + d + d = ... = a + d + d , +... + d (5)
J 11 221 3321 mm m—1 1 •
В формуле (5) аппроксимации at и детали di с увеличением уровня аппроксимации i стремятся к постоянным числам, как низкочастотные функции, а детали di с малыми значениями i пренебрежимо малы как высокочастотные функции. По сути ai и dt - это временные ряды той или иной степени гладкости.
Идея прогнозирования сложного временного сигнала состоит в том, чтобы на начальном этапе прогнозировать более гладкие аппроксимации и детали функции f (t), прогнозированные значения которых отличны от нуля или отличны от постоянного числа. Затем, на следующем этапе, в соответствии с формулами (5) вычислить комбинации аппроксимаций at и
деталей d i , получая прогнозированное значение временного ряда с точностью отброшенных деталей. 6
6
Анализ и прогнозирование циклических временных рядов с использованием вейвлетов...
В работе для решения задачи прогнозирования временного ряда температуры в качестве
вейвлетной функции w
(t - b \
V а у
использована функция “coif5”; в формулах (2) и (3) парамет-
ры, целые значения j0 и j, заданы в интервале (1 ^ 10) .
Прогнозирование декомпозированных составляющих атмосферной температуры основано на использовании адаптивной гибридной системы нейро-нечеткого вывода ANFIS [9].
5. Адаптивная гибридная сеть нейро-нечеткого вывода и ее применение
Гибридная адаптивная система нейро-нечеткого вывода ANFIS представляет собой многослойную нейронную сеть специальной структуры без обратных связей, в которой используются обычные (не нечеткие) сигналы, веса и функции активации. В нейронной сети выполняются операции суммирования
П
5 = £ w ■ х + b; y = f (5), (6)
i=1
где w, - вес синапса (i = 1,..., П ); b - значение смещения; 5 - результат суммирования; у -выходной сигнал нейрона; f- функция активации нейрона (нелинейное преобразование). В формуле (6) единственный выход у и несколько входов х,, i = 1,..., к, представляют собой
нечеткие лингвистические переменные. При этом термы входных лингвистических переменных в MATLAB, использованные нами, описываются восемью стандартными функциями принадлежности, а термы выходной переменной представляются линейной или постоянной функцией принадлежности (нами использована линейная модель). Суммирование в (6) основано на использовании нечетких операторов, таких как T - нормы, T - конормы или некоторой другой непрерывной нечеткой операции. При этом значения входов, выходов и весов гибридной нейронной сети представляют собой вещественные числа из отрезка [0, 1]. Основная идея, положенная в основу модели гибридных сетей, заключается в использовании существующей выборки данных для определения параметров функций принадлежности, которые лучше всего соответствуют некоторой системе нечеткого вывода. При этом для нахождения параметров функций принадлежности используются известные процедуры обучения нейронных сетей [9].
В гибридной сети ANFIS реализован многоэтапный нечеткий вывод преобразования входных переменных в выходные на основе нечетких правил продукций.
В качестве выборки данных используются полученные на этапе вейвлет-декомпозиции исходного временного ряда атмосферной температуры аппроксимации или детали разных уровней. Каждая из этих данных преобразуется в матрицу обучающей выборки, состоящую из пяти столбцов, первые четыре из которых являются входными, а пятый - выходной (обучающий).
Для генерации структуры системы для каждой входной переменной задаются по 3 лингвистических терма, в качестве типа их функций принадлежности выбирается одна из восьми заданных в системе, в качестве типа выходной переменной задается линейная функция. При решении задач нами выбраны треугольные или гауссовы функции принадлежности. Обучение сети производится с использованием гибридного метода обучения с уровнем ошибки нуль. Гибридный метод обучения представляет собой комбинацию метода наименьших квадратов и метода убывания обратного градиента [9]. Проверка адекватности построенной модели гибридной сети и прогнозирование выбранных составляющих ряда атмосферной температуры производится с использованием стандартной процедуры “evalfis” [9].
Входным аргументом этой процедуры является вектор значений прогнозируемой составляющей атмосферной температуры и сгенерированная в процессе обучения структура сети. Адекватность построенной нечеткой модели гибридной сети производится сравнением вы-
7
В. И. Алексеев
ходных значений модели со значением температуры tm+5 и последующими ее значениями,
которые не была использованы в построении модели гибридной сети. Здесь т - число строк обучающей матрицы.
6. Результаты прогнозирования
В соответствии с алгоритмом прогнозирования, описанным в п. 4 и 5, прогнозировалась одна из аппроксимирующих составляющих функции f (t) и несколько детализирующих. На рис. 3 для примера приведены составляющие al, d4 и dl исходного временного ряда, представленного на рис.1.
Рисунок 3а - Аппроксимирующая составляющая al временного ряда f (t)
0,15
0,15
-0,2
Рисунок 3б - Детализирующая составляющая dl временного ряда f (t)
В результате прогнозирования составляющих al, dl, d6, d5, d4, d3 исходного временного ряда и используя формулу Малла (5), получено спрогнозированное значение ряда
fпрогноз (t) = al(t) + dl(t) + d6(t) + d5(t) + d4(t) + d3(t) . (l)
8
Анализ и прогнозирование циклических временных рядов с использованием вейвлетов...
В формуле (7) составляющие d 2(t) и d1( t) не использованы как малые величины. Спрогнозированное значение этого ряда представлено на рис. 4.
Рисунок 4 - Прогнозированное значение ряда (7) до 2069 года Вейвлет фазо-частотная характеристика ряда (7) представлена на рис. 5
Рисунок 5 - Вейвлет фазо-частотная характеристика ряда (7)
Среднее значение периодов мегациклов временной динамики составляет 57 лет со стандартным отклонением о = 5 лет. Из рис. 5 следует, что современная атмосферная температура находится в фазе ее глобального понижения (с 1988 по 2017 годы). Как видно из рис. 5, в интервале времени с 1850 по 2068 годы периоды похолоданий чередуются периодами потеплений с интервалами 4-5 лет.
На рис. 1б видно, что известное похолодание 40-х годов прошлого столетия приходится к интервалу похолодания с 1925 по 1952 годы, холода в 10 - 20-х годов прошлого столетия приходятся к интервалу потепления с 1899 по 1925 годы, похолодания 60-80-х годов прошлого столетия приходятся также к интервалу потепления с 1952 по 1982 годы.
Аналогичные графики спрогнозированного временного ряда атмосферной температуры и ее вейвлет фазо-частотной характеристики получены для комбинации
fпрогноз (t) = a8(t) + d8(t) + d7(t) + d6(t) + d5(t) + d4(t) + d3(t).
7. Заключение
В работе анализ и прогнозирование циклических временных рядов производится на примере анализа и прогнозирования атмосферной температуры, полученной инструментальным путем с 1850 по 2009 годы. Этот процесс включает следующие этапы:
1) декомпозиция ряда на совокупность кратномасштабных функций с использованием заданной комплексной вейвлет-функции и вычисление вейвлет фазо-частотной характеристики ряда с установлением закономерностей роста и понижения функции; при этом ус-
9
В. И. Алексеев
тановлен бифуркационный характер изменения вейвлет фазо-чатотной характеристики ряда;
2) декомпозиция ряда на аппроксимирующие и детализирующие составляющие с последующим их прогнозированием в среде сети гибридной нейро-нечетких выводов и прогнозирование временного ряда с использованием алгоритма Малла;
3) построение и анализ вейвлет фазо-частотной характеристики прогнозированного временного ряда, позволяющая предсказать периоды роста и понижения атмосферной температуры до 2068 года с некоторой точностью.
Следует заметить, что предложенные метод анализа и прогнозирования сложных по структуре временных рядов динамики могут быть использованы в самых разных областях знаний, науки и техники.
ЛИТЕРАТУРА
1. Федулов, К. В. Структура климатических изменений (по палеоданным и данным инструментальной эпохи) [Электронный ресурс] / К. В. Федулов, Н. Н. Астафьева // ИКИ. -Режим доступа : www.iki.rssi.ru
2. Глобальное потепление [Электронный ресурс]. - Режим доступа : http://ru.wikipedia.org/wiki/Глобальное_потепление
3. Шерстюков, Б. Г. Изменения, изменчивость и колебания климата [Текст] / Б. Г. Шерстю-ков. - Обнинск : ФГБУ «ВНИИГМИ-МЦД», 2011. - 293 с.
4. Котляков, В. М. История климата Земли по данным глубокого бурения в Антарктиде [Текст] / В. М. Котляков // Природа. - 2012. - № 5. - С. 3-9.
5. Тишин, Д. В. Дендроэкология (методы древесно-кольцевого анализа) [Текст] / Д. В. Тишин. - Казань : Казанский университет, 2011. - 33 с.
6. Блаттер, К. Вейвлет-анализ. Основы теории [Текст] / К. Блаттер. - М. : Техносфера, 2004. - 280 с.
7. Дьяконов, В. П. Вейвлеты. От теории к практике [Текст] / В. П. Дьяконов. - 2-е изд. -М. : СОЛОН-Пресс, 2004. - 400 с.
8. Малла, С. Вейвлеты в обработке сигналов [Текст] / С. Малла. - М. : Мир, 2005. -671 с.
9. Леоненков, А. В. Нечеткое моделирование в среде MATLAB fuzzy и TECH [Текст] / А. В. Леоненков. - СПб. : БХВ-Петербург, 2003. - 736 с.
10. Гонсалес, Р. Цифровая обработка изображений [Текст] / Р. Гонсалес, Р. Вуде. - М. : Техносфера, 2005. - 1072 с.
11. Пелюхова, Е. Б. Синергетика в физических процессах: самоорганизация физических систем : учеб. пособие [Текст] / Е. Б Пелюхова, Э. Е. Фрадкин. - 2-е изд. - СПб. : Изд-во «Лань», 2011. - 320 с.
10