Псевдокод оценки максимального правдоподобия - PullRequest
24 голосов
/ 11 октября 2011

Мне нужно кодировать Оценщик максимального правдоподобия, чтобы оценить среднее значение и дисперсию некоторых игрушечных данных.У меня есть вектор из 100 сэмплов, созданный с помощью numpy.random.randn(100).Данные должны иметь нулевое среднее и единичную дисперсию гауссовского распределения.

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

ЕстьЕсть ли псевдокод для оценки максимального правдоподобия?Я понимаю интуицию MLE, но не могу понять, с чего начать кодирование.

Вики говорит, что принимает argmax логарифмического правдоподобия.Что я понимаю: мне нужно рассчитать логарифмическую вероятность, используя разные параметры, а затем я возьму параметры, которые дали максимальную вероятность.Чего я не понимаю: где я найду параметры в первую очередь?Если я случайно выберу другое среднее значение и дисперсию, чтобы получить высокую вероятность, когда я должен прекратить попытки?

Ответы [ 4 ]

36 голосов
/ 21 августа 2013

Я только что сталкивался с этим, и я знаю, что он старый, но я надеюсь, что кто-то еще выиграет от этого.Хотя предыдущие комментарии давали довольно хорошее описание того, что такое оптимизация ML, никто не давал псевдокод для ее реализации.Python имеет минимизатор в Scipy, который сделает это.Вот псевдокод для линейной регрессии.

# import the packages
import numpy as np
from scipy.optimize import minimize
import scipy.stats as stats
import time

# Set up your x values
x = np.linspace(0, 100, num=100)

# Set up your observed y values with a known slope (2.4), intercept (5), and sd (4)
yObs = 5 + 2.4*x + np.random.normal(0, 4, 100)

# Define the likelihood function where params is a list of initial parameter estimates
def regressLL(params):
    # Resave the initial parameter guesses
    b0 = params[0]
    b1 = params[1]
    sd = params[2]

    # Calculate the predicted values from the initial parameter guesses
    yPred = b0 + b1*x

    # Calculate the negative log-likelihood as the negative sum of the log of a normal
    # PDF where the observed values are normally distributed around the mean (yPred)
    # with a standard deviation of sd
    logLik = -np.sum( stats.norm.logpdf(yObs, loc=yPred, scale=sd) )

    # Tell the function to return the NLL (this is what will be minimized)
    return(logLik)

# Make a list of initial parameter guesses (b0, b1, sd)    
initParams = [1, 1, 1]

# Run the minimizer
results = minimize(regressLL, initParams, method='nelder-mead')

# Print the results. They should be really close to your actual values
print results.x

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

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

18 голосов
/ 11 октября 2011

Если вы выполняете вычисления максимального правдоподобия, первый шаг, который вам нужно сделать, заключается в следующем: Предположите распределение, которое зависит от некоторых параметров. Поскольку вы generate свои данные (вы даже знаете свои параметры), вы "говорите" своей программе принять гауссово распределение. Однако вы не указываете своей программе свои параметры (0 и 1), а оставляете их априори неизвестными и затем вычисляете их.

Теперь у вас есть примерный вектор (назовем его x, его элементы от x[0] до x[100]), и вы должны его обработать. Для этого необходимо вычислить следующее (f обозначает функцию плотности вероятности гауссовского распределения ):

f(x[0]) * ... * f(x[100])

Как вы можете видеть в приведенной мной ссылке, f использует два параметра (греческие буквы µ и σ). Вы теперь должны вычислить значения для µ и σ таким образом, чтобы f(x[0]) * ... * f(x[100]) принял максимально возможное значение.

Когда вы это сделаете, µ - это максимальное значение вероятности для среднего значения, а σ - максимальное значение вероятности для стандартного отклонения.

Обратите внимание, что я не говорю вам явно , как вычислять значения для µ и σ, поскольку это довольно математическая процедура, которой у меня нет под рукой (и, вероятно, я бы ее не понял ); Я просто расскажу вам технику получения значений, которая может быть применена и к любым другим дистрибутивам.

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

Если вы действительно хотите рассчитать его, вы можете сделать несколько упрощений, которые приведут к следующему термину (надеюсь, я ничего не испортил):

enter image description here

Теперь вы должны найти значения для µ и σ такие, чтобы вышеупомянутый зверь был максимальным. Это очень нетривиальная задача, называемая нелинейной оптимизацией.

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

4 голосов
/ 11 октября 2011

Вам нужна процедура численной оптимизации. Не уверен, что что-то реализовано в Python, но если это так, то это будет в numpy или scipy и друзья.

Ищите такие вещи, как «алгоритм Nelder-Mead» или «BFGS». Если ничего не помогает, используйте Rpy и вызовите функцию R 'optim ()'.

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

0 голосов
/ 17 июля 2013

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

В случае нормального распределения вы должны получить логарифмическую правдоподобие по отношению к среднему (mu) и затем вывести по дисперсии (сигма ^ 2), чтобы получить два уравнения, оба равные нулю.Решив уравнения для mu и sigma ^ 2, вы получите выборочное среднее значение и выборочную дисперсию в качестве ответов.

Подробнее см. На странице википедии .

...