Блог

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

 def Maximum_Likelihood(param, pmf):
    
    i = symbols('i', positive=True)
    n = symbols('n', positive=True)
    
    Likelihood_function = Product(pmf, (i, 1, n)) 
    # calculate partial derivative for parameter (p for Bernoulli)
    deriv = diff(Likelihood_function, param) 
    equation_to_solve = Eq(deriv,0) # equate with 0
    
    # solve above equation and return parameter (p for Bernoulli)
    return solve(equation_to_solve, param) 
 

Param означает параметр, для которого я хочу знать оценку, а pmf — это функция массы вероятности.

И, например, я хочу получить оценку для параметра p в распределении Бернулли. Как должна выглядеть максимальная вероятность:

Мой код. Импорт:

 import numpy as np
import sympy as sym

from sympy.solvers import solve

from sympy import Product, Function, oo, IndexedBase, diff, Eq, symbols
 

Теперь, используя Sympy, я определил его:

 def Maximum_Likelihood(param, pmf):
    
    i = symbols('i', positive=True)
    n = symbols('n', positive=True)
    
    Likelihood_function = Product(pmf, (i, 1, n)) 
    deriv = diff(Likelihood_function, param) 
    equation_to_solve = Eq(deriv,0) 
    
    return solve(equation_to_solve, param) 
 

и пример Бернулли:

 x = IndexedBase('x')
i = symbols('i', positive=True)
n = symbols('n', positive=True)
formula = (p**x[i])*((1-p)**(1-x[i]))

Likelihood_function = Product(formula, (i, 1, n))
Likelihood_function
 

Когда я хочу получить результат maximum_likelibility(param, pmf):

 param = p 
pmf = formula
print(Maximum_Likelihood(param, pmf))
 

Я получаю «[]».
Я хочу получить оценку p, которая должна выглядеть так:

Не могли бы вы взглянуть на это и посоветовать, что я делаю не так. Спасибо!

По какой-то причине diff продукта фактически не оценивает производную, но вы можете использовать doit для принудительной оценки:

 In [14]: solve(Eq(diff(Product(p**x[i]*(1 - p)**(1 - x[i]), (i, 1, n)), p), 0).doit(), p)
Out[14]: 
⎡  n       ⎤
⎢ ___      ⎥
⎢ ╲        ⎥
⎢  ╲       ⎥
⎢  ╱   x[i]⎥
⎢ ╱        ⎥
⎢ ‾‾‾      ⎥
⎢i = 1     ⎥
⎢──────────⎥
⎣    n     ⎦
 

Это говорит о том, что оценка MLE для p — это всего лишь выборочное среднее значение данных, которые я предполагаю.