Блог

Я определил следующее нормальное распределение N. Здесь r-случайная величина (вы можете думать о r как о «возрасте»), в то время как среднее значение N задается функцией P(d), которая (в качестве параметра) фиксирует N каждый раз (вы можете думать о d как о «росте»):

 def p(d, a, b):
    return a-b*d

def N(r, d, a, b, s):
    return (1/(s*sqrt(2*pi)))*exp(-(1/2)*((r-p(d, a, b))/s)**2)
 

Другими словами, для разных значений d (рост) N становится другим PDF-файлом (в форме a, b и s), который описывает случайную величину r (возраст).

У меня есть много (18 миллионов) пар d, r, и я хотел бы разместить PDF-файл на этих данных, найдя оптимальные a, b и s.

Как я могу это сделать?

Итак, вы хотите найти параметры a, b и s, которые максимизируют вероятность получения данных? Таким образом, я могу предположить, что ваша функция потерь будет произведением N(r, d, a, b, s) с учетом r и d в ваших данных. Существует множество методов оптимизации, учитывая, что эти функции дифференцируемы, вы даже можете использовать фреймворк autograd, такой как Tensorflow или PyTorch. Но для простоты я буду использовать scipy, как вы его отметили, это должно быть нормально, если ваши данные невелики (

 import numpy as np
import scipy.optimize
from numpy import pi, sqrt, exp, log

def p(d, a, b):
    return a-b*d

def N(r, d, a, b, s): # Writen as numpy-friendly (accepts numpy arrays as inputs)
    return (1/(s*sqrt(2*pi)))*exp(-(1/2)*((r-p(d, a, b))/s)**2)

def minus_log_likelihood(p): # params, a, b, s. Log sum is equivalent to product
    return -np.sum(log(N(dataset[:, 0], dataset[:, 1], p[0], p[1], p[2])))

dataset = np.random.uniform(size=(100, 2)) # 100 points with d and r values
res = scipy.optimize.minimize(minus_log_likelihood, [0, 0, 1])