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



На правах рукописи

Беляев Михаил Геннадьевич

Построение аппроксимации многомерных зависимостей по выборкам с факторным

планом

05.13.17 — Теоретические основы информатики

АВТОРЕФЕРАТ диссертации на соискание учёной степени кандидата физико-математических наук

"ОЯ2ШЗ

Москва — 2013

005538988

Работа выполнена в Федеральном государственном бюджетном учреждении науки Институте проблем передачи информации им.

A.A. Харкевича Российской академии наук (ИППИ РАН)

Научный руководитель:

доктор технических наук, академик РАН Кулешов Александр Петрович Официальные оппоненты:

доктор физико-математических наук, профессор Магарил-Ильяев Георгий Георгиевич, профессор кафедры Общих проблем управления механико-математического факультета Московского государственного университета имени М.В. Ломоносова

кандидат физико-математических наук, доцент Гасников Александр Владимирович, доцент кафедры Математических основ управления Московского физико-технического института (государственного университета) Ведущая организация:

Федеральное государственное бюджетное учреждение науки Институт прикладной математики им. М.В. Келдыша Российской академии наук

Защита состоится 16 декабря 2013 г. в 16.30 на заседании диссертационного совета Д 212.156.04 при Московском физико-техническом институте (ГУ) по адресу: 141700, г. Долгопрудный, Московская обл., Институтский пер., д. 9, ауд. 204 Нового корпуса.

С диссертацией можно ознакомиться в библиотеке Московского физико-технического института (государственного университета). Автореферат разослан 14 ноября 2013 года.

Ученый секретарь

диссертационного совета Д 212.156.04

к.ф.-м.н.

Стрыгин Л. В.

Общая характеристика работы

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

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

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

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

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

P. Dierckx2 поставил задачу построения бикубической аппроксимации с явным контролем гладкости. Задача свелась к решению системы с матрицей вида (M®N+A®I+I®B)a; = у, которая эффективно не решалась (символ <g> обозначает произведение Кронекера, т.е. тензорное произведение матриц). Была предложена модифицированная штрафная функция, которая теряла содержательный смысл, но позволяла строить модель эффективно. Использование такого штрафа приводило к системе вида (М + А) ® (N В):;: = у, которая легко решалась с помощью использования базовых свойств произведения Кронекера.

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

Наконец, L. Xiao с соавторами4 предложили подход, позволяющий строить аппроксимацию для многомерной задачи с одномерными факторами и явно контролировать изменчивость модели с помощью вектора параметров регуляризации. Из недостатков следует

1De Boor С. Efficient computer manipulation of tensor products // ACM Transactions on Mathematical Software. 1979. V. 5, no. 2. P. 173-182.

2 Dierckx P. A fast algorithm for smoothing data on a rectangular grid while using spline functions // SIAM Journal on Numerical Analysis. 1982. V. 19, no. 6. P. 1286-1304.

3Marx B. D., Eilers P. H. C. Multidimensional penalized signal regression // Technometrics. 2005. V. 47, no. 1. P. 13—22.

4Xiao L., Li Y., Ruppert D. Fast bivariate P-splines: the sandwich smoother // Journal of the Royal Statistical Society: Series В (Statistical Methodology). 2013. V. 75, no. 3. P. 577-599.

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

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

Результатов по неполному факторному плану существенно меньше. В частности, G. Pisinger и А. Zimmermann6 рассмотрели задачу построения бикубической интерполяции по неполной двумерной решетке (без какого-либо штрафа на изменчивость модели) и предложили решение, основанное на использовании техники псевдообращения. Метод легко обобщается на многомерный случай, однако помимо отсутствия штрафа у него есть еще один существенный недостаток: необходимо решать систему из М уравнений (М — число пропущенных точек). Эта система не имеет какой-то специальной структуры, что не позволяет использовать этот метод при большом числе пропущенных точек.

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

5R. Paulo Default Priors for Gaussian Processes // The Annals of Statistics. 2005. V. 33, no. 2. R 556-582.

