Подгоните 1000 образцов к нормальному распределению, используя Максимальное правдоподобие в C++
У меня есть образец 5000 парных как
образец = {1,23, -4,67, 0,17, 1,25, 6,89, -2,03, ...}
и хотим подогнать данные к параметрическим распределениям, таким как N(mu, sigma) или обобщенный студент t(loc, scale, DoF)...
У меня уже есть PDF-файлы этих распределений PDF_normal(mu, sigma)(x) и PDF_t(loc, scale, DoF)(x), и я могу рассчитать сумму логарифмов PDF-файлов для 5000 образцов для фиксированных параметров распределения.
Теперь я хочу использовать некоторые алгоритмы C++ для решения задач нелинейной оптимизации, чтобы найти параметры (mu_max, sigma_max) или (loc_max, scale_max, DoF_max), которые дадут мне значения максимального логарифмического правдоподобия.
Проект R для статистических вычислений решает проблему в пакете MASS следующим образом: .. прямая оптимизация логарифмического правдоподобия выполняется с помощью optim. Расчетные стандартные ошибки взяты из наблюдаемой информационной матрицы, рассчитанной численным приближением. Для одномерных задач используется метод Нелдера-Мида, а для многомерных задач - метод BFGS...
К сожалению, я не могу использовать решение R, но мне нужно найти решение в Microsoft VS2010 C++, и я не хочу самостоятельно писать код оптимизации, и я не хочу смотреть на исходный код R и переписывать его для C++...
Какие-нибудь предложения, где я могу найти быструю и хорошо проверенную реализацию BFGS (или подобную) для C++?
Есть ли что-нибудь доступное в Boost, Intel MKL и т. Д.?
Спасибо за вашу помощь, Мэтт
2 ответа
Хорошо, мне не нужна оптимизация для MLE нормального распределения, потому что это может быть решено в закрытой форме, см. 1: http://de.wikipedia.org/wiki/Maximum-Likelihood-Methode
Но я хочу решить эту проблему для разных семейств дистрибутивов, где я знаю только PDF-файлы. Таким образом, мне все еще нужна хорошая реализация C++ нелинейного решателя...
Все упомянутые вами дистрибутивы имеют очень мало параметров. Возможно, было бы более полезно сделать прямой метод Ньютона. Вы можете определить функцию логарифмического правдоподобия и ее градиент и гессиан (в отношении параметров, а не данных).
Когда вы решаете систему в конце, вы должны использовать факторизацию Холецкого или ЛПНП-факторизацию вместо гауссовой элиминации. Вы будете точно знать все факторы с точностью до sqrt(машинный эпсилон), что в большинстве систем составляет 10-7.
Более надежный способ получить факторизацию Холески по Гессиану, если она является суммой набора матриц ранга один, состоит в том, чтобы включить каждую новую матрицу ранга один в факторизацию Холецкого, используя последовательность вращений Гивенса. Это дает вам точность, близкую к машинному эпсилону, но примерно в два раза медленнее, чем формирование гессиана и последующая обработка факторизации Холецкого обычным способом.
Технология здесь достаточно проста, так что вы, вероятно, можете написать ее самостоятельно.
Мне сложно предложить вам универсальный программный пакет. Это отчасти потому, что у меня не было положительного опыта с библиотеками оптимизации C++ других людей. В основном, тем не менее, написание метода Ньютона не так уж много работы, когда у вас есть код для оценки градиента и Гессиана, так что я действительно не чувствовал необходимости.
Тем не менее, вы можете посмотреть и посмотреть, что может предложить проект COIN-OR. Возможно, вы сможете сэкономить время, используя инструменты автоматического дифференцирования --- (у COIN-OR есть один). Полагаю, я должен упомянуть, что я никогда не использовал для нелинейной оптимизации что-то из COIN-OR или инструмент автоматического дифференцирования.