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

Вестник РУДН Серия Математика. Информатика. Физика. № 2. 2009. С. 66-73

УДК 519.254.1

Моделирование и прогноз развития лесопожарных ситуаций с использованием векторных авторегрессионных процессов

А. В. Парамонов

Центр Математических Исследований МГТУ «СТАНКИН» Вадковский пер. д. 3А, Москва, Россия, 117936

Проведено исследование статистических свойств совместной выборки ежедневных значений количества очагов лесных пожаров и комплексного метеорологического показателя в Иркутской области за 1969—1988 гг. Построена математическая модель стохастического процесса развития лесопожарных ситуаций с использованием векторного авторегрессионного процесса ВАР(1). Разработан вычислительный алгоритм моделирования будущих значений количества очагов возгорания. Показано, как он может быть использован для прогноза развития лесопожарной ситуации с высокой степенью надёжности на 1—2 суток вперёд.

Ключевые слова: лесные пожары, комплексный показатель Нестерова, векторный авторегрессионный процесс, величина риска.

1. Введение

Лесные пожары представляют собой масштабное стихийное природное явление, которое регулярно наносит серьёзный урон российской экономике. Процесс развития лесных пожаров практически не поддаётся контролю со стороны человека.

В то же время существуют и активно используются практические методики профилактики и противодействия лесным пожарам. Эффективность применения этих методик самым непосредственным образом зависит от умения моделировать и предсказывать развитие лесопожарных ситуаций. Точный и своевременный прогноз даёт возможность сосредоточить ограниченные ресурсы на тех мероприятиях, которые позволят минимизировать потенциальные последствия пожара.

2. Показатель Нестерова

Общие требования по мониторингу и прогнозированию лесных пожаров и чрезвычайных лесопожарных ситуаций регулируются ГОСТ Р 22.1.09-99 «Безопасность в чрезвычайных ситуациях. Мониторинг и прогнозирование лесных пожаров. Общие требования» от 25 мая 1999 г. [1]. Степень пожарной опасности в лесу по условиям погоды должна определяться по принятому в лесном хозяйстве комплексному показателю В. Г. Нестерова, который вычисляется на основе данных о температуре воздуха, температуре точки росы и количестве выпавших осадков:

п

КП = £ Т(Т - То),

¿=1

где Т — температура воздуха в градусах Цельсия, Т^ — температура точки росы в градусах Цельсия, п — число дней после последнего дождя (осадки до 2, 5 мм в расчёт не принимаются) [2]. В настоящей работе расчёт значений точки росы произведён на основе имеющихся значений температуры и влажности воздуха по формулам:

т = Ь1(Т,КН) Л а - 7 (Т, ЯН),

Статья поступила в редакцию 24 ноября 2008 г.

7 (T,RH)

аТ b + Т

+ 1п(ЛЯ/100),

(2)

где ЯН — относительная влажность в процентах, а = 17, 27, Ь = 237, 7°С. Формулы (1)—(2) обладают погрешностью не более 0,4 °С в следующем диапазоне значений [3]:

0°С <Т< 60°С, 1% <ПН< 100%, 0°С <Та < 50°С.

3. Модель МАР(1)С

Проведено исследование статистических свойств совместной выборки ежедневных значений количества очагов возгорания (далее — ОВ) и комплексного показателя В. Г. Нестерова (далее — КП) в Иркутской области за 1969-1988 гг.:

- (кпГ,) , к = ркр

Xt, —

где р = 20 — число сезонов (периодов), кр = 214 — число наблюдений за сезон (наблюдения проводились с 1 апреля по 31 октября)1 .

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

m

о ft

се *

о В о Ч о Н

О

00

Jll^l!^ Mil likiiMl

1л liulilU LuiiU

"1 I г

1970 1975 1980

о о о

се m о ft

К

ч и К

1985

Рис. 1. Исходные данные

К полученным значениям сезонной компоненты применена сплайн-аппроксимация. Это позволило рассчитать значения сезонной компоненты для дат, отсутствующих в исходной выборке (до 1 апреля и после 31 октября).

Стохастическая компонента получена путём вычитания сезонной компоненты из логарифмированных данных:

к mod кр

Ytk = 1n№fc) - s(Tk), тк =---р.

Кр

