Планирование экспериментов для дискриминации регрессионных моделей тема диссертации и автореферата по ВАК РФ 01.01.07, кандидат наук Гученко Роман Александрович

• Гученко Роман Александрович
• кандидат наук
• 2018

Оглавление диссертации кандидат наук Гученко Роман Александрович

Введение

Глава 1. Постановка задачи планирования дискриминационных

экспериментов и обзор известных результатов

1.1 Задачи анализа и планирования в контексте восстановления регрессионных зависимостей

1.2 Планы эксперимента

1.3 Планирование экспериментов для дискриминации регрессионных моделей

1.4 Связь задачи дискриминации и задачи по оцениванию параметров

1.5 Критерий Т-оптимальности и его статистическая интерпретация

1.6 Теорема эквивалентности для критерия Т-оптимальности

1.7 Связь критерия Т-оптимальности с теорией аппроксимации

1.8 Алгоритм для численного построения Т-оптимальных планов

1.9 Различные обобщения критерия Т-оптимальности

Глава 2. Т-оптимальные планы для дискриминации

дробно-рациональных и полиномиальных моделей

2.1 Простейшая дробь и полином

2.2 Другие дробно-рациональные модели

2.3 ЕМАХ-модель и квадратичная модель

2.3.1 Локально-оптимальный план

2.3.2 Робастный план

Глава 3. Построение байесовских Тр-оптимальных планов

3.1 Различные аналоги критерия Т-оптимальности

3.2 Сведение байесовских планов к локально оптимальным

3.3 Теорема эквивалентности для Тр-оптимальных планов

3.4 Численные алгоритмы для построения Тр-оптимальных планов

3.4.1 Квадратичное программирование

3.4.2 Градиентный метод

3.5 Примеры

3.5.1 Байесовские Тр-оптимальные планы для дискриминации экспоненциальных моделей

3.5.2 Байесовские Тр-оптимальные планы для дискриминации моделей зависимости эффективности препарата от дозы

Глава 4. Построение байесовских ^Ьр-оптимальных планов

4.1 Обобщение критерия Т-оптимальности на случай произвольно распределенных ошибок наблюдения

4.2 Байесовские ^Ьр-оптимальные планы и их сведение к локально оптимальным

4.3 Теорема эквивалентности для ^Ьр-оптимальных планов

4.4 Численные алгоритмы для построения ^Ьр-оптимальных планов

4.4.1 Квадратичное программирование

4.4.2 Градиентный метод

4.5 Примеры

4.5.1 ^Ь-оптимальные планы для дискриминации ЕМАХ-модели и модели Михаэлиса Ментен

4.5.2 Байесовские ^Ьр-оптимальные планы для дискриминации экспоненциальных моделей

4.5.3 Байесовские ^Ьр-оптимальные планы для дискриминации моделей зависимости эффективности препарата от дозы

Глава 5. Построение полу-параметрических оптимальных планов

5.1 Полу-параметрические оптимальные планы

5.2 Теоремы эквивалентности для полу-параметрических оптимальных планов

5.3 Связь полу-параметрических критериев с критерием Т-оптимальности

5.4 Примеры

5.4.1 Полу-параметрические оптимальные планы для дискриминации ЕМАХ-модели и модели Михаэлиса Ментен

5.4.2 Полу-параметрические оптимальные планы для дискриминации экспоненциальной модели и модели Михаэлиса Ментен

5.4.3 Сравнение подходов

Заключение

Список литературы

Список рисунков

Список таблиц

Введение диссертации (часть автореферата) на тему «Планирование экспериментов для дискриминации регрессионных моделей»

Введение

Актуальность темы. Математическая теория планирования эксперимента предоставляет методологию для оптимального выбора условий проведения наблюдений при практических исследованиях, которая позволяет получать результаты, касающиеся числовых характеристик регрессионных моделей для изучаемых явлений, с заданной степенью статистической погрешности при наименьших затратах. Название этому направлению в науке дала книга «The Design of Experiments» знаменитого английского статистика Р. Фишера. Базовый аппарат теории был разработан во второй половине прошлого века в трудах Дж. Элфвинга, Дж. Кифера, Дж. Вольфовица, К. Pao, В. Стаддена, В.В. Федорова, Г. Уинни и других авторов. Целью планирования является нахождение (приближенного) плана эксперимента — дискретной вероятностной меры, заданной на множестве всех возможных условий проведения измерений, оптимальной с точки зрения некоторого заранее заданного критерия.

Большинство работ по оптимальному планированию экспериментов посвящено планам для оценивания неизвестных параметров регрессионных моделей. В этом случае предполагается, что функция регрессии задана с точностью до неизвестных параметров, которые необходимо оценить с некоторой точки зрения оптимально посредством выбора условий для проведения измерений. Примерами оптимальных планов для оценивания параметров являются D-оптпмальные планы, минимизирующие объем доверительного эллипсоида для оценок неизвестных параметров в предположении о независимости, гомос-кедастичности и нормальной распределенности случайных ошибок измерения, или ^-оптимальные планы, минимизирующие в том же эллипсоиде длину максимальной оси. Однако при проведении прикладных исследований в различных областях знаний нередко возникает ситуация, когда вид регрессионной модели с точностью до параметров не известен a priori до проведения эксперимента, но тем не менее у экспертов есть несколько гипотез о возможном виде модели. В этом случае проводят эксперимент специального вида — дискриминационный эксперимент, планируемый таким образом, чтобы по его результатам можно было оптимально относительно некоторого критерия проверить гипотезу об истинном виде исследуемой регрессионной модели. Самыми популярными критериями, используемыми для решения задач дискриминации, являются критерий

Т-оптимальности, введенный в работах Аткинсона и Федорова, и различные его обобщения.

Критерий Т-оптимальности для дискриминации двух конкурирующих моделей предполагает наличие априорно заданного фиксированного значения параметров для одной из моделей. Это его свойство называется локальностью. Современная статистическая практика требует нахождения планов, устойчивых относительно некорректного выбора фиксированного вектора параметров, таких как байесовские оптимальные планы, где вместо точечного фиксированного значения параметров берется некоторое распределение. Явный вид опорных точек и весов даже для локального Т-оптимального плана получен в литературе только в случае дискриминации полиномиальных моделей, отличающихся на один или на два порядка. Отыскание же байесовских планов в явном виде обычно не представляется возможным. Для численного нахождения дискриминационных планов обычно используют различные варианты алгоритма, предложенного Аткинсоном и Федоровым. В случае байесовских планов применение этого алгоритма становится проблематичным. Данная диссертационная работа посвящена разработке эффективных численных алгоритмов для нахождения дискриминационных планов, а также нахождению Т-оптимальных планов в явном виде для некоторых специальных пар регрессионных моделей.

