УДК 519.226
Е. В. Бурнаев, П. Д. Ерофеев, П. В. Приходъко
Московский физико-технический институт (государственный университет) Институт проблем передачи информации им. A.A. Харкевича ООО «ДАТАДВАНС»
Выделение главных направлений в задаче аппроксимации на основе гауссовских процессов
Для решения множества прикладных задач необходимо уметь восстанавливать неизвестную зависимость по выборке ее значений. Одним из наиболее эффективных методов восстановления зависимостей является метод построения аппроксимации на основе гауссовских процессов. Зачастую ковариационную функцию гауссовского процесса моделируют с помощью гауссовской ковариационной функции на основе взвешенного евклидова расстояния. Однако такое априорное предположение о ковариационной структуре источника данных означает, что направления наиболее значимого изменения аппроксимируемой функции и направления координатных осей в пространстве данных совпадают, что редко выполняется на практике. В работе предложен эффективный способ моделирования ковариационной структуры данных для более общего случая. Разработан эффективный алгоритм настройки параметров ковариационной функции. Показано, что разработанный алгоритм также позволяет решать задачу эффективного снижения размерности. Предложен статистический тест, позволяющий оценить сниженную размерность данных. Все основные утверждения подкреплены экспериментальными результатами на тестовых задачах. Общая эффективность разработанного подхода продемонстрирована на примере решения реальной задачи из области самолетостроения.
Ключевые слова: эффективное снижение размерности, гауссовские процессы, расстояние Махалонобиса.
1. Введение
Для решения многих задач инженерного проектирования в различных областях промышленности необходимо уметь строить высокоточные аппроксимации многомерных зависимостей на основе данных, собранных в результате проведения натурных экспериментов и/или симуляционных экспериментов с физическими моделями.
В настоящее время для восстановления зависимостей широко распространено применение методов, основанных на гауссовских процессах. Этот подход предполагает, что неизвестная зависимость порождена гауссовским полем с неизвестной ковариационной функцией, принадлежащей некоторому параметрическому семейству. В рамках такой модели задача построения аппроксимации сводится к оценке неизвестных параметров ковариационной функции по имеющемуся набору данных с помощью максимизации соответствующей функции правдоподобия.
Стандартным параметрическим семейством, используемым для моделирования ковариационной функции, является гауссовская ковариационная функция (на основе взвешенного евклидова расстояния). Использование такого семейства фактически неявно означает, что верно предположение о совпадении направлений наиболее значимого изменения аппроксимируемой функции и направлений координатных осей в пространстве данных. Очевидно, что это предположение является достаточно ограничивающим и не всегда выполняется в реальных задачах. Для того чтобы снять указанное ограничение и сделать модель данных более общей, было предложено использовать ковариационную функцию на основе расстояния Махалонобиса [7]. Однако такой подход не позволяет работать с данными высокой размерности («проклятье высокой размерности», см. [4]) в связи с существенно большим
числом параметров ковариационной функции, которые необходимо оценивать. Таким образом, актуальна задача разработки алгоритма настройки параметров ковариационной функции на основе расстояния Махаланобиса, который бы позволял проводить эффективную настройку параметров и в случае данных высокой размерности.
В данной работе предложен новый алгоритм подстройки параметров ковариационной функции на основе расстояния Махаланобиса, позволяющий эффективно решить поставленную выше задачу. Кроме того, разработанный подход позволяет линейно снизить размерность данных в некоторых задачах, выделив главные направления изменения аппроксимируемой зависимости.
Работа устроена следующим образом. В разделе 2 приведены основные обозначения и кратко изложена суть метода аппроксимации на основе гауссовских процессов. В разделе 3 предложен новый метод оценки параметров ковариационной функции гауссовского процесса на основе расстояния Махаланобиса. В разделе 4 приведен подход к решению задачи эффективного снижения размерности на основе предложенного метода выделения главных направлений изменчивости аппроксимируемой зависимости, описан статистический тест для оценки размерности сжатия. В разделе 5 приведены результаты вычислительных экспериментов на искусственных и реальных данных.
2. Регрессия на основе гауссовских процессов
Задача аппроксимации многомерной зависимости ставится следующим образом. Пусть у = f (х) — некоторая неизвестная функция со входом х £ X С Мт и выходом у € М1, Огеат = = {(хг, у^ = /(х^)) ,г = 1, Ж} — обучающая выборка. Задача состоит в
построении аппроксимации у = /(х) = /(х|^еагга) для исходной зависимости у = /(х) по обучающей выборке В1еагп.
Если для всех х € X (не только для х € 01еагп) имеет место примерное равенство
/(х) * /(х), (1)
то считается, что аппроксимация хорошо воспроизводит исходную зависимость. Выполнение соотношения (1) провкряют на независимых тестовых данных = = {(х^, у^- = £(х^-)),] = 1,Х}, оценивая качество построенной
аппроксимации с помощью среднеквадратичной ошибки:
е(/|А„1) =
\
Т м*
тт £(у - / (х; ))2. (2)
.7 = 1
2.1. Гауссовские процессы
Гауссовский процесс является одним из возможных способов задания распределения на пространстве функций. Гауссовский процесс / (х) полностью определяется своей функцией среднего т(х) = Е[/(х)] и ковариационной функцией соь(у, у') = к(х, х') = Е[(/(х) — т(х))(/(х') — ш(х'))]. Если положить функцию среднего нулевой т(х) = Е[/(х)] = 0, а ковариационную функцию считать известной, то функция апостериорного (для заданной обучающей выборки) среднего значения гауссовского процесса в точках контрольной Выборга X имеет вид [6] /(X) = К*К-1У, где Кч = К(х, х) = [к(ъ, ху), г = 1ГЖ, 3 = 17х], К = К(X, X) = [к(хг, х), г, ] = Т^].
Обычно предполагается, что данные наблюдаются с шумом: у(х) = /(х) + е(х), где е(х) ~ М(0, <г2) — белый шум. В таком случае наблюдения у(х) являются реализацией гауссовского процесса с нулевым средним и ковариационной функцией соу(у(х), у(х')) = к(х, х') + а25(х — х'), где 5(х) — дельта-функция. Таким образом, функция апостериорного (для заданной обучающей выборки) среднего значения гауссовского
процесса f (x) в точках контрольной выборки X * принимает вид
f (Х*) = К* (К + a2l) lY,
(3)
Где / — единичная матрица размера N х N.
Заметим, что наличие в формуле (3) дисперсии шума а2 фактически приводит к регуляризации, что позволяет улучшить обобщающую способность аппроксиматора. При этом апостериорная ковариационная функция гауссовского процесса в точках контрольной выборки имеет вид
У[Х*] = К(X*,Х*) + а21* — К*(К + а21)~ 1К^, (4)
где К(Х*,Х*) = [к(х.г, х), г,] = 1,..., Ж*], I* — единичная матрица размера Ы* х Ы*.
Дисперсии гауссовского процесса в точках контрольной выборки могут быть использованы как оценки ожидаемой ошибки аппроксимации в этих точках. Заметим, что для этого нет необходимости вычислять по формуле (4) всю матрицу У [X*], а достаточно вычислить только элементы ее главной диагонали, которые и являются искомыми дисперсиями.
Кроме того, зная среднее и ковариационную функцию, можно также получить апостериорную оценку среднего и дисперсии производной гауссовского процесса в точках. Если
9/(х)
9(xo) =
dx
x=xo
ТО
Law(g(xo)l(X,Y)) = N(JT(К + a2l)-1Y, В - JT(К + a2!)-1 J), где
J
т
dk(xo - xi) dk(xo - xn)
d xo
dxo
В =
cov(gi(xo),gi(xo)) . . . cov(gi(xo),gm(xo))
cov(g (xo),gi(xo)) .
cov(di,dj) = i-я компонента вектора градиента g.
cov(gm(xo),gm(xo))_
d2k(xo, xo) dxldxi
2.2. Оценка параметров гауссовского процесса
При работе с реальными данными ковариационная функция породившего их гауссовского процесса, как правило, неизвестна, поэтому необходимо уметь идентифицировать ее по данным.
Предположим, что ковариационная функция гауссовского процесса принадлежит некоторому параметрическому семейству к(х, х') = к(х, х'|а), где а € М^ — вектор параметров ковариационной функции. Семейство к(х, х'|а) обычно берется из класса так называемых стационарных ковариационных функций, т.е. функций, значения которых зависят только от разности значений аргументов к(х, х'|а) = к(х — х'|а). Стандартный способ оценки значений параметров а по обучающей выборке 01еагп состоит в использовании принципа максимума правдоподобия. Логарифм правдоподобия гауссовского процесса в точках обучающей выборки [6] при этом имеет вид
log p(YIX, а, а) = -2YT(К + a2I)-iY - 1 log + а211 - у log 2ж,
(5)
где IK + а211 — детерминант матрицы К + а21.
Кроме параметров а ковариационной функции, параметром функционала (5) является также значение дисперсии шума наблюдений а2, которое также можно настраивать по обучающей выборке. Таким образом, нахождение оптимальных значений параметров сводится к отысканию максимума правдоподобия по параметрам log p(YIX, а, а) ^ maxa,a-
9
2.3. Ковариационная функция
Выбор конкретного семейства ковариационных функций к(х, х'|а) обычно продиктован соображениями удобства, а также априорными представлениями о свойствах аппроксимируемой зависимости.
Для моделирования ковариационной функции часто используется экспоненциальное семейство Ко(х, х') = а¡2 ехр(—((х, х')), где (1(х, х') — расстояние между точками х и х' в некоторой метрике. Наиболее распространены два типа расстояния:
• Взвешенное евклидово расстояние
т
((х, х') = ^ е2(х* — (ж'Л2, х =(х1, ■■■ ,Хт), (6)
г=1
где вг € М,г = 1,т. При использовании такого расстояния возникает необходимость настроить т параметров в1,..., вт.
• Расстояние Махалонобиса
((х, х') = (х — х' )ТА(х — х'), (7)
где А € Мтхт — некоторая положительно определенная матрица. Для параметри-
А
разложение Холецкого А = ЬТ Ь, где Ь — верхняя треугольная матрица с положительными элементами на главной диагонали. В таком случае необходимо настроить
т(т+1) Т
—^2—- параметров, задающих матрицу Ь.
В работе [7] показано, что использование ковариационной функции на основе расстояния Махаланобиса обеспечивает инвариантность точности аппроксимации по отношению к степени коллинеарности направлений наиболее значимого изменения аппроксимируемой функции и направлений координатных осей в пространстве данных. В то же время при использовании ковариационной функции на основе взвешенного евклидова расстояния получаются аппроксимации, точность которых существенно зависит от взаимного расположения направлений наибольшей изменчивости аппроксимируемой зависимости и координатных осей в пространстве данных. Однако качество построенных моделей, испольующих ковариационную функцию на основе расстояния Махаланобиса, существенно ухудшается с ростом размерности. Это связано с квадратичным ростом числа настраиваемых параметров при увеличении размерности задачи.
3. Выделение значимых направлений через аппроксимацию градиента
В этом разделе описан подход к оценке параметров ковариационной функции на основе расстояния Махаланобиса, который позволяет проводить эффективную настройку параметров и в случае данных высокой размерности.
При определении ковариационной функции типа (7) на основе расстояния Махаланобиса можно использовать следующую параметризацию1:
((х, х') = (х — х')тВт ЛВ(х — х'),
где В — ортогональная матрица т х т, Л — диагональная положительно определенная матрица размера т х т. При такой параметризации матрица В осуществляет поворот осей
Л
осей.
1 Матрица А из расстояния Махаланобиса (7) может быть приведена к такому виду через сингулярное разложение.
В
/(х) = д(хВ), где д — аппроксимация па основе гауссовских процессов, использующая функцию ковариации на основе взвешенного евклидова расстояния (6). В таком случае
Л
онной функции соответствующего гауссовского процесса.
Л В
Г
■ „ „ N
df(x) df(x)
dx1
dxr-
x=x;
/ i=\
x = (x\...,xm).
В таком случае мы можем развернуть координатные оси так, чтобы производные вдоль
В
Гт Г.
В Л
шагов:
learn■
• Строим аппроксимацию f на основе гауссовских процессов по выборке Die
• Вычисляем градиенты аппроксимации Г = {al(x') ,..., 1 в точках
I ax X=Xi ax x=Xi) i=1
обучающей выборки.
• Вычисляем матрицу ковариации градиентов Х^ = ГТГ.
• Вычисляем матрицу В собственных век торов Хр.
• «Поворачиваем» входные векторы выборки Diearn, получаем новую выборку Dleam = (Z, Y) = {(z, = ХгВТ, уг ),i =T,N }. * * *
• Настраиваем веса во взвешенном евклидовом расстоянии с помощью максимизации правдоподобия выборки Diearn для соответствующего гауссовского процесса. Полученные значения весов используются в качестве оценки матрицы Л.
Возможно дополнительно улучшить качество итоговой аппроксимации через итератив-В
• Зададим t = 1, x(0) = x, (i = T,..., N), Во = I, где I — единичная матрица.
• По выборке данных |x(i 1), y }N=1 оцениваем матрицу проектирования В, используя базовый алгоритм.
• Пересчитываем итоговую матрицу проектирования согласно формуле Bt = В ■ Bt-1.
(t) Т
• Проектируем xj из текущего базиса в базис, заданный В¿: x^ = х^ВТ ,i = 1,..., N.
• Если мы достигли последней итерации, то останавливаемся; иначе t = t + 1, и переходим к шагу 2.
В таком случае итоговый поворот можно записать как В = ПгВг. Построение аппроксимации на основе новой выборки Diearn = (Z,Y) = {(zj = x^Т, yi) ,i = 1,N} с «повернутыми» векторами входов зачастую существенно увеличивает точность аппроксимации.
4. Эффективное снижение размерности на основе гауссовских процессов
x=x
В этом разделе приведена постановка задачи эффективного снижения размерности, а также показано, что алгоритм, предложенный в предыдущем разделе, позволяет ее решить.
4.1. Задача эффективного снижения размерности
Пусть процесс, порождающий данные, имеет вид у = /(х) + е(х), где предполагается, что /(х) = д(хВт), £ — случайный шум с Ее = 0 и В € МЛхт, й <т, ВВТ = 1аха-
Под задачей эффективного снижения размерности [19] понимают задачу нахождения выпуклой оболочки 5 = ,врап{В}, которая задает центральное среднее подпространство регрессии.
Пусть В — оценка матрицы В, тогда должно быть выполнено следующее равенство:
/(хВтВ) « /(х). (8)
Если В содержит в себе центральное среднее подпространство, т.е. В € врап(В), то равенство (8) выполняется в точности.
Пусть {Д, р2,..., рт} - собственные векторы, а {А1 > Х2 > ■ ■ ■ > Хт} — соответствующие им собственные числа матрицы Х^ = ГТГ, полученные в алгоритме поворота координатных осей, описанном в разделе 3.
Если для некоторого к выполнено, что Хг ~ 0,1 > к, то рассматриваемая аппроксимация имеет почти нулевую производную вдоль направлений > к, т.е. значение функции /(х) слабо меняется вдоль направлений > к. Значит, возможно снизить размерность входных данных, спроектировав данные на линейное подпространство В, составленное из собственных векторов матрицы Х^, соответствующих наибольшим собственным числам. При этом «потеря информации» об аппроксимируемой функции будет минимальной.
4.2. Оценка эффективной размерности данных
В случае применения предложенного подхода эффективная размерность данных (то есть размерность центрального среднего подпространства) определяется рангом матри-
Г
процедуры снижения размерности, сформулированной в предыдущем разделе, сводится к определению размерности линейного подпространства, натянутого на столбцы матрицы градиентов.
В литературе предложено большое количество подходов к оценке эффективной размерности данных для рассматриваемой постановки задачи.
При использовании метода главных компонент часто в качестве оценки требуемой размерности используют число компонент, объясняющих заданную долю дисперсии данных [5]. Оценка доли объясненной дисперсии в свою очередь равна доле суммы наибольших
у^ д.
к собственных чисел, т.е. 1к = Чт-т1 Д. Достаточное значение 1к обычно задается внешним
2^г= 1 д
для процедуры образом и никак не учитывает особенностей конкретных данных. На практике часто выбирают значения между 70% и 99%.
В [14] описывается тест на сферичность. Тест основан на том, что если вся значимая дисперсия объяснена выбранными главными компонентами, то в дополнении к подпространству, задаваемому матрицей В, дисперсия изотропна и разброс данных в этом подпространстве может быть представлен многомерной сферой. При заданном к проводится тест на равенство между собой собственных значений.
Близкий подход предлагается в [20]. Здесь используется то обстоятельство, что если присутствуют два одинаковых собственных числа, то соответствующее им двумерное подпространство будет устойчиво к возмущениям выборки, однако главные направления в нем будут нестабильны. С помощью бутстрепа можно сгенерировать набор случайных реализаций выборки. Размерность многомерной сферы, содержащей оставшуюся дисперсию, может быть найдена из тех соображений, что оценка этого подпространства по различным бутстреп-реализациям будет наиболее стабильной.
В [10,18] предлагается универсальный подход, позволяющий оценить ранг матрицы В линейного проектирования, если ее оценка В оказывается асимптотически нормально рас-
пределенной и несмещенной (например, в работе приведен явный вид теста на размерность данных для метода эффективного снижения размерности SIR [8]).
Далее опишем статистический тест для предложенного подхода к эффективному снижению размерности, позволяющий оценить размерность искомого подпространства.
Так как данные являются некоторой реализацией нормального поля, то для оценки ранга можно применить общие результаты из [10,18], которые позволяют проверять гипотезу о ранге случайной матрицы.
Пусть параметры случайного гауссовского поля заданы, то есть известно совместное распределение Y ~ М(ц, K), Y = (y1,..., yN)• В таком случае матрица градиентов апостериорных средних рассматриваемой зависимости /(x) имеет безусловное распределение:
Г
= г Mx)
I dx1
df (x)
dxr
x=xl J г=1
1 N
( JK-V, £),
=1
где £ = J1 K-1 J, матрица J имеет вид J = [JT , JT,..., J^ ]T , J¿ - матрица размера m xN вида Ji = ^É^
. к(х) = {К(х, х
Выпишем сингулярное разложение матрицы Г = иИУ, где и — ортогональная матрица размера N х^ V — ортогональная матрица размера т х т, И — диагональная матрица размера N х т.
Если исходная матрица В на самом деде имеет ранг к, то можно представить матрицу И в виде
"А 0 0 И2
D
где Di — невырожденная диагональная матрица размера к x к, D2 (N — к) x (m — к), элементы которой должны равняться нулю. Аналогично U =
матрица размера U1
U2
, где
Ui и U2 имеют размер к х N и к х (N — к), V = [Vi V2], матрицы Vi и V2 имеют размер т х к и т х (т — к).
Используя результаты [10,18], можно показать, что если I = vес(UTГ V2), вде vес() — вектор, составленный из элементов соответствующей матрицы, и Q = (V2T®UТ)Х(^ ),
где ® — произведение Кронекера, то lTQl ~ X2mm((m-k)(N-k),rank(s))-
Q
диент не во всех N точках, а достаточно рассмотреть только N % т < N * < N точек. В
2
таком случае процедура оценки сводится к последовательной проверке критерия % для величины lTQl и различных размерностей к начиная с единичной.
5. Экспериментальные результаты 5.1. Устойчивость к повороту осей
Продемонстрируем предложенный подход на примере функции Mystery (рис. 1):
у = 2 + 0.01(х2 — х1)2 + (1 — xi)2 + 2(2 — х2)2 + 7sin(0.5xi) sin(0.7xix2),
где x1 е [0, x2 G [0, 5].
Пусть задана некоторая случайная выборка данных Diearn = |x¿, y¿}8:= ная рассматриваемой функций. Будем строить аппроксимации по этой выборке в виде У = /(x VT), где
У = cos (f — sin ф sin (f cos (f
x=x
есть матрица поворота осей координат на угол
Рис. 1. Функция Mystery-
Рис. 2. Чувствительность аппроксимации к повороту координатных осей
.0 0.5 1.0 1.5 2.0 2.5 3.0
Угол поворота координатных осей
— VEGA
Взвешенное Евклидово расстояни Расстояние Махаланобиса
На рис. 2 изображена зависимость качества модели от угла поворота осей координат <р для различных способов задания расстояния в ковариационной функции2. Под качеством модели понимается среднеквадратичная ошибка аппроксимации на тестовой выборке (2).
Оценка гиперпараметров аппроксимации, использующей ковариационную функцию на основе расстояния Махаланобиса, с помощью стандартного и предложенного подходов позволяет добиться независимости качества модели от поворота осей координат. В приведенном двумерном случае стандартный и предложенный в данной работе алгоритмы настройки гиперпараметров дают одинаковую точность, однако, с ростом размерности стандартный метод обучения модели с ковариационной функции, использующей расстояние Махаланобиса, существенно проигрывает в точности по сравнению с предложенным подходом поворота координатных осей в связи с ростом числа гиперпараметров, которые необходимо оценивать в процессе обучения (см. также результаты тестирования на многомерных данных в [1|).
5.2. Сравнение точности методов эффективного снижения размерности
В данном разделе протестирована точность предложенного метода эффективного снижения размерности в том случае, когда зависимость, порождающая данные, имеет вид
f(x)=g(xBT). (9)
Для иллюстрации качества решения задачи эффективного снижения размерности с помощью разработанного алгоритма было проведено сравнение предложенной процедуры с популярными подходами: Sliced Inverse Regression [8,9], Structural Adaptation via Maximum Minimization [19], Minimum Average Variance Estimation [12], Outer Product of Gradients [12], Partial Lest Squares [11].
Точность процедуры определялась близостью пространств, задаваемых настоящей матрицей В и ее оценкой В. Мера близости вычислялась согласно формуле Е = (I(/ — ВтВ)ВТВК. Размерность сжатия при этом считалась известной.
Для тестирования использовался набор из 27 гладких функций, которые обычно применяются для тестирования алгоритмов оптимизации [21]. Размер обучающей выборки составлял 300 точек, исходное пространство имело размерность 10, матрица В размер 5 х 10. Для осуществления сжатия генерировалась ортогональная матрица Bfuu размера 10 х 10. Из ее 5 строк составлялась матрица В. Далее для функции / в точке х £ R10 значение вычислялось согласно формуле у = f(x.BT).
"Предложенный метод здесь и далее обозначен как VEGA Variable Extraction via Gradient. Approximation.
Рис. 3. Кривые Долана Мора для набора тестовых функций по ошибкам сжатия/восстановления
Результаты сравнения (с помощью кривых Долана Мора, см. [3]) методов приведены на рис. 3. Предложенная процедура оказалась наиболее эффективной в сравнении со стандартными подходами.
5.3. Эксперименты по оценке размерности
В этом разделе рассмотрена эффективность предложенного метода оценки размерности в случае выполнения модели (9). Далее приведены результаты вычислительных экспериментов для случая пространства размерности 10 (x £ [-1, 1]10) и случайной ортогональной матрицы В = , } размер а 5 х 10. В экспериментах использовались следующие искусственные функции:
• Сумма квадратов (SqSum): /1(x) = ^®=1 i(fyx)2.
• Функция Griewank: Д(х) = Ei=1 ТОО^ - П®=1 cos (^^/р) + 1.
• функция Trid: Д(х) = ^5=1 (^¿х)2 - ^5=1(^ix)(^i-1x).
• функция Perm: /1 (x) = EiU^ + (^ - 1) )2 .
В таблице 1 приведены оценки размерности для рассмотренных функций в зависимости от размера обучающей выборки.
Т а б л и ц а 1
Результаты оценки размерности с помощью предложенного
статистического теста
Function Sample size
50 60 70 80 90 100 200 300 400
SqSum 2 2 7 5 5 5 5 5 5
Griewank 1 7 7 2 3 5 5 2 5
Trid 1 4 8 1 6 5 5 5 5
Perm 1 1 7 7 7 6 6 5 5
С ростом выборки оценка размерности с помощью предложенного статистического теста постепенно сходится к настоящему значению. Теми сходимости определяется сложностью функции.
Для функции Trid приведен график, показывающий разброс результатов оценки от размера выборки (см. рис. 4). Для каждого размера выборки, рассмотренного в предыдущем тесте (см. таблицу 1), было сгенерировано по 5 случайных наборов данных. По оси ординат отложена разность оцененной размерности и настоящей, то есть значение 0 соответствует правильной оценке размерности. На этом примере видно, что разброс по запускам уменьшается с ростом выборки и оценка сходится к правильному значению. Аналогичный результат для более сложной функции Perm изображен на рис. 5. В этом случае для рассмотренных размеров выборки не удается получить такую же хорошую картину сходимости, однако в целом с ростом выборки оценка также постепенно сходится к настоящему значению.
6
Л
- Медиана оценки Настоящая размерность Диапазон между 5% и 95% квантилями оценки
7
Я
о
400 600
Размер выборки
400 600
Размер выборки
Рис. 4. Оценка размерности на функции Trid для Рис. 5. Оценка размерности на функции Perm нескольких запусков для нескольких запусков
5.4. Эксперименты на реальных данных
0.6
0.5
0.4
0.3
Рис. 6. Ошибка аппроксимации для задачи предсказания энергозатратности маневров вертолета
Для демонстрации эффективности предложенного подхода была рассмотрена задача аппроксимации энергозатрат по параметрам маневра вертолета, включающим в себя 51 величину: масса, моменты вращения, углы, высота и т.п. Были рассмотрены аппроксимации на основе гауееовеких процессов с тремя подходами к подстройке параметров ковариационной функции: 1) стандартный подход с использованием взвешенного евклидова расстояния; 2) предложенный в статье метод на основе расстояния Махаланобиса; 3) пред-
ложенный в статье метод в комбинации со снижением размерности на основе выделения направлений с наиболее значимой изменчивостью аппроксимируемой функции. Ошибка аппроксимации оценивалась во всех случаях с помощью десятикратного скользящего контроля. Эффективная размерность данных, оцененная с помощью статистического теста, составила от 12 до 15 параметров в зависимости от номера итерации процесса скользящего контроля. Результаты, приведенные на рис. 6, говорят о целесообразности применения предложенного метода на данной задаче. Более того, сокращение числа параметров в ~ 4 раза на основе выделения наиболее значимых направлений позволяет сохранить качество аппроксимационной модели, существенно сократив сложность задачи.
6. Вывод
В работе описывается процедура построения регрессии на основе гауссовских процессов, позволяющая осуществить линейное снижение размерности данных. Кроме того, предлагается процедура оценки размерности сжатия данных. Эксперименты показали высокое качество предложенных процедур.
Работа выполнена при поддержке Лаборатории структурных методов анализа данных в предсказательном моделировании МФТИ, грант правительства РФ дог. 11.G34.31.0073.
Литература
1. Burnaev Е., Erofeev P., Prikhodko P. Extraction of the most sensitive directions via Gaussian process modelling // Uncertainty in Computer Models Conference. — 2012. — Sheffield, UK.
2. Ferrari S. and Stengel R.F. Smooth function approximation using neural networks // IEEE Transaction on neural network. - 2005. - V. 16. - P. 24-38.
3. Dolan E., More J. Benchmarking optimization software with performance profiles // Mathematical Programmming. Ser. A 91. — 2002. — P. 201-213.
4. Donoho D.L. Aide-Memoire. High-Dimensional Data Analysis: The Curses and Blessings of Dimensionality. — 2000.
5. Hastie Т., Tibshirani R., Friedman J. The elements of statistical learning: data mining, inference, and prediction. — Berlin: Springer, 2008. - 763 p.
6. Rasmussen C.E., Williams C.K.I. Gaussian Processes for Machine Learning. — MIT Press, 2006.
7. Vivarelli F. and Williams С. K. Discovering hidden features with Gaussian processes regression // Advances in Neural Information Processing Systems 11, ed. by Kearns M. S. and Solla S. A., and Cohn D. A. - MIT Press, 1999.
8. Li Ker-Chau. Sliced Inverse Regression for Dimension Reduction // Journal of the American Statistical Association. — 1991. — V. 86, N 414. — P. 316-327. Published by: American Statistical Association.
http: / / www.jstor.org/stable/2290563.
9. Yehua Li and Tailen Using. Deciding the dimension of effective dimension reduction space for functional and high-dimensional data // Ann. Statist. — 2010. — V. 38, N 5. — P. 30283062.
10. Burn E., Yang J. Dimension Estimation in Sufficient Dimension Reduction: A Unifying Approach // Journal of Multivariate Analysis. — 2011.
11. Tenenhaus M., Vinzia V.E., Chatelinc Y-M., Laurob C. PLS path modeling, Computational Statistics and Data Analysis. — 2005. — V. 48. — P. 159 205.
12. Yingcun Xia, Howell Tong, W. К. Li, and Li-Xing Zhu. An adaptive estimation of dimension reduction space // Journal of the Royal Statistical Society: Series В (Statistical Methodology). - 2002. - V. 64(3). - P. 363-410.
13. Louis Ferre. Determining the dimension in Sliced Inverse Regression and related methods // Journal of the American Statistical Association. - 1997.
14. Peres-Neto P., Jackson D., Somers K. How many principal components stopping rules for determining the number of non-trivial axes revisited. Computational Statistics and Data Analysis 49. - 2005. - P. 971 997.
15. Ruben Daniel and Pedro Valero-Mora. Determining the Number of Factors to Retain in EFA: An easy-to-use computer program for carrying out Parallel Analysis Ledesma. — 2007.
16. Buja A., Eyuboglu, N. Remarks on parallel analysis // Mult. Beh. Res. 27. - 1992. — P. 509^540.
17. Seghouane A. and Cichocki A. Bavesian estimation of the number of principal components. Signal Processing, - 2007. - V. 87(3). - P. 562-568.
18. Ratsimalahelo Z. Strongly Consistent Determination of the Rank of Matrix // Econometrics. — 2003.
19. Dalalyan A. S., Juditsky A., Spokoiny V. A new algorithm for estimating the effective dimension reduction subspace // Journal of Machine Learning Research. — 2008. — V. 9. - P. 1647-1678.
20. Zhu Y. and Zeng P. Fourier methods for estimating the central subspace and the central mean subspace in regression // Journal of the American Statistical Association. — 2006. _ у. loi. _ R 1638-1651.
21. GoProblems
http://www-optima.amp i.kyoto-u.ac.jp/member/student/hedar/Hedar_files/TestGO .htm
Поступила в редакцию 11.12.2012.