Алгоритм параметризации для модели Ваганова-Шашкина формирования годичных колец древесных растений

Модель Ваганова Шашкина формирования годичных колец древесных растений [1; 2; 3] (далее VS-модель) описывает влияние климатических условий на клеточную структуру годичных колец. Основное предназначение модели заключается в использовании ее как инструмента, позволяющего для произвольной индивидуальной или обобщенной древесно-кольцевой хронологии (ДКХ) в течение набора лет с доступными метеостанционными данными (см. ниже описание входных данных VS-модели): выделить климатически обусловленную компоненту рассматриваемой ДКХ; указать, если качество моделирования рассматриваемой ДКХ VS-моделью удовлетворительно, для произвольных суток данного набора лет, какой из двух факторов – температура воздуха или осад- ки, лимитировал рост древесного растения, соответствующего этой ДКХ. ДКХ является временным рядом значений некоторой числовой характеристики соответствующего этой ДКХ дерева/деревьев [4].

VS-модель может служить тестером качества индивидуальных и обобщенных ДКХ, которые планируется использовать в роли посредников, несущих информацию о климатических данных [2]. ДКХ полагается пригодной для использования в роли такого посредника тогда и только тогда, когда выделенная VS-моделью ее климатически обусловленная компонента достаточно удовлетворительно согласуется с самой ДКХ.

VS-модель является детерминированной динамической имитационной моделью. Ее входные данные состоят из двух блоков.

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

Второй блок представляет собой набор значений параметров VS-модели. Модель имеет 42 параметра. Все параметры, кроме двух, являются вещественными однозначными переменными, причем три параметра – целочисленные однозначные переменные. Два параметра являются векторами равной размерности. Один вектор содержит вещественные однозначные переменные, второй – булевы однозначные переменные. Размерность векторов является параметром VS-модели.

Значение каждого параметра обеспечивает VS-модель информацией либо о моделируемом дереве, либо о его местообитании. VS-модель имеет 10 опций. Набор значений опций задает вариант VS-модели, который будет использован при моделировании. Опции VS-модели определяют совокупность ее модулей, запускаемых при моделировании, точность вычислений при моделировании и т. п.

Выходные данные VS-модели состоят из двух блоков. Данные первого блока являются численными характеристиками, отражающими динамику во времени, временной шаг отслеживания которой не превосходит одних суток, некоторых процессов, протекающих в смоделированных дереве и его местообитании (влажность почвы, величина транспирации дерева, количество камбиальных клеток и т. п.). Содержимое этого блока зависит от используемого набора значений опций VS-модели.

Второй блок выходных данных VS-модели является моделью ДКХ. VS-модель предполагает, что моделируемая индивидуальная ДКХ является временным рядом значений ширины сформировавшегося годичного кольца или количества клеток древесины в таком кольце. VS-модель предполагает, что все индивидуальные ДКХ, использованные при построении моделируемой обобщенной ДКХ, являются такими, как описано выше. Для каждого года климатических входных данных VS-модели второй блок содержит значение смоделированной ДКХ. Если моделировалась обобщенная ДКХ, то это значение не имеет единиц измерения. Если моделировалась ин- дивидуальная ДКХ, то это значение либо не имеет единиц измерения, либо является количеством клеток древесины в годичном кольце (имеющий место вариант определяется используемым набором значений опций VS-модели).

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

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

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

Проблема параметризации для VS-модели – это проблема выбора значений параметров VS-модели, для которых известны лишь диапазоны их значений. Выбор должен осуществляться внутри этих диапазонов таким образом, чтобы полученный набор значений всех параметров VS-модели обладал необходимой внутренней структурой. Как, руководствуясь чем, осуществлять этот выбор? – сущность проблемы параметризации. Проблема параметризации существующей математической модели – это проблема выбора входных данных модели в условиях недостатка информации.

Существуют различные применяемые на практике подходы к решению проблемы параметризации произвольной математической параметрической модели [5]. Представленное ниже решение проблемы параметризации для VS-модели использует имеющийся произвол в выборе значений параметров, для которых отсутствует необходимая для определения их точного значения информация, целенаправленно. Выбираются те значения параметров, которые согласуются с целью исследовате- ля, формализованной в виде достижения оптимума целевой функции, оптимальным, на области определения целевой функции, образом. Предложенное решение проблемы параметризации для VS-модели применимо только в случаях, когда доступна моделируемая сущность, в частности, при использовании VS-модели по ее основному назначению.