Современное состояние исследований. Первые работы по планированию дискриминационных экспериментов появились в начале 70-х годов прошлого века и были связаны с дискриминацией вложенных моделей, когда одна модель является частным случаем другой при определенных значениях параметров. Так в работе [1] для дискриминации между двумя вложенными полиномиальными моделями в предположении о независимости, нормальности и гомоскедастичности ошибок было предложено искать план, доставляющий минимум объему доверительного эллипсоида для тех параметров более общей модели, которые не входят в менее общую. В рамках этого подхода также написаны следующие работы: [2], где планы для оценивания старших коэффициентов полиномиальной модели были найдены в явном виде; [3], где определяется план, оптимальный для выбора количества слагаемых в модели Фурье; [4], где рассматриваются одновременно задача оценивания степени полинома и задача оценивания его коэффициентов; [5], где рассматривается задача дискриминации для нескольких вложенных нелинейных по параметрам моделей.

Другой критерий оптимальности, Т-критерий, применимый для решения задачи о дискриминации двух конкурирующих регрессионных моделей в случае независимых нормально распределенных гомоскедастичных ошибок, не предполагающий вложенности этих моделей, был введен в [6]. В этой же работе был предложен численный алгоритм для нахождения Т-оптимальных планов. Критерий Т-оптимальности тесно связан с задачей наилучшей чебышёвской аппроксимации. Наиболее полно соответствующие результаты изложены в [7], а в [ ] они использованы для получения в явном виде Т-оптимальных планов для дискриминации двух полиномиальных моделей, отличающихся на два порядка.

Базовый критерий Т-оптимальности имеет ряд существенных ограничений. Одно из ограничений заключается в требовании о нормальности и гомос-кедастичности ошибок наблюдения. На случай дискриминации двух моделей при нормальных гетероскедастичных ошибках Т-критерий был обобщен в [ ]. В [ ] был предложен критерий ^Т-оптимальности, основанный на расстояниях Кульбака Лейбл ера и связанный с тестом отношения правдоподобия, не требующий нормальности и пригодный для дискриминации двух произвольных конкурирующих моделей для плотностей ошибок. Критерии из работ [6] и [9] являются частными случаями ^Т-критерия. В [ ] был введен полу-парамет-рический критерий, в котором функции регрессии для всех конкурирующих моделей и закон распределения ошибок для одной из них считаются известными, а закон распределения для оставшейся модели получается как решение специальной задачи вариационного исчисления. Другое ограничение Т-критерия состоит в количестве сравниваемых моделей. Авторы оригинальной работы предложили в [12] вариант его обобщения на случай дискриминации произвольного количества конкурирующих моделей. Симметричная версия Т-критерия для дискриминации многих моделей, Тр-критерий, а также численный алгоритм для нахождения оптимальных планов, основанный на многомерной чебышёвской аппроксимации, были введены в [ ]. Еще одно ограничение Т-критерия — зависимость от априорных значений для параметров одной из моделей, может быть компенсировано с помощью байесовского подхода, обсуждавшегося еще в [6].

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

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

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

2. Предложить метод построения байесовских Тр-оптимальных планов для дискриминации нескольких моделей, преодолевающий недостатки метода Аткинсона и Федорова.

3. Обобщить этот метод на случай байесовских ^Ьр-оптимальиых планов.

4. Исследовать полу-параметрические критерии оптимальности и их связь с классическими критериями, упростить их численное нахождение.

Научная новизна. Все результаты, представленные в диссертационной работе, являются новыми.

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

Методология и методы исследования. В работе применяются методы теории аппроксимации, функционального анализа, математической статистики, вариационного исчисления и общие методы теории планирования эксперимента. Численные примеры выполнены с использованием статистического пакета Р.

Положения, выносимые на защиту.

1. В явном виде получены Т-оптимальные планы для дискриминации полиномиальных моделей и аналогичных моделей, содержащих дополнительное дробно-рациональное слагаемое.

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

сто используемого в литературе алгоритма с новым, выявившее значительное преимущество последнего.

3. Сформулирован байесовский критерий ^Ьр-оптимальности и теорема эквивалентности для него. Результаты, касающиеся байесовскихТр-оптималь-ных планов, обобщены на случай байесовских ^Ьр-оптимальных планов.

4. Предложен упрощенный метод численного нахождения полу-параметриче-ских оптимальных планов. Доказаны две теоремы, связывающие полу-пара-метрические критерии с критерием Т-оптимальности.

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

Основные результаты диссертационной работы докладывались на семинаре кафедры статистического моделирования математико-механического факультета СПбГУ. Исследования по теме диссертации были частично поддержаны грантами СПбГУ (проект 6.38.435.2015) и РФФИ (проект 17-01-00161).

Публикации на тему диссертации. Основные результаты диссертации изложены в трех работах [14], [15], [16], которые опубликованы в журналах, индексируемых системой Зсорин.

Вклад диссертанта в совместные работы. Вторая глава основана на работе [14], в которой все результаты и основной текст принадлежат диссертанту. Постановка задачи и текст вводной секции принадлежит соавтору. Третья глава основана на работе [15], в которой метод сведения байесовских планов к локально оптимальным, алгоритм 3.2, секция про метод оптимизации по весам, основанный на квадратичном программировании, теорема 3.3 и все численные результаты принадлежат диссертанту. Постановка задачи, итоговый английский текст и остальные теоремы принадлежат соавторам. Четвертая глава основана на работе [16], в которой адаптация алгоритма из предыдущей работы на случай ^Ьр-оптимальных планов и все численные результаты принадлежат диссертанту. Постановка задачи, итоговый английский текст и все остальные результаты принадлежат соавторам. Пятая глава основана на работе, выполненной во время визита диссертанта в Рурский университет Бохума

совместно с X. Детте, В.Б. Медасом и В.К. Вонгом. Теоремы 15, 17 и 18, лемма 4, а также все численные результаты принадлежат диссертанту.

Структура работы. Диссертация состоит из введения, пяти глав и заключения.

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

В главе 1 сформулирована задача планирования эксперимента для дискриминации регрессионных моделей, описано ее место в математической статистике и в теории планирования эксперимента, а также дан обзор имеющихся в литературе результатов по данной теме. В § 1.1 приведены постановки задач анализа эмпирических данных и оптимального планирования экспериментов в контексте восстановления регрессионных зависимостей, описана разница между этими двумя постановками. В § 1.2 приведены понятия точного и приближенного планов эксперимента, описано соотношение между ними. Задача планирования эксперимента для дискриминации конкурирующих регрессионных моделей в общем виде сформулирована в § 1.3. В § 1.4 описано то, как связаны планы эксперимента для оценивания неизвестных параметров одной регрессионной модели и планы эксперимента для дискриминации нескольких регрессионных моделей в случае вложенных конкурирующих регрессионных моделей. Наконец, в § приведен критерий Т-оптимальности для планов дискриминации (предложенный в работах [6; 12]), исследованию различных обобщений которого посвящена содержательная часть работы. В том же параграфе дана статистическая интерпретация Т-критерия для случая точных планов эксперимента. В § для Т-критерия приведена соответствующая теорема эквивалентности, дающая необходимые и достаточные условия для оптимальности плана. В сущности, все последующие теоремы эквивалентности для обобщений Т-критерия, которые будут рассматриваться в данной работе, имеют схожую структуру и доказываются аналогично теореме из этого параграфа. В § 1.7 приведены теоремы, связывающие задачу нахождения Т-оптимального плана с задачей наилучшей чебышёвской аппроксимации. Используя эти теоремы в главе мы получим в явном виде Т-оптимальные планы для дискриминации полиномиальных моделей и моделей, отличающихся от полиномиальных добавлением дробно-рационального слагаемого. В § 1.8 описан самый распространенный

