У меня есть набор наблюдаемых данных, и я пытаюсь подогнать к нему кривую. Однако эта кривая генерируется системой дифференциальных уравнений (она не имеет аналитической формы). Установив начальные условия, я должен изменить параметры, чтобы соответствовать ему.
Моя попытка вычислить среднеквадратичную ошибку и минимизировать. Я знаю, что мои данные поступают очень точно c раз, поэтому я ввожу параметры и численно решаю дифференциальное уравнение. Я выбираю свои решения в определенные c раз из данных и вычисляю MSE между обоими наборами данных. Я использовал numpy .minimize, чтобы минимизировать эту функцию, но я думаю, что она работает неправильно. Даже когда я генерирую данные из указанных c параметров, он не может найти параметры (или делает это с слишком большой ошибкой). Мне также нужен гессиан, но это не имеет смысла с моими ожидаемыми результатами. Я публикую свой код:
def loss_function(x):
trial_sol = odeint(f, U0, t, args = (x[0], x[1]))
I = trial_sol[:,1]
t_sample = np.linspace(0,999,20)
std_des = np.std(i_noise)
i_trial_sample = [I[int(t)] for t in t_sample]
mse = mean_squared_error(i_trial_sample,i_noise)
return mse
res = minimize(loss_function,[0.8,0.07], method = 'TNC', bounds = ((0,None),(0,None)))
minima = res.x
res
## calculation of the Hessian
import numdifftools as nd
hess = nd.Hessian(loss_function, step = 1e-8,method = 'central')
h =hess(res.x)
cov_matrix = np.linalg.inv(h)
cov_matrix
Что-то важное, это то, что гессиан отправляет предупреждения о делении на ноль, когда он работает по умолчанию. Кроме того, метод минимизации имеет тенденцию останавливаться перед нахождением минимума. Буду признателен, если кто-нибудь поведет меня в правильном направлении!