Понятие критерия отличия моделируемой ДКХ от ее модели. Введем понятие «критерий отличия моделируемой ДКХ от ее модели», представляющее ценность и вне контекста проблемы параметризации для VS-модели, чтобы изложить предлагаемое в этой статье решение последней. Критерий отличия является однозначной вещественной неотрицательной функцией двух аргументов, обозначим ее DC ( trCr , mCr ), определенной на множестве пар вещественных векторов размерности, равной мощности множества Yrs , с неотрицательными компонентами. Множество Yrs является пересечением двух множеств – набора лет, для которых осуществлялось моделирование, и набора лет, в которые доступны значения моделируемой ДКХ.

Первый аргумент trCr должен быть вектором относящихся к годам Yrs значений моделируемой ДКХ, которые расположены в векторе согласно годам, упорядоченным по возрастанию. Второй аргумент mCr должен быть вектором соответствующих значений ДКХ, смоделированной посредством VS-модели.

Значение критерия отличия характеризует удаленность модели ДКХ, полученной посредством VS-модели, от этой моделируемой ДКХ. При фиксированной моделируемой ДКХ меньшему значению критерия отличия соответствует более удачная модель этой ДКХ. Критерий отличия вводит на множестве моделей некоторой фиксированной ДКХ отношение эквивалентности и упорядочивает классы эквивалентности линейным образом.

Критерий отличия моделируемой обобщенной ДКХ от ее модели должен обладать двумя дополнительными свойствами: свойством симметрии на его области определения и свойством рефлексивности ( trCr = mCr влечет DC ( trCr , mCr ) = 0). Эти два свойства не налагаются на критерий отличия моделируемой индивидуальной ДКХ от ее модели, поскольку семантика значений моделируемой ДКХ (включая единицы измерения) совпадает с семантикой значений модели ДКХ, полученной посредством VS-модели, гарантированно лишь для обобщенных ДКХ.

Необходимые свойства критерия отличия не гарантируют его непрерывности. Многозначные критерии отличия не рассматриваются в этой статье.

Апробация нижеизложенного алгоритма параметризации для VS-модели выполнена со следующими двумя критериями отличия, которые предлагается использовать как критерии отличия по умолчанию. Критерий отличия DC _ITRC индивидуальной ДКХ от ее модели и критерий отличия DC _GTRC обобщенной ДКХ от ее модели определены формулами

DC ITRC ( trCr , mCr ) = 1 - crln( trCr , mCr ) + + 1,25 ( 1 - sync( trCr , mCr ) ) ²,

DCGTRC (trCr, mCr) = DCITRC(trCr, mCr) + +0,4(maxItrCri - mCrI)2, x i e Yrs где crln(trCr, mCr) и sync(trCr, mCr) – коэффициенты корреляции Пирсона и синхронности между векторами trCr и mCr, trCri (mCri) – соответствующая году i-я компонента вектора trCr (mCr).

Выбор таких критериев отличия обусловлен желанием получать расчетные ДКХ как можно более положительно коррелированными с соответствующими им моделируемыми ДКХ; не иметь низких значений коэффициента синхронности между моделируемой ДКХ и ее моделью; не иметь серьезных визуальных отличий между ломаными, представляющими моделируемую обобщенную ДКХ и ее модель. Критерии отличия DC _ITRC и DC _GTRC апробированы на обширном материале, совместимы с нижеизложенным алгоритмом параметризации для VS-модели и отражают взгляд большинства исследователей на понятие близости двух ДКХ.

Алгоритм параметризации. Введем определения и обозначения. Обозначим через p = (p1, …, pn) вектор используемых при моделировании параметров (ИпМП) VS-модели. Вектор p однозначно определяется используемым при моделировании набором значений опций VS-модели. Обозначим номера рассматриваемых n параметров VS- модели, для которых известны конкретные числовые значения, через i1, …, ik; известные конкретные числовые значения этих параметров обозначим ci1 ,..., cik. Диапазон значений i-го параметра VS-модели обозначим [ai ; bi]. Для удобства изложения считаем a. = bij = ci. для j = 1, ..., k.

Назовем пространством оптимизации S подмножество ( n – k ) мерного параллелепипеда

P = { p e R n : a i <  p i <  b i для 1 <  i <  n } , лежащего в R ⁿ , определенное условием: элемент P , рассматриваемый как набор значений ИпМП, принадлежит S тогда и только тогда, когда он обладает необходимой внутренней структурой (см. выше). Набор значений ИпМП, принадлежащий S , назовем допустимым. Считаем S + 0 . Структура S обусловлена используемыми при моделировании взаимосвязями между ИпМП. S напоминает прямоугольный кусок сыра. Типична ситуация: из 100 000 значений непрерывной n мерной случайной величины, равномерно распределенной на P , только одно принадлежит S .

Определим предикатную функцию pFail(p ): P ^ {«TRUE», «FALSE»}, принимающую логическое значение «TRUE» только если p П S . Как правило, для вычисления pFail ( p ) необходимы выходные данные VS-модели, полученные для набора p .

Назовем точностью вдоль i -й оси, 1 <  i <  n , величину