численный алгоритм для нахождения Т-оитимальных планов. Альтернативный этому новый численный алгоритм для нахождения байесовскихТ-оптимальных планов описан в главе 3. В § 1.9 дан обзор критериев, обобщающих критерий Т-оптимальности.

В главе в явном виде получены Т-оптимальные планы для дискриминации полиномиальных моделей и моделей, отличающихся от полиномиальных добавлением дробно-рационального слагаемого. Нахождение аналитических представлений для точек и весов Т-оптимальных планов является трудной минимаксной задачей, поэтому для получения Т-оптимальных планов обычно используют численные алгоритмы. До недавнего времени аналитическое решение было известно только для случая дискриминации двух полиномиальных моделей, отличающихся на один порядок. Известно (смотри, например, [7]), что задача нахождения Т-оптимального плана тесно связана с задачей наилучшей чебышёвской аппроксимации в том смысле, что опорные точки оптимального плана являются точками максимума функции, представляющей из себя квадрат разности между регрессионной моделью из нулевой гипотезы и ее наилучшей чебышёвской аппроксимации альтернативной моделью. На основе этого факта в работе [ ] были аналитически получены Т-оптимальные планы для дискриминации полиномиальных моделей, отличающихся на два порядка. В § 2.1 с использованием известных фактов из теории аппроксимации доказано, что опорные точки Т-оптимального плана для дискриминации простейшей дробно-рациональной модели и полиномиальной модели являются корнями полинома специального вида. Корни этого полинома, а значит и опорные точки соответствующего плана, могут быть выписаны для небольших степеней полиномиальной модели в явном виде. В § 2.2 аналогичный результат доказан для двух других дробно-рациональных моделей. В § 2.3 при помощи теоремы из § аналитически найден Т-оптимальный план для дискриминации ЕМАХ-модели и квадратичной модели, а также численно построен план, оптимальный с точки зрения байесовского стандартизированного критерия Т-оптимальности. Изложение в главе 2 опирается на материал из работы [14].

В главе предложен метод построения байесовских Тр-оптимальпых планов. У стандартного критерия Т-оптимальности есть ряд существенных недостатков: он предназначен для сравнения только двух моделей и он является локальным, то есть зависит от фиксированного априорного значения парамет-

ров одной из моделей. В § 3.1 перечислены различные обобщения критерия Т-оптимальности, призванные устранить эти недостатки. Для дискриминации нескольких моделей в работе [ ] был предложен критерий Тр-оптимальности, представляющий из себя взвешенную сумму из попарных квадратов разностей между конкурирующими моделями, где параметры одной из моделей в паре фиксированы, а параметры другой подбираются по методу наименьших квадратов. Стандартным же способом компенсации локальности любого критерия является введение его байесовской версии, где вместо фиксированного значения параметров используется фиксированное априорное распределение. В § 3.2 вводится байесовский критерий Тр-оптимальности и показывается, что он сводится к обычному критерию Тр-оптимальности с большим числом попарных сравнений в случае, когда априорные распределения для параметров дискретны. В § сформулирована теорема эквивалентности для Тр-критерия. В § предложен новый двухэтапный численный метод для нахождения Тр-оптимальных планов, суть которого состоит в чередовании обновления носителя промежуточного плана и оптимизации по его весам. Далее доказывается сходимость этого алгоритма, а также рассматриваются две процедуры для выполнения эффективной оптимизации по весам промежуточного плана. В § 3.5 приведены примеры полученных численно байесовских Тр-оптимальных планов с дискретными априорными распределениями. Изложение в главе 3 опирается на материал из работы [15]. Алгоритмы из § 3.4 реализованы в качестве пакета [17] для языка Р.

В главе все результаты, полученные для байесовских Тр-оптимальных планов, обобщены на случай байесовских ^Тр-оптимальных планов. Критерий Т-оптимальности предполагает гомоскедастичность и нормальную распределенность случайных ошибок наблюдения. В [ ] был предложен критерий ^Т-оп-тималыюсти, основанный на расстояниях Кульбака Лейблера, не требующий нормальности и пригодный для дискриминации двух произвольных конкурирующих моделей для плотностей ошибок. Этот критерий и его статистическая интерпретация описаны в § . Байесовский критерий ^Тр-оптимальности вводится в § по аналогии с байесовским критерием Тр-оптимальности. Соответствующая теорема эквивалентности сформулирована в § 4.3. Численные алгоритмы из главы 3 обобщены на случай нового критерия и описаны в § 4.4. Заключительный § 4.5 посвящен численным примерам. Изложение в данной

главе опирается на материал из работы [16]. Алгоритмы, описанные в § 4.4, также реализованы в пакете [17] для языка Р.

Глава посвящена так называемым полу-параметрическим .^Т-оптималь-ным планам. Рассматривается только случай дискриминации двух конкурирующих регрессионных моделей. В работах по планированию дискриминационных экспериментов, как правило, исследуется ситуация, когда конкурирующие модели, определяемые парой из функции регрессии и закона распределения для случайных ошибок, известны с точностью до параметров. В работе [11] был предложен полу-параметрический критерий оптимальности для планирования дискриминационных экспериментов, обобщающий критерий КЬ-оптимальности из [10], суть которого заключается в том, что функции регрессии для всех конкурирующих моделей и закон распределения ошибок для одной из них считаются известными, а закон распределения для оставшейся модели получается как решение специальной задачи вариационного исчисления с ограничениями. В § 5.1 рассмотрены два варианта полу-параметрического критерия из [11], которые в данной работе обозначаются как 8КЬ^ау и б^Т^-критерии, а также доказана теорема, упрощающая нахождение планов, оптимальных с точки зрения этих критериев. Упрощение заключается в том, что предложенный в теореме подход требует решения одного нелинейного уравнения на заданном интервале для вычисления расстояния Кульбака-Лейблера в критерии оптимальности вместо решения системы из двух нелинейных уравнений. В § 5.2 формулируется теорема эквивалентности для и б^Т^-критериев. В § доказаны две теоремы, связывающие критерии и б^Т^-оптимальиости с критерием Т-оптимальиости. Первая теорема утверждает, что б^Т^-оптимальиые планы совпадают с Т-оптимальными планами в случае, когда плотность случайных ошибок первой модели является усеченной и симметричной относительно условного среднего. Вторая теорема утверждает, что б^Т^-оптимальиые планы совпадают с Т-оптимальными планами, если ошибки второй модели распределены нормально. Наконец, в § 5.4 приводятся численные результаты.