6G. Pisinger, A. Zimmermann. Linear least squares problems with data over incomplete grid 11 BIT Numerical Mathematics. 2007. V. 47. P. 809-824.

""Friedman J. H. Multivariate adaptive regression splines // The Annals of statistics. 1991. V.19, по. 1. P. 1-67.

сложных зависимостей (и соответствующих больших объемов выборки) такое упрощение может существенно повлиять на качество аппроксимации.

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

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

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

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

4. Предложить методы выбора вектора параметров регуляризации для разработанного способа решения задачи аппроксимации.

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

Научная новизна работы состоит в том, что в ней впервые

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

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

• получено вычислительно-эффективное оптимальное решение задачи аппроксимации для выборок с неполным факторным планом.

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

Основные положения, выносимые на защиту:

1. Разработан новый метод решения задачи аппроксимации для случая факторных планов.

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

3. Для нахождения оптимального решения задачи аппроксимации в случае полного факторного плана эксперимента разработана явная вычислительно-эффективная формула.

4. Для нахождения оптимального решения задачи аппроксимации в случае неполного факторного плана эксперимента предложен вычислительно-эффективный итеративный метод.

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

• International workshop on advances in predictive modeling and optimization (2013, Берлин, Германия);

• 9-я Международная конференция "Интеллектуализация обработки информации" (2012. Будва, Черногория);

• 15-я Всероссийская конференция "Математические методы распознавания образов" (2011, Петрозаводск, Россия);

• Конференции молодых ученых "Информационные Технологии и Системы" (2011, Геленджик, Россия; 2012, Петрозаводск, Россия; 2013, Калининград, Россия);

• Third International Workshop on Surrogate Modelling and Space Mapping For Engineering Optimization (2012, Рейкьявик, Исландия);

• 53-я, 54-я и 55-я научные конференции Московского физико-технического института (2010-2012, Долгопрудный, Россия).

Кроме того, полученные в работе результаты обсуждались на научных семинарах лаборатории структурных методов анализа данных в предсказательном моделировании МФТИ (2012, 2013).

Публикации. Основные результаты по теме диссертации изложены в 6 печатных изданиях, 3 из которых изданы в журналах, рекомендованных ВАК ([1]-[3]), 3 — в тезисах докладов. Основные результаты изложены в работах [1], [2], написанных без соавторов, и совместной работе [3], в которой вклад диссертанта был определяющим.

Объем и структура работы. Диссертация состоит из введения, трех глав, заключения и приложения. Полный объем диссертации 98 страниц, включая 7 рисунков. Список литературы содержит 89 наименования.

Диссертационная работа была выполнена при поддержке лаборатории структурных методов анализа данных в предсказательном моделировании МФТИ, грант правительства РФ дог. 11.G34.31.0073.

Содержание работы

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

Первая глава посвящена постановке задачи аппроксимации. В первом разделе формулируется задача аппроксимации по конеч-

ному набору данных в заданном классе функций. Назовем обучающей выборкой S совокупность конечного множества точек £ = {хг е Х.Х С Rd)iLi и значений некоторой неизвестной функции д : X —» R1 в точках из этого множества:

5 = {хг еХ,№= g(xi) = {£, <?(£)}.

Задача (Аппроксимации). Пусть заданы обучающая выборка S, некоторый класс функций F и штрафная функция Р\ : F —» R\. Необходимо найти 6 F, минимизирующую сумму выборочной функции ошибки и штрафа:

г = argmin y>(x)-/(z))2 + Px(f).

Для конкретной постановки задачи аппроксимации необходимо задать класс функций F и штрафную функцию P\{f)- В следующих разделах мы выберем их, принимая во внимание особенности плана эксперимента

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

ак = {4 е Х*}?;=1, xfc с Rd\ к = 1,..., К.

Под полным факторным планом будем понимать множество точек обучающей выборки вида Е = о\ х • • ■ х ак = Е/иИ (d = ^fc)- В случае неполного факторного плана множество точек обучающей выборки S не равно множеству Е¡г1ц, а составляет его достаточно большое подмножество.