x

e R : VInpOpt1 VInpOpt2 x sup <

(| p i — p ²1 <  x ^ Q ( Output ], Output ₂

где InpOpt1 и InpOpt2 – две совокупности входных данных VS-модели и значений ее опций, отличающиеся лишь значениями pi1 и pi2 i-го ИпМП; Output1 и Output2 – выходные данные VS-модели, соответствующие InpOpt1 и InpOpt2. Предикат Q истинен тогда и только тогда, когда Output1 и Output2 отличаются настолько мало, что этим отличием можно пренебречь и оно несущественно. Обо- значим используемую при моделировании оценку точности вдоль i-й оси через hi .

Введем на S метрику d : S х S ^ [0; + да ), положив

d ( x , y ) = max

xi

-

yi

hi

-

Выбор метрики на S не является тривиальной задачей. На S введена именно такая метрика, поскольку пространство оптимизации имеет по разным осям, вообще говоря, разные физические единицы измерения; константы h, , соответствующие осям с одинаковыми физически-

ми единицами измерения, могут различаться; она является естественным аналогом метрики с_г на R n [6].

Изложим алгоритм параметризации для VS-модели. Предлагаемый алгоритм выбирает набор значений ИпМП из пространства оптимизации S , доставляющий глобальный минимум на S целевой функции F(р ). Алгоритм требует, чтобы F(р ) являлась однозначной вещественной неотрицательной функцией, определенной всюду на S . Алгоритм не требует каких-либо дополнительных свойств F(р ). Целевая функция определена равенством

F(р ) = DC ( trCr , mCr (р )), р О S , где mCr(р ) представляет собой вектор относящихся к годам Yrs значений модели ДКХ (см. выше), полученной посредством VS-модели с набором р значений ИпМП. Предлагаемый алгоритм параметризации для VS-модели применим с любым критерием отличия. Он учитывает лишь те свойства минимизируемого им критерия отличия, которые необходимо имеет любой критерий отличия.

Целевая функция F(р ) и ее область определения S индивидуальны для каждой ситуации, в которой применяется излагаемый алгоритм параметризации. Исследования задачи поиска глобального на S минимума F(р ) показали, что F(р ), как правило, имеет более одного локального минимума. Входные и выходные данные алгоритма, его структура представлены на рисунке.

Метод нулевого порядка многомерного однокритериального поиска глобального минимума, используемый алгоритмом параметризации, является методом мультистарта [7] метода покоординатного спуска [7]. В статье не затрагивается вопрос существования более эффективного метода для решения семейства задач поиска глобального на S минимума F(р ). Используемый метод покоординатного спуска отличается от классического следующим: порядок обхода координатных осей S задается пользователем; на S используется метрика, отличная от евклидовой; используемый метод одномерной оптимизации осуществляет поиск глобального оптимума; отсутствие требования непрерывности F(р ) и сложная структура S требуют модификации классического метода покоординатного спуска.

Настройки алгоритма параметризации представляют собой совокупность шести констант: m , level , maxIter , е , q , maxStart и двух векторов т и wght . Алгоритм использует одни и те же значения его настроек на протяжении всех осуществляемых им действий. Следующие настройки

ВХОДНЫЕ ДАННЫЕ АЛГОРИТМА

ПАРАМЕТРИЗАЦИИ

АЛГОРИТМ ПАРАМЕТРИЗАЦИИ

ВЫХОДНЫЕ ДАННЫЕ АЛГОРИТМА ПАРАМЕТРИЗАЦИИ

Структурная схема алгоритма параметризации для VS-модели алгоритма параметризации предлагается использовать как е"0 настройки по умолчанию: m = 3, level < infpеS F(p), maxIter = 500, ε ≤ (n – k)–1, q = 1, maxStart = 1 000, τi ≤ hi ⋅ min (1, ε) (1 ≤ i ≤ n), wghtj = 0 для j = i1, …, ik (значения по умолчанию остальных n – k положительных компонент вектора wght в статье не приводятся). Нечетность m и истинность неравенства τi ≤ hi ⋅ ε, 1 ≤ i ≤ n, необходимы для корректного функционирования алгоритма.

Метод мультистарта запускает с одними и теми же настройками метод покоординатного спуска из maxStart стартовых точек в S . Наименьший из maxStart найденных локальных минимумов является оценкой глобального минимума и выдается как результат работы метода мультистарта. Стартовые точки генерируются случайным образом и являются различными значениями непрерывной n мерной случайной величины, равномерно распределенной на P .

Веса wght задают порядок обхода координатных осей S методом покоординатного спуска. Каждая итерация последнего заключается в выполнении одномерных минимизаций вдоль координатных осей S с положительными весами; в конце итерации проверяется условие останова. Одномерные минимизации вдоль осей, имеющих максимальный вес, выполняются в первую очередь. Оси с одинаковыми положительными весами обходятся в порядке, индуцированном их номерами, – чем меньше номер, тем раньше осуществляется минимизация.

Минимизация вдоль i -й оси, 1 ≤ i ≤ n , означает поиск глобального минимума функции f ( x ) = F (…, p_{i –1}, x , p_i+1 , …) вещественной переменной x на множестве, лежащем внутри отрезка [ a_i ; b_i ] и задаваемом предикатной функцией xFail ( x ) = pFail (…, p_{i –1}, x , p_i+1 , …). После его осуществления значение x ^*, доставляющее глобальный минимум f ( x ), записывается в i -ю координату текущей точки p в S . При входе в первую итерацию текущая точка p совпадает со стартовой, выбранной методом мультистарта. Координаты p с нулевыми весами не изменяются на протяжении всех осуществляемых методом покоординатного спуска действий.

Условие останова выполняется в случае истинности хотя бы одного из условий: количество итераций достигло порога maxIter ; F ( p ) <  level , где p является текущей точкой в конце итерации (после выполнения всех одномерных минимизаций данной итерации); d ( p^N , p^N–i ) <  ε для всех i = 1, …, q ( N – номер текущей итерации, p ^j – текущая точка в S в конце j -й итерации, q и ε – положительные константы); d ( p^N , p^{N –1}) = 0.

Одномерные минимизации вдоль координатных осей S выполняет метод, являющийся модификацией итерационного метода деления отрезка на несколько равных частей [7]. На каждой итерации метод разбивает текущий отрезок на m равных частей и либо делает текущим отрезком одну из них, либо завершает работу, сообщая пользователю о необходимости изменить значение m . Итерации прекращаются, когда значение минимизируемой функции на текущем отрезке становится меньше уровня level или когда длина текущего отрезка становится меньше величины д, заданной до запуска метода.

Все одномерные минимизации вдоль i -й оси, осуществляемые в ходе работы метода покоординатного спуска, выполняются с δ = τ _i (1 ≤ i ≤ n ).

Реализация и апробация алгоритма параметризации. Алгоритм параметризации для VS-модели реализован на языке программирования C++ с привлечением средств стандартной библиотеки C++. Его текстовый код удовлетворяет стандарту С++03 (ISO/IEC 14882:2003). Исполняемый код реализации для конкретной операционной системы может быть получен без изменений текстового кода реализации и оформлен в виде dll-библиотеки (если целевой платформой является ОС семейства Windows) или mex-файла (если целевой платформой является ОС семейства Unix).

Основными факторами, влияющими на продолжительность параметризации некоторой ДКХ посредством рассматриваемой реализации алгоритма параметризации, являются следующие: используемый при моделировании набор значений опций VS-модели; мощность набора лет, на которых выполняется моделирование; используемые при моделировании взаимосвязи между ИпМП; значение величины n – k и длины диапазонов значений ИпМП, для которых неизвестны конкретные числовые значения; используемый при параметризации критерий отличия; используемые настройки алгоритма параметризации; аппаратное и программное обеспечение устройства, выполняющего код реализации алгоритма параметризации. Варьируя любым из этих факторов, можно добиться многократного изменения продолжительности параметризации.

На ноутбуке Asus F3JP с двухъядерным процессором Т5300 (максимальная тактовая частота каждого ядра равна 1,73 GHz) параметризация ДКХ посредством рассматриваемой реализации алгоритма параметризации обычно длится от нескольких часов до нескольких суток.

Предложенный алгоритм параметризации для VS-модели был использован на более чем 700 различных парах ( S , F ( p )). Следствием осуществления этих параметризаций явилось то, что ряд утверждений о VS-модели, как новых, так и ранее сформулированных, подтверждены впервые вычислительным экспериментом; впервые сформулирована и подтверждена вычислительным экспериментом гипотеза о взаимосвязи моделей индивидуальных и обобщенной ДКХ, полученных посредством VS-модели; VS-модель впервые в полной мере использована как тестер качества индивидуальных и обобщенных ДКХ, которые планируется использовать в роли посредников, несущих климатическую информацию; установлены несколько свойств рассматриваемого алгоритма параметризации. Успешная апробация алгоритма параметризации и его реализация является следствием этих результатов.

Результаты данной работы востребованы преимущественно в дендроклиматологии и дендрохронологии, математическом моделировании процессов, протекающих в древесных растениях. Созданный алгоритм параметризации и его реализация существенно расширяют круг ситуаций, в которых использование VS-модели уже осуществимо на практике; позволяют тестировать и анализировать VS-модель на качественно новом уровне; позволяют решать на практике обратные задачи восстановления некоторых условий произрастания древесного растения по имеющейся его ДКХ.