В заключении сформулированы основные результаты работы.

Глава 1. Постановка задачи планирования дискриминационных экспериментов и обзор известных результатов

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

1.1 Задачи анализа и планирования в контексте восстановления

регрессионных зависимостей

Пусть зависимость результатов измерений у\,... , у^ е К от условий их проведения ... , х^ описывается следующим уравнением:

уг = ^(ж<,0) + £г, ж, ех ,(9 е 0, г = 1,...д, (1.1)

где N — общее число проводимых измерений, — вещественнозначная

функция, в = (^1,... , в^)т — ее неизвестные пара метры, £1,..., _ некоррелированные случайные ошибки с нулевым математическим ожиданием и заданной дисперсией. От требования о некоррелированности случайных ошибок можно отказаться (смотри, например, [18]), но в рамках данного текста этот случай рассматриваться не будет.

Функцию принято называть регрессионной моделью, а уравне-

ние (1.1) уравнением регрессии (смотри, например, книгу [19] с.192). Регрессионная модель называется линейной (по параметрам), если она пред ставима в виде

= вт/(Ж), /(х) = [/1 (ж),... Мх)]т, (1.2)

гДе Л(х),... , fd{x) есть набор вещественнозначных линейно независимых функций. В противном случае модель называется нелинейной.

В этом параграфе рассматриваются две задачи, связанные с уравнением (1.1): задача анализа эмпирических данных и задача оптимального планирования эксперимента. Постановки этих задач связаны с двумя соответствующими типами прикладных исследований: наблюдательными исследованиями (observational studies), где сам исследователь никак не вмешивается в процесс сбора данных и работает с уже готовой выборкой, и экспериментальными исследованиями (experimental studies), где исследователь может непосредственно влиять на процесс сбора данных посредством изменения условий проведения эксперимента.

Задача анализа эмпирических данных в контексте уравнения (1.1) состоит в построении функции регрессии по имеющейся выборке

(Ж1,Ш),..., (xn,vn)

из совместного распределения х и у7 то есть при определенных предположениях о природе случайных ошибок (например, что £1,..., являются нормальными случайными величинами с постоянной дисперсией) необходимо найти такую функцию ^(х,0) (и оценить при этом ее неизвестные параметры), которая приближала бы имеющиеся выборочные данные достаточно хорошо относительно какого-нибудь заранее заданного критерия качества ^-(0), например, суммы квадратов разностей реальных наблюдений и предсказаний

1 N 1 N

п(0) = ^ Е с [^(х^)] = ^ Е^ - ^(х^)]2>

г=1 г=1

имея при этом также хорошую обобщающую способность (грубо говоря, функция, которую мы ищем, не должна быть слишком сложной). Более подробное описание этой задачи можно найти, например, в соответствующих главах книг [20] и [21].

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

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

предполагается известной a priori с точностью до параметров в. Задача планирования эксперимента состоит в том, чтобы с точки зрения некоторого заранее заданного критерия оптимальности рационально выбрать условия для проведения измерений, то есть точки ... , xn- Так, например, объем доверительного эллипсоида для оценок неизвестных параметров в по методу наименьших квадратов в случае линейной по параметрам регрессионной модели вида (1.2) при нормальных ошибках с одинаковой дисперсией является монотонно возрастающей функцией от величины

Список литературы диссертационного исследования кандидат наук Гученко Роман Александрович, 2018 год

У Р -

У Рг,П [1г(Х,вг) - Щ (Х,ви ^ ((1вг,3 )

тах

Х£Х '-' I Р* (£*)

г,з=1 (^ )

принимает минимальное значение. Определим также функцию

У Р -

= У] Рг,П [Пг(Х,вг) - Щ (Х,вг,3 Д 2¡Ц3 ((10^) . (3.13)

г,,=1 (а*)

Теорема 8. Пусть вы,полнено предположение . Если план £ не является Тр-оптимальным, тогда существует точках £ X, такая что

) >тР (С).

Доказательство теоремы 8 также можно найти в [15]. Рассмотрим

Предположение 5. Для любых % и ], таких что р^ = 0, и, для любого плана, для которого Тр(£) > 0 ^в^^в ( ) достигается в единственной точке = (£) во внутренности множества О,-.

Если кроме предположения 4 выполнено также предположение 5, то функция

3.4 Численные алгоритмы для построения Тр-оптимальных планов

В работе [13] был предложен численный алгоритм для построения локальных Тр-оптимальных планов, основанный на многомерной чебышевской аппроксимации. Этот алгоритм достаточно сложен как с точки зрения описания, так и с точки зрения реализации, поэтому в этом параграфе он не приводится. Также он требует существенных вычислительных ресурсов (смотри в [13] параграф 7 с численными результатами), из-за чего применять его можно только в случае небольшего количества попарных сравнений между конкурирующими моделями в критерии (3.5). В этом параграфе мы опишем более эффективный метод, способный вычислить план даже в случае большего числа попарных сравнений. Как уже отмечалось в параграфе 3.2, этот алгоритм также можно будет использовать для построения байесовских планов при дискретных априорных распределениях.

Если выполнено предположение , то по теореме существует точкам £ X 7 такая что

для любого плана который не является локальным Тр-оптимальным, где функция ) определяется в ( ), а Тр(£) в ( ). Алгоритм, предложен-

ный в [6], использует похожее свойство для того, чтобы построить последовательность планов, сходящуюся к Т-оптимальному. Приведем аналог этого алгоритма для случая Тр-оптимальных планов.

Алгоритм 2 (аналог алгоритма из [ ]). Пусть ^о ~ некоторый начальный план и (а3)с^=0 — последовательность положительных вещественных чисел,

Более того,

) = Е ЫХ,°г) - Щ (хАз )]2.

(3.14)

) >Тр (С),

удовлетворяющая условиямИш8а = 07 ^= ТО < то. Тогда

для, в = 0,1,... поочередно выполняем два шага:

1. Находим, ж5+1 = а^шахж€д> Ф(ж,<^)7

2. Берем, £8+1 = (1 - а8)+ (ж^).

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

Теорема 9. Пусть вы,полнено предположение 4 и пусть {^}в=0,1,... — это последовательность планов, получаемая, с помощью алгоритма 2. Тогда

НшТр(&) = тр(£ *),

в ^то

где — локальный ТР-оптимальный план.

Доказательство этой теоремы аналогично доказательству теоремы 3 из главы 1.

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

^ = а^т£/ [^г(жД) - ^(ж, 0г^)]2£(¿¿ж) (3.15)

е^Зх

в критерии (3.5). Также в работе [13] было показано, что алгоритм 2 может не сходится в том случае, если в критерии (3.5) слишком много попарных сравнений.

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

критерия по опорным точкам (шаг 1) и по весам (шаг 2). Первый шаг не представляет затруднений и подробно мы его описывать не будем. Для второго шага предложены несколько вариантов оптимизационных процедур, которые мы опишем ниже по тексту в параграфах 3.4.1 и 3.4.2.