1 Данные по метеорологическим показателям (температура воздуха, относительная влажность, количество осадков) предоставлены Государственным учреждением «Всероссийский научно-исследовательский институт гидрометеорологической информации — Мировой центр данных». Данные по количеству очагов возгорания предоставлены предприятием «Авиалесоохрана».

н

о &

се *

о В

о Ч

св т

о &

К

Ч м К

1970

1975

1980

1985

Рис. 2. Стохастическая компонента (логарифмированные данные за вычетом

сезонной компоненты)

Полученные значения характеризуют, на сколько порядков наблюдаемая величина отличается от ожидаемой (средней за все годы наблюдений):

Х1к = ехр(У1к )т(тк), т(тк) = ехр(з(тк)).

Эмпирические авто- и кросскорреляционные функции совместного распределения стохастической компоненты существенно отклоняются от нуля и имеют вид затухающих экспонент. Это указывает на наличие тренда, который может быть описан векторной авторегрессионной моделью ВАР(1) [4]:

Уьъ = Ау1к_г + е1к

(3)

где у1к е К2 4, к > 1, А

вектор значений стохастической компоненты в момент времени — матрица параметров модели, е^ е К2 — вектор остаточ-

2x2

ных значений в момент времени (случайная компонента). Численные оценки значений параметров модели получены методом Юла-Уолкера [5]:

А =

ОВ КП

ОВ КП [ 0,69 0, 06 А ^0,06 0, 79) ■

Авто- и кросскорреляционные функции совместного распределения остаточных членов (случайной компоненты) не отклоняются существенно (на уровне значимости 0, 05) от нуля (рис. 3).

Частное распределение значений случайной компоненты ОВ хорошо описывается нормальным законом, особенно для положительных значений. Статистические тесты Андерсона—Дарлинга, Крамера—фон—Мизеса и Шапиро—Франсиа требуют отвергнуть гипотезу об отклонении эмпирического распределения от нормального даже на уровне значимости 0,1. Наблюдаемое расхождение для отрицательных значений обусловлено, по-видимому, дискретностью исходных данных.

Частное распределение значений случайной компоненты КП является асимметричным и, в отличие от распределения значений случайной компоненты ОВ, не может быть описано нормальным законом. Как видно, распределения случайной компоненты ОВ и КП существенно различаются по своим статистическим свойствам. Это подтверждает мнения некоторых исследователей о том, что индекс Нестерова в том виде, как он описан в ГОСТ, не может считаться вполне адекватной характеристикой пожароопасности [6,7].

Число пожаров

Число пожаров и инд. Нестерова

0,00 0,05 0,10 0,15

(я) АКФ количества очагов возгорания

0,00 0,05 0,10 0,15

(Ь) ККФ количества очагов возгорания — значений комплексного показателя

'-' Н м '-'

1-1-1 4—1 '— -

-0,15 -0,10 -0,05 0,00

(с) ККФ значений комплексного показателя — количества очагов возгорания

0,00

0,05

0,10

0,15

