Материал из MachineLearning.
Выбор регрессионных моделей и гипотез порождения данных
Общий подход к сравнению нелинейных моделей заключается в следующем. Рассмотрим набор конкурирующих моделей 
. Априорная вероятность модели
определена как )
. При появлении данных 
апостериорная вероятность модели)
может быть найдена по теореме Байеса,
=\frac{P(f_i)p(D|f_i)}{\sum_{\iota=1}^Mp(D|f_\iota)P(f_\iota)},)
где )
- функция соответствия модели данным. Знаменатель дроби обеспечивает выполнение условия=1)
.
Вероятности моделей 
и 
, параметры которых идентифицированы по данным , сравнимы как
}{P(f_2|D)}=\frac{P(f_1)p(D|f_1)}{P(f_2)p(D|f_2)}.)
Отношение }{p(D|f_2)})
- есть отношение правдоподобия моделей. Отношение }{P(f_2)})
является априорной оценкой предпочтения одной модели другой. При моделировании отдается предпочтение наиболее простым и устойчивым моделям. Если априорные оценки моделей одинаковы, то есть, нет причины предпочитать одну модель другой, то их необходимо сравнивать по значениям :
=\int{}p(D|\w,f_i)p(\w|f_i)d\w.)
Апостериорная плотность распределения параметров 
функции при заданной выборке равна
=\frac{p(D|\w,f_i)p(\w|f_i)}{p(D|f_i)},)
где )
- априорно заданная плотность вероятности параметров начального приближения, )
- функция правдоподобия параметров модели, а знаменатель обеспечивает выполнение условия d\w=1)
. Он задан интегралом в пространстве параметров p(\w'|f_i)d\w'.)
Вышеприведенные формулы называются формулами Байесовского вывода первого и второго уровня.
Рассмотрим регрессию +\nu)
с аддитивным Гауссовским шумом с дисперсией 
и с нулевым матожиданием. Тогда плотность вероятности появления данных
\triangleq{}p(D|\w,\beta,f)=\frac{\exp(-\beta{}E_D(D|\w,f_i))}{Z_D(\beta)},)
где 
. Нормирующий множитель )
задан выражением
=\left(\frac{2\pi}{\beta}\right)^{\frac{N}{2}})
и взвешенный функционал ошибки в пространстве данных
-y_n)^2.)
Введем регуляризующий параметр 
, который отвечает за то, насколько хорошо модель должна соответствовать зашумленным данным. Функция плотности вероятности параметров с заданным гиперпараметром имеет вид
=\frac{\exp(-\alpha{}E_W(\w|f_i))}{Z_W(\alpha)},)
где --- обратная дисперсии распределения параметров, 
, а нормирующая константа 
определена дисперсией распределения параметров как
=\left(\frac{2\pi}{\alpha}\right)^\frac{W}{2}.)
Требование к малым значениям параметров предполагает Гауссовское априорное распределение с нулевым средним:
=\frac{1}{Z_W}\exp{(-\frac{\alpha}{2}||\w||^2)}.)
Так как переменные и 
являются параметрами распределения параметров модели, в дальнейшем будем называть их гиперпараметрами. Исключая нормирующую константу , которая не зависит от параметров и логарифмируя, получаем:
Эта ошибка регуляризирует параметры, начисляя штраф за их чрезмерно большие значения.
При заданных значениях гиперпараметров и для фиксированной функции получаем:
=\frac{p(D|\w,\beta)p(\w|\alpha)}{p(D|\alpha,\beta)}.)
Записывая функцию ошибки в виде =\alpha{}E_W+\beta{}E_D)
, получаем
=\frac{\exp(-S(\w|f_i))}{Z_S(\alpha,\beta)},)
где 
- нормирующий множитель.
Нахождение параметров модели
Рассмотрим итеративный алгоритм для определения оптимальных параметров и гиперпараметров 
при заданной модели . Корректный подход заключается в интегрировании всех неизвестных параметров и гиперпараметров. Апостериорное распределение параметров определяется как
=\iint{}p(\w,\alpha,\beta|D)d\alpha{}d\beta=\iint{}p(\w|\alpha,\beta,D)p(\alpha,\beta|D)d\alpha{}d\beta,)
что требует выполнить интегрирование апостериорного распределения параметров )
по пространству, размерность которого равна количеству параметров. Вычислительная сложность этого интегрирования весьма велика. Интеграл может быть упрощен при подходящем выборе начальных значений гиперпараметров.
Приближение интеграла заключается в том, что апостериорная плотность распределения гиперпараметров )
имеет выраженный пик в окрестности наиболее правдоподобных значений гиперпараметров 
. Это приближение известно как аппроксимация Лапласа. При таком допущении интеграл упрощается до
\approx{}p(\w|\alpha^\m,\beta^\m,D)\iint{p(\alpha,\beta|D)}d\alpha{}d\beta\approx{}p(\w|\alpha^\m,\beta^\m,D). )
Необходимо найти значения гиперпараметров, которые оптимизируют апостериорную плотность вероятности параметров, а затем выполнить все остальные расчеты, включающие )
при фиксированных значениях гиперпараметров.
Для нахождения функционала , который использует апостериорное распределение параметров, рассмотрим аппроксимацию ошибки )
на основе рядов Тейлора второго порядка:
\approx{}S(\w^\m)+\frac{1}{2}(\w-\w^\m)^TA(\w-\w^\m).)
В выражении нет слагаемого первого порядка, так как предполагается, что 
определяет локальный минимум функции ошибки, то есть }{\partial{}w_\xi}=0)
для всех значений 
. Матрица 
- это матрица Гессе функции ошибок:
=\beta\nabla^2{}E_D(\w^\m)+\alpha{}I. )
Обозначим первое слагаемое правой части через 
, тогда 
. Подставив полученное приближенное значение в уравнение и обозначив 
, получим
=\frac{1}{\hat{Z}_S}\exp\left(-S(\w^\m)-\frac{1}{2}\Delta\w^TA\Delta\w\right), )
Оценим нормирующую константу 
, необходимую для аппроксимации кривой Гаусса, как
)(2\pi)^\frac{W}{2}(\det{}A)^{-\frac{1}{2}}.)
Максимизируем функцию )
, изменяя значения гиперпараметров и . Это можно выполнить, интегрируя функцию плотности вероятности данных по пространству параметров :
=\int{}p(D|\w,\alpha,\beta)p(\w|\alpha,\beta)d\w=\int{}p(D|\w,\alpha,\beta)p(\w|\alpha)d\w, )
где второй интеграл справедлив по причине того, что распределение параметров не зависит от дисперсии шума в силу гипотезы о Гауссовском распределении шума. Для упрощения вычислений мы допускаем, что распределение )
является равномерным.
Запишем вероятность данных через функцию ошибки:
=\frac{1}{Z_D(\beta)}\frac{1}{Z_D(\alpha)}\int\exp(-S(\w)))d\w. )
Из значений )
, и предыдущего выражения получим:
=-\alpha{}E_W^\m-\beta{}E_D^\m-\frac{1}{2}\ln|A|+\frac{W}{2}\ln{\alpha}+\frac{N}{2}\ln\beta-\frac{N}{2}\ln{}(2\pi).)
Для того, чтобы оптимизировать это выражение относительно , найдем производную:
)
В этом выражении 
- собственные значения матрицы . Так как функция ошибки на данных не является квадратичной функцией параметров, как при линейной или RBF регрессии, то непосредственно оптимизировать величину невозможно, гессиан 
не является константой, а зависит от параметров . Так как мы принимаем 
для вектора , который зависит от выбора , то собственные значения косвенным образом зависят от . Таким образом, формула (\ref{dda}) игнорирует параметры модели.
С использованием этого приближения, производная с учетом равна
=-E_W^\m-\frac{1}{2}\sum_{j=1}^W\frac{1}{\lambda_j+\alpha}+\frac{W}{2\alpha}. )
Приравнивая последнее выражение к нулю и преобразовывая, получаем выражение для
Обозначим вычитаемое правой части через 
Те компоненты суммы, в которых 
привносят вклад, близкий к единице, а те компоненты суммы, в которых 
, привносят вклад, близкий к нулю. Таким образом, может быть интерпретирована как мера числа хорошо обусловленных параметров модели.
Для нахождения гиперпараметра рассмотрим задачу оптимизации. Обозначим через 
собственные значение матрицы 
. Так как 
, то 
и следовательно, 
Отсюда,
=\frac{1}{\beta}\sum_{j=1}^W\frac{\lambda_j}{\lambda_j+\alpha}.)
Дифференцируя, как и в случае нахождения , мы находим, что оптимальное значение определено как
Таким образом, на каждом шаге для модели определены гиперпараметры 
. При этом определена для всей модели, а и для каждой функции из суперпозиции. Так как оптимизация параметров 
дает нам положительно определенную форму гессиана, его собственные значения 
больше нуля, и, таким образом, меньше 
.
Мы имеем следующий итерационный процесс пересчета и :
})
Процесс сходится, так как увеличение ведет к увеличению , что на следующем шаге ведет к уменьшению .
Код, считающий гиперпараметры:
for m=1:limit
gamma(m)=0;
for i=1:size(Model.wFound,2)
gamma(m)=gamma(m)+max(alpha)/(lambda(i)+max(alpha));
end
for i=1:size(Model.wFound,2)
alpha(i)=(size(Model.wFound,2)-gamma(m))/Model.wFound(i)^2
end
beta(m)=(size(y,1)-gamma(m))/Model.errTest;
end