Алгоритм 3. Пусть £0 — некоторый начальный план. Тогда для в = 0,1,... поочередно выполняем два, шага:

1. Пусть — носитель плана Найдем множество Е^] всех локальных максимумов функции на X и определим = 5[5] иЕ^.

2. Определим £ = {Б^+^ш} как план с носителем и вектором, весов ш. Теперь вычислим локальный, Тр-оптимальный план в классе всех планов с носителем то есть найдем такой вектор ш^+ц, которым доставляет максимум

где — это вес в точке х £ Б^+ц. Уберем из множества те точки, которым в полученном оптимальном векторе весов ш^+ц соответствуют нули, и, получим, новый носитель, которым также обозначим, Наконец, определим план^а+1 как план с носителем

и соответствующими ненулевым и весами изш^+ц.

На первом шаге каждой итерации все локальные максимумы функции добавляются к множеству возможных опорных точек плана Накапливания точек в носителе не происходит благодаря следующему важному свойству функции ):

Лемма 3. Пусть вы,полнены, предположения, 4 и, 5. В конце каждой, итерации, алгоритма функция, ^(х£8+1) принимает одно и то же значение во всех опорных точках плана

Доказательство леммы 3 можно найти в [15]. Из этой леммы сразу следует

Теорема 10. Пусть вы,полнены, предположения, 4 и 5. Тогда в конце каждой, итерации алгоритма число точек в носителе план а ^^ не превосхо-

д(щ) = Тр({5[8+1],^},01,... А)

dum максимального числа корней уравнения

6+1) = с, се [0, тахФ^ £s+i)].

В программной реализации на каждом шаге из носителя промежуточного плана исключаются те точки, итоговый вес которых после решения оптимизационной задачи по весам w меньше, чем т0'25, где т = 2.2 х 10-16 — самое маленькое число с плавающей запятой в R.

Для алгоритма 3 верен следующий результат о (ходимости:

Теорема 11. Пусть выполнено предположение 4 и пусть {s=0;1;... — это последовательность планов, получаемая с помощью алгоритма 3. Тогда

НшТр (&+1)=ТР К *),

в ^то

где — локальный Тр-оптимальный план. Доказательство. Очевидно, что неравенство

ТР ({^[ф^Н}) ^ ТР ({^[5+1],^[5+1]})

выполнено для всех в, так как С 5[5+1], то есть последовательность Тр(^)5 в = 0,1,... — неубывающая. Более того, эта последовательность также ограничена сверху величиной Тр(£*), а значит она имеет предел. Обозначим этот предел как Тр*. Существует подпоследовательпость планов ] = 1,2,..., которая сходится к некоторому плану £**. Функция Тр является полу-непрерыв-ной сверху как инфимум непрерывной функции, поэтому Тр(£**) = Тр*. Теперь предположим, что

трК**) <трк*). ** Тр

ществует константа 5 > 0, такая что

supQ(C**) -ТРК**) = 2Ö,

где функция Q определяется в формуле ( ). Поэтому при достаточно больших j, таких что j > N для некоторого заранее фиксированного N, получаем (используя тот факт, что функция sup CQ(() является полунепрерывной сни-

supQ(c, ) -тР() >а,

с

при всех j > N. Так как по построению последовательность Тр (£s), s = 0,1,... не убывает, верно неравенство

тр) - ТР() > ТР+1) - ТР). (3.16)

Чтобы оценить правую часть этого неравенства, определим при j > N и a £ [0,1] план

4+1(a) = (1 -a) + a0,

где Q — это мера, на которой функция Q(£,) достигает своего максимального значения на множестве всех планов, сосредоточенных в точках локальных максимумов функции Ф(ж,), и

as = argmax^(£ (a)).

0<a<1

По построению план +1 — это лучший план, сосредоточенный на множестве точек supp() U supр(Q), поэтому

Тр(60+1)) > ТР+1) > Тр(¿,+1 (aSj+i)). (3.17)

Введем обозначение

ВД,а)=ТР (L (a)),

и заметим, что

<9ТР (¿,+1 (a))

5a

= Q(0,) -тр() = supQ(C, ) - Тр() > £

a=° С

Теперь, используя разложение Тейлора, получаем

+ 1,о^.+1) - ни + 1,0) = тах [Тр(£,+1 (а)) - Тр^ (0))]

а£[0,1]

> тах

«£[0,1]

> тах

«£[0,1]

■ дТр( (а)) а—

да

а6 - -а2К 2

- -а2К

а=0 2

2

2 К

К

ния (3.17) получаем

тр(с^) - Тр(с) > Тр(с+1) - Тр(с)

>7р(^.+1) -Тр(с )

2

= к(з + 1,а81.1) - + 1,0) >

2 К

Теперь просуммируем последнее выражение по ] от N до Ь - 1 > N при Ь > N + 1 и получим

Ь-1 62 Е [Тр(^(,+1)) - Тр(^)] = ТР() - ТР(и) > [Ь - N1 2К.

Левая часть последнего неравенства стремится к конечному значению Т(£**) -Т), а правая часть стремится к бесконечности при Ь ^ ж. Мы пришли к противоречию, а значит наше исходное предположение о том, чтоТр(^**) < Тр(£*) неверно, что завершает доказательство теоремы . □

С практической точки зрения эффективное решение оптимизационной задачи на втором шаге алгоритма 3, то есть нахождение максимума функции д(ш), является самым важным. Далее описаны два метода оптимизации по весам. Будем считать, что выполнены предположения 4 и 5.

3.4.1 Квадратичное программирование

Пусть = {х1,... ,хп} обозначает множество, получаемое на первом

шаге й-й итерации алгоритма . Определим £ как план с носителем 5^+1] и весами ... ,шп, которые и предстоит определить на втором шаге итерации

алгоритма посредством максимизации функции

V п

) = Е ^¿^^[^(жА; ,0 г) - ^ (жА Д,?)] ,

г,_7=1 а=1

Где = определяются в (3.15). Для этой цели предлагается линеаризо-

вать функции т^- (ж*,0г,^) в окрестности точек 0г,р Рассмотрим функцию

5М = £ шл £ [*(ж* Д) - „ (ж* Д„) - ^

■ . , а ,

г^=1 ^ А=1

V

д

. . ац

г^=1

д щ (жг, 6>г^)

где ^ — это размерность множества возможных параметров j-й модели О^

^ = diag(.l,... ,.п), а матрицы J^ £ Мпх^'5 Кг^ е М пхи векторы ^ £ Мп задаются как

^,3 V д=0^/ Г=1

^ = ([Д) - V, (жГ Д, )] й . )

