Е:
Ч Б. Нармурадов
~ и номов на эле-: - Ouiо известно определяется
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ 2005 год, том 17, номер 9, стр. 43-52
ГЛОБАЛЬНЫЕ ПОКАЗАТЕЛИ ЧУВСТВИТЕЛЬНОСТИ
ДЛЯ ИЗУЧЕНИЯ НЕЛИНЕЙНЫХ МАТЕМАТИЧЕСКИХ МОДЕЛЕЙ
И М. Соболь
I
■ : раз больше, - " чов, должно гоцессе рас-ерии расче-яце, задан-
:тода (СМ) фально-се-- з алгебраи-*о< = два раза.
Таблица 3
- Новосибирск:
(Spectral)
для Ново-
LPl4.
Институт математического моделирования РАН, г.Москва, 125047, Миусская пл., д.4А, e-mail: hq@imamod.ru
Работа выполнена при поддержке гранта РФФИ № 04-01-00367
Статья представляет собой обзор теории глобальных показателей чувствительности, которые позволяют численно исследовать структуру нелинейной функции, заданной аналитически или как «черный ящик».
GLOBAL SENSITIVITY INDICES FOR THE INVESTIGATION OF NONLINEAR MATHEMATICAL MODELS
I.M. Sobol'
Institute for mathematical modeling RAS, Moscow, e-mail: hq@imamod.ru
The paper is a review of the theory of global sensitivity indices, that allow to investigate numerically the structure of a nonlinear function defined analytically or as a "black box".
1. Введение
Глобальные показатели чувствительности для изучения нелинейных функций были введены в 1990-м году [1]. Перевод этой статьи на английский язык был опубликован в 1993-м году. С тех пор появилось несколько десятков статей на английском языке, посвященных теории и приложениям этих показателей, и всего одна статья на русском языке [2].
Цель настоящей работы - изложить современное состояние этой теории. С приложениями (в технике, экологии, медицине, финансах и др.) можно познакомиться по трудам конференций SAMO (Sensitivity Analysis of Model Output) [3-5] и сборнику [6]. Предыстория вопроса имеется в [2] и [6].
Предположим, что исследуемая математическая модель описывается функцией
где х=(х1,...гх„) - входные данные, а и - выход (скаляр). Традиционный анализ чувствительности, который можно назвать локальным, рассматривает какое-нибудь решение и =j[x ) и оценивает чувствительность и по отношению к х,; она зависит от частной производной
Глобальный анализ чувствительности не предполагает выделения индивидуальных решений. Функция fix) рассматривается во всей области определения: изучается влияние отдельных переменных х, и их групп, выделяются существенные и несущественные переменные, выясняется структура^) и возможности аппроксимации^*) более простыми функциями.
2. Разложение ANOVA
Рассмотрим функцию/(*). где х=(хь.. ,гх„), определенную и интегрируемую с квадратом в единичном n-мерном кубе К",
К" = {О <<1,...,0<*„ <1}.
В дальнейшем символ I без указания пределов интегрирования означает, что интегрирование осуществляется от 0 до 1 по всем переменным; dx=dx1...dx„. Определение. Формула
/w = /o+Z X fK.....(D
где внутренняя сумма берется по всевозможным i\,...,is, удовлетворяющим неравенствам 1 < ¡1 <... < is < п, называется ANOVAразложением функцииJ[x), если
/о = \f{x)dx (2)
и
I
\fhr..,hdxip = 0 при1<р<5. (3)
о
Развернутая запись формулы (1)
/М = /о + Z fi (*< ) + S fij (Xi'Xj ) ■+ ••• ■+ /l,2,...,n (*1 .*2, )•
i i<j
В [1] ив предшествующих работах автора использовался термин «разложение на разноразмерные слагаемые». Термин ANOVA происходит от слов Analysis of Variances [7]. Разложение (1) можно построить с помощью любой полной ортонормированной системы функций <Ро(х),ч>|(х),...,ч>*(х),..., содержащей функцию (p0(x)sl: достаточно построить разложение/(*) в ряд Фурье по системе {<р*(*)} и сгруппировать слагаемые различных размерностей.
Нетрудно доказать, что условия (2) и (3) однозначно определяют все члены разложения (1). Во-первых, проинтегрируем (1) по всем переменным, кроме х}:
\f{x)Y\dxp =fn+Mx,). p*i
Отсюда находим одномерные слагаемые /¡(х,). Во-вторых, проинтегрируем (1) по всем переменным, кроме х, и х/.
J/M П dxP = Л + f,(*г) + /;(*;) + fij(*,•>*;)•
P*i,j
Отсюда находим двумерные слагаемые ¡¡¡(х^х,). И так далее. Последнее слагаемое
/12.....„(*|,х2,...,*„) определяется из соотношения (1).
Нетрудно также доказать, что все слагаемые в (1) ортогональны (в силу (3)): если (ii,...,i,)*Ui,--J'i), то ,...,;,/;„..=
3. Дисперсии
Величины
ИМ.Соболь
Глобальные показатели чувствительности для изучения нелинейных математических моделей
45
-чую с квадратом в
интегрирование
щ*м
(1)
неравенствам
(2)
(3)
ние на разно-7]. Разложе-»1 функций жение /(х) в
разложения
> всем перемен-
слагаемое
называются дисперсиями, а
/)= \/2(х)сЬ-/02
- полной дисперсией. Происхождение этих названий очевидно: если бы точка х была случайной точкой, равномерно распределенной в К", то значения Дх) и ^ (х^^х^) оказались бы случайными величинами, а/) и /),• ^ - их дисперсиями.
Возводя (1) в квадрат и принимая во внимание ортогональность слагаемых, получим равенство
Д = 1 I А.....
(4)
Величина Б характеризует изменение/(дс) в К", а величина определяет вклад /......, в это
изменение. Для кусочно-непрерывных функций равенство Ц ^ =0 означает, что Д,...^ - 0.
Стоит отметить, что для более или менее сложных функций Дх) вычислить многомерную компоненту / ( (л^ ) практически невозможно, а постоянные Д ^ , как будет показано ниже, достаточно просто вычисляются по значениям Дх) в специальным образом выбранных точках.
4. Глобальные показатели чувствительности
Определение. Глобальными показателями чувствительности называются отноше-ния дисперсий
5, , =3,- , т.
Очевидно, все ^ нормирована:
неотрицательны, не превосходят единицы и, в силу (4), сумма их
X I .....,,=••
(5)
$ = 1 I, <...<!,
В приложениях чаще всего используют одномерные показатели Например, с их помощью ранжируют переменные х{. чем больше тем более влиятельна переменная х,. Впрочем, для полного учета влияния х, следует использовать полные глобальные показатели Я!0', которые определены в следующем разделе. Здесь мы отметим один вывод, важный для приложений: Для того, чтобы Дх) представляла собой сумму одномерных слагаемых, необходимо и доста-
п
точно, чтобы =1.
/=1
5. Глобальные показатели чувствительности для множества переменных
Рассмотрим произвольную группу переменных хк] ,...,хк , где \<к\<...<кт<п, 1 <т<п-\.
Условимся обозначать их одной буквой у = (хк> ,...,хк ); и пусть г - совокупность всех п-т остальных переменных. Таким образом, х=(у,2).
Через М обозначим совокупность индексов (ки...,кт). Для множества >> введем два типа глобальных показателей чувствительности:
где сумма берется по всем группам ¿1,...,/„ в которых все ¡реМ, и
.....
где сумма берется по всем группам таким, что хотя бы один индекс Очевидно,
Показатели предложил ввести А.Салтелли [2].
В случае множества;)', состоящего из одной переменной у=(х,), значение 5^=5,, а 8'°', которое уместно обозначить 5'°', есть сумма всех показателей 5,- , , один из индексов которых равен г.
Пример. ПустьЕслиу=(х,), то 6^=5,, а Б1"' = + 512 + + 5|23. Еслиу=(хьх2), то ¿^=51+52+512. В этом случае г=(х}) и
Другой подход к определению и связан с преобразованием формулы (1). Обозначим через £•](>>) сумму всех слагаемых, зависящих только от переменных множества у, и пусть g2(z) - сумма всех слагаемых, зависящих только от переменных множества 2. Тогда вместо (1) получим равенство
о <sy<s'y01 <1.
Наиболее информативны крайние значения:
A) Sу = S'°' = 1 означает, что flx) зависит только от_у;
Б) S у = S'"' = 0 означает, что flx) зависит только от z. Из определения и равенства (5) следует, что
S'y01 = 1-SZ, S'z0'=\-Sy.
/to = /0 + g, (у) + g2(z) + gi2(y,z). Справедливы соотношения, аналогичные (3):
(6)
Поэтому, возведя (6) в квадрат и проинтегрировав по dx=dydz, получим, что
D = Dy + Dz + D]
уг>
где
Dy = \g]dy, Dz = \g\dz, Dyz=\ g*2dx.
Обозначим
D'y01 =Dy+Dyz, D'°'=D2+Dj
Тогда нетрудно проверить, что
\tOt
У
D '
и аналогично выражаются S2 и S'"'.
=\f. Очевидно,
ж S=S„ a S'y°', ко-
кз индексов кото-
-5,
123 ■
чулы (1). Обозна-
: жества у, и пусть Тогда вместо (1)
(6)
Большое достоинство излагаемой теории заключается в том, что не надо суммировать ,2 , чтобы найти Бу и 5'°'. Дисперсии Д, и можно вычислять непосредственно по значениям функции^*).
6. Интегральные представления для Оу и
Рассмотрим две точки х и х' из К" и будем считать, что х=(у,2), а х'=(у',гг). Теорема 1.
йу = ¡Дх)Ду,г^1'-/102. (7)
Доказательство. Интеграл в (7) можно преобразовать:
Iйу \/(у, ¡Ду,= ¡Оу [ \Ду,.
Подставив (6) во внутренний интеграл и проинтегрировав по ¿г, мы получим /о +g^(y)■ Возведя это выражение в квадрат и проинтегрировав по с1у получим Оу + /02. о
Теорема 2.
О'? = X-\[f{x)-f{y\z)-fdxdy'.
Доказательство. Выражение, стоящее в (8) справа, равняется
(8)
Используя теорему 1, получим
¡/2(х)<1х-(В2+/2) = В-В2. о
В ходе доказательства теоремы 2, мы использовали теорему 1 применительно к множест-
ву z:
D, = \f(x)f(y\z)dxdy'~f2
(9)
и затем вычислили В'°' =D-DZ. Таким образом, формула (9) может заменить формулу (8) для нахождения В'°'. Вполне аналогично, записав формулу (8) применительно к D'°', можно вычислить Dy = D-D'°'.
Алгоритмы Монте-Карло, соответствующие этим представлениям, сравнивались в статье [8]. Оказалось, что алгоритмы, соответствующие представлениям (7) и (8), обладают меньшими дисперсиями (или, что то же, большей точностью при одинаковом количестве испытаний),
7. Алгоритм Монте-Карло
Испытание номер к, где k=\,2,...,N, начинается с выбора 2п стандартных случайных чисел. Из этих чисел формируются две случайные величины тц и г\'к , равномерно распределенные в 1С, и две случайные величины Q и Çk , равномерно распределенные в Затем математическая модельДх) вычисляется в трех точках:Д%,^, /(л*,Ci) и после чего вычисляются четыре оценки:
Ф*=/(Л Ф*, V* =Ф*/(Л*>^) и X* = -/Oli.Ç*)]2-
После N независимых испытаний, получим
1 к
— Е<Р* ——>/о,
14 к=I
" *=I
Вполне аналогично строится алгоритм квази-Монте-Карло. Для реализации испытания номер к, выбирается точка Qk = (д,* ,...,<?*„) последовательности, равномерно распределенной в К2". Тогда
Л* =(^+и+|.->?2(|)-
Замечание. Если значение /02 велико (по сравнению с Эу), то при использовании
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.