В третьем разделе рассматривается вопрос выбора класса функций F. Построим словарь функций Atcnsor с помощью тензорного произведения заданных словарей Д/; = функций меньшей размерности:

A tensor = {Ф)г {jk = 1, • • • , Рк)к=1 } •

Будем искать решение задачи в классе функций Fiensor, который состоит из всех возможных линейных комбинаций функций из Ate„SOr: fix) - EJiJ2,...jk ОЛЛ,...^^ {х1Щ2{х2) ... ^fK{xK). Пе-реиишем это выражение с помощью тензорных операций:

f(x) = А®^1^1) ®2тР2{х2) ... ®кфк(хк).

где х = (х1 .х2,... ,хк), 1[)к(хк) — это вектор-столбец значений функций из словаря Д*. в точке хк, а А — тензор коэффициентов разложения по словарю A tensor размера п\ х • • - х пц- Символ ¡^используется для обозначения умножения тензора на матрицу вдоль направления8. Для двумерных тензоров (т.е. матриц) эта операция сводится к умножению слева или справа: Л®^ = ВТА, = АВ.

Четвертый раздел посвящен выбору штрафной функции P\(f), которая позволила бы контролировать гладкость модели относительно различных групп переменных. Рассмотрим сначала способ ограничения изменчивости модели вдоль срезов.

Определение. Будем называть срезами модели (функции) f(x) вдоль переменных х1 следующее множество функций

{fi2,...,ik(xl) = г2 = 1,п2, • ■ • , гк = 1,пк}.

где /(ж1, х?2,..., xf) функция от х1 при зафиксированных переменных х2 — xf ,..., х" = х? .

12' ' гк

Определение естественным образом обобщается на срезы вдоль переменных из нескольких факторов.

Мы будем использовать штраф, основанный на ограничении изменчивости модели вдоль срезов. Рассмотрим в качестве примера штраф для срезов вдоль ж1. Под изменчивостью будем понимать интегральную норму второй производной, что является весьма распространенным подходом и используется, например, в thin-plate сплайнах9. Введем обозначения для переменных, составляющих этот фактор: х1 = (х1'1,.... х1''1'). Тогда определим штраф

8см., например, Kolda Т. G., Bader В. W. Tensor decompositions and applications // SIAM Review. 2009. V. 51. №. 3. P. 455-500.

9Wood S. N. Thin plate regression splines / / Journal of the Royal Statistical Society: Series В (Statistical Methodology). 2003. V. 65, no. 1. P. 95-114.

рРл (/) как сумму по всем срезам и всем возможным частным производным порядка 2 (индекс е\ описывается ниже):

£¿1 (¿1

рм = Е ЕЕ

¿2,9=1 г=1

д2П2.....ь(х1)

дхдх1<г

/V \ их-™ их'

,,гк (?=1 г=1

где е\ = (1,0,.. -, 0) € в является бинарным вектором длины К,

© {0,1}^, а ре1(1) — скаляр, поскольку срезы не зависят от х2,... хк, а по х1 проводится интегрирование.

Штраф можно обобщить на случай среза вдоль переменных из нескольких факторов. Введем бинарный вектор в € 6 \ 0, ненулевые компоненты которого будут определять штрафуемые факторы. Возьмем в качестве меры изменчивости все возможные суммы смешанных частных производных порядка 2 в к (по 2 переменные из каждого штрафуемого фактора и в к таких факторов). В примере выше штрафуется изменчивость только по переменным из х1, поэтому в — е\ = (1, 0,... ,0). Можно показать, что для произвольного в штраф ро(1) является квадратичным и представим в виде

Ро(Л = (Л + (1 - б1)Ф?Ф1) ... + (1 - вк)х

К (1) х ф£фк, Л)) = уес(Л)т0(№ + (1 - 0*)Ф£Ф*)тес(.Д),

к= 1

где

• каждая матрица Пк (размера Рк*Рк) — это сумма интегральных норм произведения вторых частных производных функций словаря А к,

• (Л,В) = ~ скалярное произведение тензоров,

• ®к=1 используется для сокращенной записи произведения Кронекера К матриц,

• vec (Л) означает переупорядочивание элементов тензора в вектор-столбец.

Рассмотрим класс штрафных функций, который состоит из взвешенных сумм Ро(1) для всех возможным в € Э \ 0, причем веса 70(А) слагаемых этой суммы параметризованы вектором параметров регуляризации А > 0 длины К:

Ра = {ВД):Л(/)= £ 7в(АЫ/),7в(А)>0}. (2) ее(е\о)