\ д(7г. 1 Рг.1 = Ргл'/ 5

= ([??г(жг Д) - Щ (жг А^Р^

соответственно. Нетрудно видеть, что минимум по достигается при

/ Т1 \ — 1 Т1

=1..... п

^ = ^,У Г=1,...,п

откуда получаем, что

#(.) = -. +

где

= > Рг.г^г.г иТ'.Шг.,-)

2

Матрица Q(w) зависит от ш, но если игнорировать эту зависимость и взять матрицу

Q = diag(^i,... ,шп)

фиксированной, то получим следующую задачу квадратичного программированы^) = — шTQ(w) ш + bTш ^ max, (3.18)

ш

п

шк = 1; шк > 0, к = 1,... ,п.

k=i

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

Отметим, что получение оптимальных значений (3.15) требует в алгоритме 3 гораздо меньше вычислительных ресурсов, чем в алгоритме 2, так как количество опорных точек в промежуточном плане, который находится в ходе итерации, существенно меньше за счет удаления точек с нулевыми весами.

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

При практической реализации предложенной процедуры для шага 2 рекомендуется проводить только первые несколько итераций и прерывать процедуру, как только будет достигнута существенная разница между значениями критерия при начальном для шага 2 плане и плане, полученном в результате итерационного процесса. При этом общая скорость сходимости алгоритма 3 существенно возрастает, однако равенство значений функции Ф для опорных точек вычисленного на шаге 2 плана (как указано в лемме 3) может и не достигаться.

3.4.2 Градиентный метод

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

V п

д(ш ) = Е Рьз^шк[г1г(хк ,вг) - ц3 (хк ,6i,j )] , (3.19)

i,j=1 к=1

где Qij = 6ij(ш) определяются в (3.15). Введем функции

V

vk (ш) = Pi,j [пг(хк, 0г) - Vj (хк,0г,з (ш))] , к = 1,...,п,

и посчитаем последовательность векторов (ш(7))7=од;... по следующему алгоритму. Вначале возьмем какой-нибудь вектор ш(о) = ш (например, с равными весами). Если ш(7) = (ш(7),1,... ,ш(1),п) задан, будем производить вычисления для 7 = 0,1,... по следующей схеме. Найдем индексы кик, такие что

к = arg max Vk (ш(7)), к = arg min Vk (ш(7)).

1<k<n 1<к<п

Обозначим

а* = arg max д(ш (7 )(а)), (3.20)

где вектор ш(7)(а) = (ш(7у (а),... ,ш(7),п(а)) определяется как

|ш(7у + а, если i = к;

— а, если i = к; ш(7),г, иначе.

Вектор ш(7+1) на следующей итерации определяется как ш>(7+1) = ш(7)(а*). Следующая теорема показывает, что сгенерированная таким образом последова-

ш*

Теорема 12. Последовательность (ш(7))7ем сходится кш* G argmaxg(ix>).

Доказательство теоремы можно найти в [15].

Стоит отметить, что одномерная оптимизационная задача (3.20) является достаточно дорогой с вычислительной точки зрения. В программной реализации этого метода мы используем линеаризацию подобно той, что была описана в параграфе 3.4.1.

3.5 Примеры

Теперь проиллюстрируем предлагаемые алгоритмы с помощью двух примеров. Алгоритм 3 с двумя способами оптимизации по весам из параграфов 3.4.1 и был реализован в качестве пакета [ ] для языка К Примеры из этого параграфа также можно найти в пакете. В программе пользователь задает начальный план £о, метод решения оптимизационной задачи на шаге 2, весарг,^ и априорную информацию о параметрах #г. Таблицу

Р =

Р1,1 Р1,2 ... Р1, !/-1 Рм

... !/-1

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

Тр

Тр

кретных априорных распределений для параметров может быть сведена к зада-

Тр

рующих моделей и попарных сравнений, мы опишем соответствующую таблицу сравнения для байесовского критерия. Для наглядности в первом примере мы рассмотрим случай двух конкурирующих моделей и = 2. Байесовский крите-

Тр

распределение на О1 имеет носит ель А1,...и вес а Т1,... т^. Этот критерий может быть переписан в виде

1+1 » 2 Тр(0= Е ^ /пи Ьг(жД) - Щ (ж, №), (3.21)

г ¿=1 ^ * I- ]

где таблица сравнения задается как

0 0 ... 0 п 0 0 ... 0 т2

P = (Pij) ¿J=1,...,l+1

G Rl+1xl+1, (3.22)

0 0 ... 0 П 0 0 ... 0 0

т]г(ж,0г) = ^(ж, Аг), г = 1,..., £ и т]1+1(ж,0^-) = т]2(ж,6'г,1+1). Таким образом, байесовский критерий для дискретного априорного распределения предусматривает

Построение таблиц сравнения для случая дискриминации более чем двух моделей делается по аналогии.

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

3.5.1 Байесовские Тр-оптимальные планы для дискриминации

экспоненциальных моделей

Рассмотрим задачу дискриминации между двумя конкурирующими экспоненциальными регрессионными моделями

где множество планирования — это интервал [0,10]. Экспоненциальные модели вида (3.23) широко используются в различных приложениях. Например, модель 2

альное название (Mitscherlich growth law), при описании зависимости урожайности от количества использованного удобрения. В исследованиях в области рыболовства эта же модель носит название закона роста Берталанффи и используется при описании зависимости длины рыбы от ее возраста (смотри [33]). Оптимальные планы для экспоненциальных моделей обсуждались, например, в

щ(ж,01) = 01,1 - 01,2exp(-01>эж01-4), ^(ж,02) = $2,1 - 02,2 exp(-02>3Ж),

(3.23)

работе [ ]. Нетрудно показать, что Тр-оптимальный план для дискриминации моделей ( ) не зависит от линейных параметров модели щ, поэтому можно выбрать их произвольными, например, $1,1 = 2 и 6^,2 = 1. Для параметров $13 и $1,4 были взяты независимые априорные распределения, сосредоточенные в точках

д + г = 1,... ,5 ; j = 3,4 , (3.24)

где д3 = 0.8, д4 = 1.5. Веса в этих точках пропорциональны (в обоих случаях) величинам

, 1 expf - ( ~ 3) ); г = 1,..., 5 . (3.25)

Р\ 8 ); , , К }

25

Тр

минации моделей (3.23), имеющие эффективность

EffтР(0 = ТРТ\ ) > 0.999,

были найдены тремя способами: 1) алгоритмом 3 со вторым шагом из § 3.4.1, 2) алгоритмом 3 со вторым шагом из § 3.4.2, и 3) алгоритмом 2. В таблице 4 представлены планы, полученные первым способом, при различных значениях параметра а2, а в таблице — сравнительное время работы в секундах до получения плана с нужной эффективностью для всех трех алгоритмов. В качестве

11

0,1, . . . ,10

Вычисления проводились на персональном компьютере с центральным процессором intel core i7-4790K. Из таблицы 5 видно, что предложенные в настоящей главе алгоритмы значительно превосходят стандартный алгоритм по времени работы.

Отметим, что в случае небольшой дисперсии а2, оптимальные планы имеют четыре опорные точки, а при а2 > 0.285 байесовские Т^-оптимальные планы сосредоточены в пяти точках. Графики функции Ф из теоремы эквивалентности показаны на рисунке при разных значениях дисперсии а2.