(^ АКФ значений комплексного показателя

Рис. 3. Эмпирические авто- и кросскорреляционные функции совместного распределения остаточных членов модели ВАР(1). Пунктиром обозначен 95% доверительный интервал в предположении о стохастической независимости

Как бы то ни было, распределения обладают выраженной статистической связью: коэффициент корреляции совместного распределения остаточных членов равен 0, 28. В качестве параметрической модели совместное распределения случайной компоненты использовано скошенное нормальное распределение [8]:

в1к ~ ЗЫ2(аПат,а), т.е. р(егк) = 202 (е^,аПат)Ф(ате1к),

где ф2(егк, 0"Пстт) — плотность 2-мерного нормального распределения с ковариационной матрицей аО,ат Е К2х2, а £ К2 — параметр формы. Численные оценки значений параметров модели получены методом максимального правдоподобия:

ОВ КП

_ ОВ / 0,0А о_ОВ {0, 79 0,17Л а = КП у -2, 9), " = КП Ц 17 0, 26 у!.

Таким образом, построена математическая модель МАР(1)С (Мультипликативная АР(1) с Сезонностью) стохастического процесса лесопожарных ситуаций, включающая в себя следующие компоненты:

— мультипликативная сезонная компонента з(т);

— тренд, описываемый векторной авторегрессионной моделью ВАР(1);

— случайная компонента е^, описываемая двумерным скошенным нормальным распределением:

хгк =ехр(у^ + в(пк)),

угк = Ауг-1 + е1к, (4)

е1к ~ ЗИ2(аПат,а).

4. Предсказание и прогноз

В соответствии с пунктом 5.2.6 ГОСТ Р 22.1.09-99, количество лесных пожаров прогнозируют, исходя из:

— степени пожарной опасности в лесу по условиям погоды;

— класса пожарной опасности лесных участков на рассматриваемой территории;

— количества потенциальных источников огня;

— количества пожаров в ретроспективе в аналогичных условиях;

— теоретических законов распределения случайных событий, которыми удовлетворительно описываются (аппроксимируются) некоторые ряды распределения лесных пожаров во времени и по территории.

В стандарте приведена простая методика прогнозирования пожарной опасности с использованием комплексного показателя Нестерова. В соответствии с этой методикой, значение КП определяет так называемый класс пожарной опасности, который характеризует ожидаемый масштаб развития лесопожарной ситуации. Существует общероссийская шкала соответствия значений КП уровню пожарной опасности.2

Описанный способ прогнозирования с использованием значений КП хотя и является очень простым, но не позволяет получать точные и надёжные оценки пожарной опасности. Так, корреляция между значениями КП и ОВ в наиболее пожароопасный период (с 1 мая по 30 сентября в сезон 1990 г.) не превышает 0, 5, а в сезон 1989 г. не превышает 0, 25 (см. табл. 1).

В настоящей работе для прогнозирования будущих значений ОВ использована предложенная модель МАР(1)С. Значения метеорологических показателей и рассчитываемого на их основе комплексного показателя Нестерова могут быть получены на несколько суток вперёд. Математические модели и методы прогнозирования метеорологических показателей выходят за рамки настоящей работы. Следует отметить, что обычно значения комплексного показателя доступны с приемлемой степенью надёжности не более чем на 3 суток [9]. Поэтому прогнозирование развития лесопожарных ситуаций на больший срок является затруднительным.

Разработан функциональный вычислительный алгоритм предсказания значений ОВ при заданных значениях КП на основе предложенной модели МАР(1)С. Прогнозирование уровня пожарной опасности осуществлено методом стохастического моделирования.

Алгоритм позволяет получить не только ожидаемое значение ОВ, но и, что особенно важно, оценить потенциальный риск развития лесопожарной ситуации, т.е. расчитать доверительный интервал для заданного уровня достоверности. В настоящей работе для обозначения верхней границы доверительного интервала используется термин Величина Риска (Value at Risk):

VaRp = inf {¡e|F(x) > p},

где F(x) — функция распределения случайной величины х. Концепция VaR широко используется в финансовой математике для расчёта потенциальных убытков. В случае моделирования лесопожарных ситуаций использование этой концепции также выглядит обоснованным, поскольку экономический урон, наносимый лесными пожарами, также как и затраты на его профилактику и ликвидацию, являются монотонными функциями количества очагов возгорания.

Получены ожидаемые значения уровня пожарной опасности, а также значения VaR для уровней достоверности 90% и 95% в наиболее пожароопасный период с 1 мая по 30 сентября, в сезоны 1989 и 1990 гг. Прогноз производился на 1, 2 и 3 суток вперёд. Значения КП в этот период считались известными заранее. Результаты для п = 2 суток приведены на рис. 4.

2Необходимо отметить, что эта шкала не является универсальной и должна быть скорректирована с учётом местных особенностей, причём не вполне ясны принципы корректировки.

Сезон 1989 г.

Июл

---Действительное значение

- Предсказанное ожидаемое значение

Предсказанный 90% уровень риска Предсказанный 95% уровень риска

Сезон 1990 г.

---Действительное значение

- Предсказанное ожидаемое значение

Предсказанный 90% уровень риска Предсказанный 95% уровень риска

Рис. 4. Прогноз количества пожаров на 2 суток вперёд

Рассчитаны значения ранговой корреляции Спирмена между предсказанными и действительными значениями ОВ. Результаты расчёта приведены в табл. 1. Как видно, прогноз, полученный с использованием модели МАР(1)С, является существенно более близким к действительным значениям, чем прогноз с использованием только показателя Нестерова.

Таблица 1

Численные оценки коэффициента ранговой корреляции Спирмена между действительнымы и предсказанными значениями количества очагов возгорания

Метод прогноза 1989 г. 1990 г.

Инд. Нестерова 0,21 0, 39

МАР(1)С, 1 сутки 0, 76 0, 81

МАР(1)С, 2 суток 0,56 0, 66

МАР(1)С, 3 суток 0,41 0, 56

С целью установления степени адекватности рассчитанных значений УаБ, проведена проверка статистической гипотезы Но о равенстве действительного и теоретического количества превышений значения УаИдо%. Для прогноза на 1 и 2 суток гипотеза Но не может быть отвергнута на уровне значимости 0, 01. Для прогноза на 3 суток и более гипотеза Но должна быть отвергнута на любом разумном уровне значимости. Таким образом, рекомендуемый максимальный горизонт прогноза количества очагов возгорания с помощью модели МАР(1)С составляет 2 суток.

5. Заключение

В настоящей работе проведен анализ статистических свойств совместной выборки значений КП и ОВ за 1969-1990 гг. Построена трёхуровневая математическая модель МАР(1)С, связывающая исследуемые показатели. Разработан вычислительный алгоритм моделирования будущих значений ОВ, который может быть использован для прогноза развития лесопожарной ситуации на 1-2 суток вперёд.

Предложенная модель МАР(1)С использует значения комплексного показателя Нестерова, но также учитывает динамические свойства процесса лесопо-жарных ситуаций. Благодаря этому она позволяет получать существенно более точные и надёжные оценки пожарной опасности.

Необходимо отметить, что существует определённый потенциал для совершенствования методики расчёта комплексного метеорологического показателя. Это тем более важно и актуально, поскольку точность оценок и максимальный горизонт прогноза модели МАР(1)С во многом определяются свойствами используемого комплексного показателя.

Литература

1. ГОСТ Р 22.1.09-99 «Безопасность в чрезвычайных ситуациях. Мониторинг и прогнозирование лесных пожаров. Общие требования».

2. Нестеров В. Г. Горимость леса и методы её определения. — М.: Гослесбумиздат, 1949.

3. MET4 and MET4A Calculation of Dew Point. — http://www.paroscientific. com/dewpoint.htm.

4. Sims C. A. Macroeconomics and Reality // Econometrica. — 1980. — № 48. — С. 1-48.

5. Brookwell P. J., Davis R. A. Time Series and Forecasting Methods. Second edition. — New York: Springer, 1991.

6. Барановский Н. В. Методика прогнозирования лесной пожарной опасности как основа нового государственного стандарта // Пожарная безопасность. — 2007. — № 4. — С. 80-84.

7. Зденева М. Я., Виноградова М. В. Метод среднесрочного прогноза степени пожарной опасности в лесах по метеорологичеким условиям / Методический кабинет Гидрометцентра России. — 2008. — http://method.hydromet. ru/methods/fire/fire.html.

8. Azzalini A., Capitanio A. Statistical Applications of the Multivariate Skew-Normal Distribution // Journal of the Royal Statistical Society, B. — 1999. — № 61. — С. 579-602.

9. Толстых М. А. Прогноз осадков по станциям регионов внетропической зоны Северного полушария на срок до 72 ч на основе полулагранжевой модели с постоянным разрешением ПЛАВ-2005 / Методический кабинет Гидрометцентра России. — 2008. — http://method.hydromet.ru/model/tolstih2/ tolstih2.html.

UDC 519.254.1

Modelling and Forecasting the Evolution of Forest Fire Using Vector Autoregression Processes

A. V. Paramonov

Center of Mathematics Research, MSUT «STANKIN» Vadkovskiy lane, 3A, Moscow, Russia, 117936

The joint sample of daily values of the number of forest fire seats and the complex meteorological index in Irkutsk region, seasons 1969-1988, is investigated. It appears that the evolution of forest fire is well described by a vector autoregression process based model. Computer algorithm of simulation of the future numbers of fire seats is presented. It is shown how the algorithm can be used for reliable forecasting of the evolution of forest fire up to 2 days ahead.