В алгоритме мы воспользуемся штрафной функцией Рд(/), в которой 7о(А) = А^ для случая 0к = 1, = 1 и 7о(А) = 0 для

К d, di

fc=l ti.i^fc 9=1 r=l

Такой выбор 7e (А) дает возможность с помощью вектора регуляризации А явно определять вклад изменчивости по переменным из каждого фактора в общую функцию ошибки R(f,H,g(£)).

Вторая глава посвящена задаче формирования словарей А Первый раздел содержит алгоритм построения аппроксимации для выборок с факторным планом. Рассмотрим этап построения словарей Дk,k = 1,... ,К по заданной обучающей выборке S. Зафиксируем словари к ф I и рассмотрим задачу выбора словаря Д;. Выборочная функция ошибки в этом случае может быть переписана следующим образом:

Q(/,S,5(S)) = trace ((У- ФгВг)Т(Уг) - Ф,В,)),

где обозначает развертку тензора вдоль направления I (переупорядочивание элементов в матрицу размера щ х Ylk^/i пк)- Задача минимизации Q по функциям {т/А = Д; (точнее, по их параметрам) может быть заменена построением аппроксимации по выборке Si = {ai,y^}. Таким образом, словари Д^, к — 1,.... К предлагается строить, решая задачу аппроксимации отдельно в каждом из

факторов. Следует отметить характерную особенность этой задачи: размерность выходных данных (равная пк) может существенно превышать не только размерность входных данных но и объем обучающей выборки щ.

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

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

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

где / (х) задается матрицей 0 параметров функций словаря и вектором а коэффициентов линейного разложения.

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

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

Заметим, что зависимость функции ошибки С} от параметров функций словаря 0 нелинейна и довольно сложна. В то же время, за-

р

б? (©, а) = {у- / (X, 0, а)) Т (у - / (X, ©,«)) + ХаТПа

висимость Q от коэффициентов разложения а квадратична. При фиксированных Э оптимальное значение а может быть найдено явно с как решение задачи обобщенной гребневой регрессии:

а (в) = (Ф(9)тф(0) + А^_1ф(в)ту. (3)

Рассмотрим новую целевую функцию R{Q) = Q (Э, а (Ö)), где а (О) вычисляется по явной формуле (3). В этом случае возникает зависимость а (0), которую нужно учитывать при нахождении производных i?e и 7?ее модифицированной функции ошибки. Обозначим вектор невязок как е = / (X, 0, а) — у, а матрицу производных модели по 0, подсчитанных в точках обучающей выборки, как J.

Предложение 1. Градиент и гессиан функции ошибки R выражаются следующими формулами:

Re = Qe — eTJ;

Ree = JTJ + /ее ®1ет - (Фе ®,ет + ЛТФ) х

х (ФТФ + ЛП)"1 (Фе в.е"1, + ЛТФ)Т ,

где Ф© ^е"1 = ^ZÜi ei"Фе fa).

Таким образом, учтена нелинейная зависимость а (0) и получены явные формулы для первых и вторых производных функции ошибки Ree■ Рассмотрим вместо гессиана Ree некоторое его приближение Н, которое позволяет строить выпуклые локально-квадратичные модели.

Предложение 2. МатрицаЯ =f ЛТЛ-(ЛТФ) (ФТФ + Ail)"1 (ЛТФ)Т неотрицательно определена для любого А > 0.

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