Таблица 4 — Байесовские Тр-оитимальные планы для дискриминации моделей (3.23)

а2 Оптимальный план

0.0 0.000 0.441 1.952 10.000 0.209 0.385 0.291 0.115

0.1 0.000 0.452 1.877 10.000 0.209 0.391 0.290 0.110

0.2 0.000 0.455 1.811 10.000 0.208 0.394 0.291 0.107

0.285 0.000 0.453 1.758 10.000 0.207 0.396 0.292 0.105

0.3 0.000 0.452 1.747 4.951 10.000 0.207 0.396 0.292 0.003 0.102

0.4 0.000 0.446 1.651 4.699 10.000 0.200 0.384 0.290 0.060 0.066

Таблица 5 Время работы различных алгоритмов в секундах

Пример а2 Алг. 3(1) Алг. 3(2) АФ

(3.23) 0.0 0.03 0.11 12.4

0.4 1.4 11.6 218.3

(3.26) 0.0 0.09 0.75 5.7

372 7.8 37.1 762.3

Примечание Алг. 3(1) и Алг. 3(2) алгоритм 3 с использованием методов оптимизации по весам из параграфов 3.4.1 и 3.4.2 соответственно; АФ алгоритм 2

о о

10

2 4 6 8

х

10

а) а2 = 0

б) а2 = 0.1

0

2

4

6

8

0

х

д) а2 = 0.3 е) а2 = 0.4

Непрерывной линией обозначен график функции Ф(ж,<^*2), задаваемой

формулой ( ), пунктирной — значение функционала Тр (<^*2), определенного в ( ), а выколотыми точками — опорные точки Тр-оптимальпого плана £*2. Рисунок 3.1 — Графики функций Ф(ж,<^*2) для планов из таблицы

3.5.2 Байесовские Тр-оптимальные планы для дискриминации моделей зависимости эффективности препарата от дозы

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

щ(х, Ох) = 0М + 0^х;

T|2{x, 02) = 02,1 + 02,2х(02,3 - х)

(3.26)

щ(х, 03) = 03,1 + 0з}2х/(03,3 + х);

щ(х, 04) = 04,1 + 04,2/(1 + ехр(04,3 - х)/64,4);

где множество планирования есть интервал X = [0,500]. В этой работе также приводились априорные значения для параметров этих моделей

01 = (60,0.56),

02 = (60, 7/2250, 600),

03 = (60, 294, 25),

04 = (49.62, 290.51,150,45.51).

Локальный Тр-оптимальный план для моделей ( ) был получен в работе [ ] в случае = 1/6 (1 < ] < ъ < 4) то есть критерий Тр-оптимальности (3.5)

6

Демонстрацию новых методов, предложенных в параграфе 3.4, мы начнем

Тр

моделей (3.26), что был представлен в работе [13]. Алгоритму 3 со вторым ша-

0.09

плана с эффективностью > 0.999. В качестве стартового использовался стартовый план из [ ]. График функции Ф(х,£1) после первой итерации алгоритма представлен на рисунке . Опорные точки плана £1 показаны на этом рисунке выколотыми точками. Отметим, что во всех опорных точках £1 функция Ф(х, £1) принимает одно и то же значение, как и утверждается в лемме . Чтобы найти план с аналогичной эффективностью, алгоритму из статьи [13] понадо-9

0 100 200 300 400 500

x

Рисунок 3.2 — График функции Ф(ж, £1) после первой итерации алгоритма

Тр

план. Для простоты мы задаем априорное распределение только для четырехмерного вектора параметров последней модели $4, в то время как 02 и 03 остаются фиксированными векторами. Это априорное распределение задано в 81 точке из К4 следующим образом:

{Ме ,е2,ез,е4 1 ех,e2,eз,е4 € {-Щ1}}, (3.27)

Меье2,ез,е4 = (М1 + в^, Д2 + Дз + 63^,^4 + в4^),

М = (М1,М2,М3,М4) = (49.62, 290.51,150,45.51).

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

1 „„„ f ll^ei,е2,ез,е4 М||2

(2тгa2)2 2а2

) , ei,е2,ез,е4 G {-1,0,1}, (3.28)

где || • ||2 есть евклидова норма. Конечный байесовский критерий ( ) состоит 246

Для разных значений параметра а2 байесовский Тр-оптимальный план снова был найден тремя способами (как в предыдущем примере). Планы, найденные с помощью алгоритма 3 со вторым шагом из § 3.4.1, представлены в

Таблица 6 — Байесовские Тр-оптимальные планы для

дискриминации моделей (3.26)

а2 Оптимальный план ^

0 0.000 78.783 241.036 500.0 0.255 0.213 0.357 0.175

202 0.000 84.467 234.134 500.0 0.257 0.225 0.351 0.167

302 0.000 91.029 225.713 500.0 0.259 0.237 0.345 0.159

332 0.000 92.692 222.735 500.0 0.260 0.240 0.344 0.156

352 0.000 91.743 129.322 221.118 500.0 0.260 0.214 0.036 0.336 0.154

372 0.000 89.881 129.590 170.306 220.191 500.0 0.260 0.170 0.091 0.019 0.310 0.150

таблице 6. Сравнительное время работы для всех трех алгоритмов можно найти в таблице 5. Мы снова видим, что предложенные в главе методы работают существенно быстрее. При небольших значениях параметра а2 байесовский Тр-оптимальный план сосредоточен в 4 точках считая граничные точки множества планирования. При увеличении параметра а2 количество опорных точек растет. Так, если а2 = 352 или 372, байесовский Тр-оптимальный план имеет 5 или 6 точек, считая граничные точки множества планирования. Графики функции Ф из теоремы эквивалентности приведены на рисунке

100

200

300

400

500

о о о

100 200 300 400 500

х

а) а2 = 0

б) а2 = 202

0

0

X

о о о

100 200 300 400 500

х

в) а2 = 302

о о о

о о ю

100

200

300

400

500

а2 = 332

0

0

х

о о о

о о о

о о о

100

200

300

400

500

100

200

300

400

500

д) а2 = 352 е) а2 = 372

Непрерывной линией обозначен график функции Ф(ж,£*2), задаваемой

формулой ( ), пунктирной — значение функционала Тр (£*2), определенного в ( ), а выколотыми точками — опорные точки Тр-оптимальпого плана £*2. Рисунок 3.3 — Графики функций Ф(ж,£*2) для планов из таблицы

0

0

х

X

Глава 4. Построение байесовских ^Тр-оптимальных планов

Критерий Т-оптимальности (как и обобщающий его критерий Тр-опти-малыюсти) предполагает гомоскедастичпость и нормальную распределенность случайных ошибок наблюдения. На случай дискриминации двух моделей при ге-тероскедастичных нормальных ошибках Т-критерий был обобщен в работе [ ]. В [ ] был предложен критерий .^Т-оптимальности, основанный на расстояниях Кульбака Лейблера, не требующий нормальности и пригодный для дискриминации двух произвольных конкурирующих моделей для плотностей ошибок, причем критерии из работ [6] и [9] являются его частными случаями.

В этой главе предложен метод построения байесовских КЬр-оптимальных планов для дискриминации произвольного количества регрессионных моделей. Байесовский критерий ^Тр-оптимальности вводится по аналогии с байесовским критерием Тр-оптимальности. Численные алгоритмы из главы обобщаются на случай нового критерия. Изложение в данной главе опирается на материал из работы [16]. Алгоритмы, описанные в параграфе 4.4, также реализованы в пакете [17] для языка Р.

4.1 Обобщение критерия Т-оптимальности на случай произвольно

где х Е а случайные ошибки £ независимы и имеют нулевое средние и заданную дисперсию. В работе [10] было предложено следующее обобщение критерия Т-оптимальности для случая произвольно распределенных ошибок наблюдения. Будем считать, что при фиксированных х и 0 величина У имеет абсолютно непрерывное распределение с плотностью /(у,ж,0), среднее и дисперсия которой равны

распределенных ошибок наблюдения

Рассмотрим стандартное уравнение регрессии:

У = ^(ж,0) +

(4.1)

Рассмотрим теперь две конкурирующие модели для плотности /(у,ж,0):

/1(2/,ж,0х), у € у, ж € *, 01 = (01,1,...,01,,1 )т € ©х; /2(2/,Ж,02), у € У, ж € *, 02 = (02,1, . . . ,$2,¿2)Т € ©2,

где У — общая область определения для обеих плотностей. Расстояние Куль-бака-Лейблера между конкурирующими плотностями/1 (у,ж,0^ и /2(у,ж,02) обозначим как

М^А,$2)= / /2(2/,жА)^^^^ (4.2)

}у Л (у ,ж,01)

Предположим теперь, что параметры первой модели фиксированы и обозначим 1

/1 (у,ж,0^ и /2(у,ж,02) в [ ] было предложено искать такой план эксперимента, который максимизирует величину

КЬ1,2Й, $1)= ^ / /1,2(ж, 01, 02)^(^). (4.3)

План, максимизирующий критерий ( ), называется ^Ь-оптимальным.

Отметим, что критерий (4.3) не является симметричным. По сути, поменяв /1 и /2 в () местами или заменив т£02€©2 на т£б»1€©1 в ( ), мы получим другие критерии (всего четыре разных варианта).

Статистическая интерпретация критерия ^Ь-оптимальности такова. Рассмотрим задачу проверки гипотезы

Яо :/(ж,0) = Л(ж,01),

против а л ьт е р и ат и в ы

Я1 :/(ж,0) = /2(ж,02), ©2.

Отношение правдоподобия для наблюдения у в точке ж имеет следующий вид:

Л(2/ 1)

ь( ¿Ш =

/2(У^ $2)'

Тест отношения правдоподобия основан на статистике

/2(2/02)

Л(0 2, &) = -2^(0 2,&) = 2^

Л(2/

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

ЕЯ1Л(02,Ы = 2 / /2(2/,х,02)log ^^ = 2/1,2(х,01, 02)^.

Мощность теста отношения правдоподобия тем больше, чем больше значение Ен1Л(02,£ж)- Теперь рассмотрим точный план£^ = (х1,..., хп; д1,..., дп), такой что ^ = т^/Ж, ?=1 = Ж. Используя независимость наблюдений получаем

ЕЯ1^(0 2 ) = 2 / /2(2/ п{ ^2/

/ 7=1 1/1(2/ ,Хг, 01 )J

= 2 Ж ¿Ф / /2(2/,х,02)log ^

/ /1(2/ ,Хг, 01)

= 2 Ж [ /1,2(Х,01, 02)^(<&). Jx

Переходя от приближенного плана к точному плану ^ и рассматривая значение Ея1Л(02,£) при наименее благоприятном 02, когда Ея1Л(02,<^) минимально, получаем критерий (4.3).

Теперь рассмотрим обобщение критерия ^Ь-оптимальности на случай дискриминации произвольного количества конкурирующих плотностей

/1(2/ ^ ^О,..., Л(У ^).

Фиксированные параметры рассматриваемых моделей обозначим как01,..., 91У.

I,, (хД, в, ) = | д (г/ ,Х, ^ )log ||;Х|у

а соответствующий критерий

КЦ,-(£Д)= / 4,(хД,%)^х).

сти происходит оптимизация. Следуя обозначениям главы 3, пусть — это неотрицательные веса, отвечающие за вклад попарного сравнения плотностей /г(у,х,0«) и (у,х,0^) в финальном критерии. По аналогии с тем, как мы делали это в предыдущей главе для критерия Т-оптимальности, введем симметричную версию критерия ХЬ-оитимальности для V > 2 конкурирующий моделей /1(у 01),..., ^ (у ^):

V V г,

КЬр(0 = £ Д) = £ ы \ /^хДД^х). (4.4)

План, максимизирующий критерий ( ), будем называть ХЬр-оптимальным планом для дискриминации моделей /1 (у,х,01),... ,/V(у,х,0г). Будем далее всегда предполагать, что множества^, в^ компактны и все плотности непрерывны по х. Отметим, что критерий ( ) зависит от неизвестных на практике параметров 01,..., 0 которые устанавливаются экспериментатором для конкурирующих моделей /1(3/,х,01),...,/V(у,х,0г). Следовательно, критерий является локальным в смысле работы [27]. Если фиксированные значения параметров 01,... ,0г далеки от истинных значений, найденный оптимальный план может иметь низкую эффективность.

4.2 Байесовские ХЬр-оптимальные планы и их сведение к

локально оптимальным

Аналогично тому, как мы делали это в предыдущей главе для Тр-опти-мальных планов, сформулируем байесовскую версию критерия ХЬр-оптималь-ности, являющуюся более устойчивой по отношению к неправильному выбору фиксированных параметров. Байесовским ХЬр-оптимальным планом будем на-

зывать дискретную вероятностную меру, доставляющую максимум величине ХЬВ (0 = ¿Я, / ^ К, (4.5)

= ^2 Рч / 1г( ж^^Ж^^да

^ = 1 М ^ X

Здесь р^ > 0 и р^ = 0 при любом г, а мера Тг для каждого г = 1,... задает априорное распределение для параметров ^ плотности такое что все интегралы в формуле (4.5) определены.

Отметим, что для случая двух конкурирующих моделей критерий (4.5) рассматривался в работе [36]. Как было показано в [15], вычисление байесовских планов, оптимальных с точки зрения критерия (4.5), можно свести к вычислению локально оптимальных планов для большего количества конкурирующих моделей. В большинстве приложений интегралы в (4.5) вычисляются численными методами, предполагающими аппроксимацию произвольных априорных распределений мерами с конечным носителем. Если априорные распределения "Рг дискретны и сосредоточены в точках Ал,... с весам и т^,... т^, то критерий (4.5) можно переписать в виде

V I „

^(0 = £ £й / (ж, А* Д,)№). (4.6)