Третья глава посвящена решению задачи аппроксимации для выборок с факторным планом эксперимента. В первом разделе обсуждается задача поиска оптимального решения в выбранном классе фуНКЦИЙ Ftensor с штрафом P\(f) € Рд. Показано, что задача сводится к минимизации квадратичной функции ошибки по элементам тензора коэффициентов разложения Л. Особенностью этой квадратичной функции является ее высокая размерность и неразреженная структура матрицы гессиана. В общем случае для минимизации функций такого рода требуется 0(Р3 + NP2) операций (напомним, что N — число наблюдений, а Р — объем словаря). В некоторых специальных случаях известны более эффективные способы.

В случае полного плана мы предложим явную формулу для оптимального решения, которая в случае Пк > Рк имеет сложность 0(NfuuJ2nk), где пк — число точек в факторе ак, Nfuu = Цпк— мощность полного факторного плана, а рк — объем словаря Ак- В случае неполного плана мы воспользуемся итеративным подходом, который позволит получить сложность 0(min(Nfuu — N + 1, N + 1 ,Р) х Nfvu'Enk).

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

Предположение (Положительной Определенности). Пусть все матрицы и Пк (к = 1,..., К) положительно определены.

Это техническое условие, которое фактически означает, что словари Д^ не вырождены.

Во втором разделе исследуется задача нахождения оптимального решения в случае полного факторного плана. Следующий результат дает достаточное условие на существование и единственность решения.

Предложение 3. Если если выполняется предположение ПО, то решение задачи аппроксимации по выборке с полным факторным планом в классе Ftensor со штрафом P\(f) £ Рд существует и единственно.

Как было указано выше, решение задачи аппроксимации соответствует минимизации квадратичной формы. В случае полного

плана она будет записана как

vec(^)Tíívec(^)+(vec(3^) —Ф vec(.4))T(vec(;y)—Ф vec(.A))j,

где матрица О определяется из (2) по аналогии с (1), то есть является взвешенной суммой произведений Кронекера.

Предлагается следующий способ нахождения решения задачи. С помощью методов линейной алгебры строится преобразование гессиана квадратичной функции ошибки R(A).

Лемма 1. Если справедливо предположение ПО, то гессиан Г2 + ФТФ представим в виде

Г2 + ФТФ = VDVT,

гдеУ = Vk, aD — положительно определенная диагональная

матрица.

Следующая теорема дает оптимальное решение задачи аппроксимации для полного факторного плана в явном виде.

Теорема 1. Если справедливо предположение ПО, то оптимальный тензор коэффициентов разложения по словарю Atensor вычисляется по формуле:

Л* = ((3> ... ®кФкУ^т) О р-1) ^Vf1... ®ку-!. (4)

где V определяется с помощью равенства vec (Т>) = diag(D), а символ О подразумевает поэлементное умножение тензоров или матриц.

Предложение 4. Сложность вычисления оптимального тензора коэффициентов разложения Л по формуле (4) равна

(к К к

nkPk + NJ2 nkJ[pi/ni к=1 к=1 2 = 1

Введем обозначение Т = N пк Hf=i Pi/ni Для второго сла-

гаемого, которое в большинстве случаев и вносит основной вклад

R(A) =

в вычислительную сложность. Т операций необходимо для последовательного умножении тензора У на матрицы Ф^У^Г Т вдоль всех направлений.

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

Пример 1. Рассмотрим случай интерполирующего словаря и плана эксперимента с одинаковым числом уровней во всех факторах (pk — Щ — N1/,K = Pl/K). Тогда сложность (4) составит 0(PY,Pk) — 0{Рх+1/к) вместо 0(Р3) для стандартного подхода. Для типичной ситуации п = 10, К = б эти оценки составят 107 и 1018 соответственно.

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

Y,fuu\H. Это позволит нам переписать расширенное множество значений в виде тензора. Обозначим этот тензор как У, а расширенную матрицу регрессоров, посчитанную во всех точках из Е¡иц. — как Ф. Теперь введем модифицированную взвешенную функцию ошибки R(A):

R{A) = уес(Д)тПуес(Л)+(уес(У)-Ф vec(.4))T W (тсс(5>)-Ф vec(.4)).

где W — диагональная матрица весов, которая формируется следующим образом:

• WM = 1, если хг С Е (г-я точка полного плана Е¡иц представлена в обучающей выборке);

• Witi = 0, если Xi Е (S/uii \ Е) (г-я точка была пропущена).

При таком выборе весовой матрицы W справедливо R(A) = R(A) при любом способе заполнения пропущенных значений в тензоре У. Таким образом, задача минимизации R(A) может быть заменена эквивалентной задачей минимизации функции R(A), которая также является квадратичной по коэффициентам разложения А.

Достаточное условие на существование и единственность минимума функции Я (Л) эквивалентно случаю полного факторного плана.

Предложение 5. Если справедливо предположение ПО, то решение задачи аппроксимации по выборке с неполным факторным планом в классе Ftensor со штрафом 6 Рд существует и единственно.

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

Лемма 2. Гессиан Нд функции fí.(Á) может быть представлен в виде

НА = V (б + UTWU) VT,

где D — положительно определенная диагональная матрица, U раскладывается в произведение Кронекера ортогональных матриц размеров rik * Pk U = Uk, а V = Vk-

Лемма 3. Пусть Л = {Л ... © Т>а ъ, где тензор V за-

дается с помощью равенства vcc = diag(D).

Определим квадратичную функцию R(A) с помощью соотношения R(A) = В,(Л). Тогда гессиан функции R(A) имеет не более чем N + 1 различных собственных чисел.

Лемма 4. Пусть А = (A ®,Vi... ®KVx) О V°-5.

Определим квадратичную функцию R(A) с помощью соотношения R(A) = R(A). Тогда гессиан функции R(A) будет иметь не более чем Njuu — N + 1 различных собственных чисел.

Леммы 2, 3 и 4 позволяют найти оптимальное решение.

Теорема 2. Если справедливо предположение ПО, то решение задачи минимизации R(A) будет найдено с помощью метода сопря-градиентов

• не более чем за mm(Nfuu — N + 1, N + 1) шагов;

• вычислительная сложность каждого шага не будет превосходить 0(Т), гдеТ = N nk Пг=1 Pi/nH

• сложность дополнительных вычислений, которые выполняются однократно, составит 0(Т + +Pfc))-

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

Следствие. Общая вычислительная сложность алгоритма поиска оптимальных коэффициентов разложения Л составляет

О (j2Pk(nk +Pk) + Tx min (NfuU -N + 1,N+1,P)).

Очевидно, что основным слагаемым в оценке сложности является 0(Т х min(Nfuu — N + 1, N + 1, Р)). Напомним, что основная часть сложности вычисления оптимального решения для полного факторного плана составляет 0(Т). Таким образом, сложность возрастает почти "непрерывно" при переходе от полного плана к неполному и затем увеличивается пропорционально количеству пропущенных точек.

Заметим, что любое множество £ может быть представлено в виде неполного факторного плана (с огромным числом пропущенных точек). Естественно, что для произвольной выборки предложенный подход не может превосходить по эффективности стандартную формулу. При постановке задачи аппроксимации для неполного факторного плана мы ограничились формулировкой "£ составляет достаточно большое подмножество Получив точные оценки для вычислительной сложности формализуем понятие "достаточно большое подмножество". Пусть а — это доля точек обучающей выборки в полном факторном плане, N = aNfuu. Рассмотрим распространенный вариант формирования словаря, при котором Р « N, a pk ~ пк.

Предложение 6. Пусть объем словаря равен числу точек выборки Р — N — aNfuii. Тогда отношение сложности предложенного

2 1 О -1

а. % -2 О)

5 -з

-4 -5

-6

- *« ■ Оценка сверху

- ■ - ■ Точное значение

-1.5 юд,„с

Рис. 1: Отношение вычислительных сложностей (СИ.) предложенного подхода и явной формулы

метода к сложности стандартного подхода будет равно

св = тт(а,1 - а) ^=1 пка1/к < Ек=1 пк ПГ=1 Пк О? Цк=1

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

Пример 2. Рассмотрим случай Пк = п, подставим константы из примера 1 (п = 10 ,К = 6) и построим график зависимости отношения сложностей от а (см. рисунок ) и ее верхней оценки, приведенной выше. Приведем некоторые характерные значения, полученные из графика:

• критическая доля точек равна асг — 0.0024;

• при а = 0.1 (90% точек пропущено) сложность предложенного метода в 520 раз меньше сложности стандартного;

• при а = 0.9 (10% точек пропущено) сложность предложенного метода в 129 000 раз меньше сложности стандартного.

Если же выбирать объем словаря Р = Nfuu (что необходимо в ситуации, когда S распределено существенно неравномерно в то предложенный способ будет эффективнее стандартного способа вычисления независимо от доли точек а.

Четвертый раздел посвящен выбору параметров регуляризации Л. Для полного факторного плана описывается вычислительно эффективный алгоритм выбора вектора Л, основанный на минимизации ошибки скользящего контроля10.

Предложение 7. Ошибка скользящего контроля в случае полного факторного плана выражается в виде

Qloo - (У -У, (У-У)®£2),

-1 (5)

С= (т-2>-1®1(и?,ои?')...®к(и£ ©и£)) ;

где У — Л ^Ф^ ... — значения модели f(A) в точках плана

^full; а I — тензор размера п\ х • ■ • х пк-

Сложность формулы (5) равна 0(Т),Т = Nfuu J2k= i nk П/=1 Pl/nl)•

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

Qloo = ({У - У) о W, (у - У) о С2).

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

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

10Fukunaga К., Hümmels D. М. Leave-one-out procedures for nonparametric error estimates // Pattern Analysis and Machine Intelligence, IEEE Transactions on. 1989. V. 11, no. 4. P. 421-423.

источников данных. В разделе показано, что решение этой задачи может быть сведено к задаче построения аппроксимации по неполному факторному плану, но с использованием весов W,^, отличных от {0,1},

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

Публикации по теме диссертации

1. Беляев М.Г. Аппроксимация данных, порожденных декартовым произведением // Труды МФТИ. 2013. Т. 5, № 3. С. 11-23.

2. Беляев М.Г. Аппроксимация многомерных зависимостей по структурированным выборкам // Искусственный интеллект и принятие решений. 2013. № 3. С. 24-39.

3. Belyaev М., Burnaev Е. Approximation of a multidimensional dependency based on a linear expansion in a dictionary of parametric functions // Informatics and its Applications. 2013. Vol. 7, no. 3. P. 114-125.

4. Беляев М.Г. Аппроксимация и интерполяция на. основе тензорного произведения параметрических словарей // Труды конференции Информационные Технологии и Системы. 2012. С. 32-40.

5. Surrogate Modeling of Stability Constraints for Optimization of Composite Structures / S. Grihon, E. Burnaev, M. Belyaev et al. // Surrogate-Based Modeling and Optimization. Engineering applications / Ed. by S. Koziel, L. Leifsson. 2013. P. 359-391.

6. Беляев М.Г. Аппроксимация зашумленных данных, имеющих структуру декартова произведения // Сборник докладов конференции ИОИ-9. 2012. С. 188-192.

Беляев Михаил Геннадьевич

Построение аппроксимации многомерных зависимостей по выборкам с факторным планом

АВТОРЕФЕРАТ

Подписано в печать 14.11.2013. Формат 60 х 84 1/16- печ- л- 1Д Тираж 100 экз. Заказ 360.

Федеральное государственное автономное образовательное учреждение высшего профессионального образования ¡Московский физико-технический институт (государственный

университет)»

Отдел оперативной полиграфии «Физтех-полиграф» 141700, Московская обл., г. Долгопрудный, Институтский пер., 9.