Подгонка параметра внутри интеграла с использованием python (или другого полезного языка) - PullRequest
2 голосов
/ 13 января 2020

У меня есть набор данных, в основном с информацией о f (x) как функции от x и самого x. Я знаю из теории проблемы, что я работаю над форматом f (x), который представлен в виде выражения ниже:

Eq. 1

По сути, я хочу использовать этот набор данных, чтобы найти параметры a и b. Моя проблема: как я могу это сделать? Какую библиотеку я должен использовать? Я хотел бы получить ответ, используя Python. Но Р или Джулия тоже были бы в порядке.

Из всего, что я до сих пор делал, я читал о функциональности, называемой fit fit , из библиотеки SciPy, но у меня возникли некоторые проблемы, в которых я мог бы выполнить код как long Моя переменная x находится в одном из пределов интеграции.

Для более эффективных способов работы с проблемой у меня также есть следующие ресурсы:

Пример набора, для которого я знаю параметры, которые я ищу. К этому набору я знаю, что a = 2 и b = 1 (и c = 3). И прежде чем возникают некоторые вопросы о том, как я знаю эти параметры: я знаю их, потому что я создал этот выборочный набор, используя эти параметры из интегрирования приведенного выше уравнения, просто чтобы использовать образец, чтобы исследовать, как я могу найти их и получить ссылку.

У меня также есть этот набор , для которого у меня есть единственная информация: c = 4 и я хочу найти a и b.

Я бы тоже хотел отметить, что:

i) сейчас у меня нет кода для публикации здесь, потому что я не имею понятия, как написать что-то, чтобы решить мою проблему. Но я был бы рад отредактировать и обновить вопрос после прочтения любого ответа или помощи, которую вы, ребята, могли бы предоставить мне.

ii) Сначала я ищу решение, в котором я не знаю a и b. Но в случае, если это слишком сложно, я был бы рад найти какое-то решение, в котором, как я полагаю, известен либо a, либо b.

РЕДАКТИРОВАТЬ 1: Я хотел бы сослаться на этот вопрос всем, кто интересуется этой проблемой, поскольку это параллельное, но и важное обсуждение проблемы, с которой здесь сталкиваются

Ответы [ 2 ]

1 голос
/ 14 января 2020

Я бы использовал чисто числовой подход c, который вы можете использовать, даже если вы не можете напрямую решить интеграл. Вот снимок для подгонки только параметра a:

import numpy as np
from scipy.optimize import curve_fit
import pandas as pd
import matplotlib.pyplot as plt

def integrand(x, a):
    b = 1
    c = 3
    return 1/(a*np.sqrt(b*(1+x)**3 + c*(1+x)**4))

def integral(x, a):
    dx = 0.001
    xx = np.arange(0, x, dx)
    arr = integrand(xx, a)
    return np.trapz(arr, dx=dx, axis=-1)

vec_integral = np.vectorize(integral)

df = pd.read_csv('data-with-known-coef-a2-b1-c3.csv')
x = df.domin.values
y = df.resultados2.values
out_mean, out_var = curve_fit(vec_integral, x, y, p0=[2])

plt.plot(x, y)
plt.plot(x, vec_integral(x, out_mean[0]))
plt.title(f'a = {out_mean[0]:.3f} +- {np.sqrt(out_var[0][0]):.3f}')
plt.show()

vec_integral = np.vectorize(integral)

enter image description here

Конечно, вы можете уменьшить значение dx, чтобы получить желаемая точность. В то время как для подгонки только a, когда вы пытаетесь также запустить b, подгонка не сходится должным образом (на мой взгляд, потому что a и b сильно коррелированы). Вот что вы получите:

def integrand(x, a, b):
    c = 3
    return 1/(a*np.sqrt(np.abs(b*(1+x)**3 + c*(1+x)**4)))

def integral(x, a, b):
    dx = 0.001
    xx = np.arange(0, x, dx)
    arr = integrand(xx, a, b)
    return np.trapz(arr, dx=dx, axis=-1)

vec_integral = np.vectorize(integral)

out_mean, out_var = sp.optimize.curve_fit(vec_integral, x, y, p0=[2,3])
plt.title(f'a = {out_mean[0]:.3f} +- {np.sqrt(out_var[0][0]):.3f}\nb = {out_mean[1]:.3f} +- {np.sqrt(out_var[1][1]):.3f}')

plt.plot(x, y, alpha=0.4)
plt.plot(x, vec_integral(x, out_mean[0], out_mean[1]), color='green', label='fitted solution')
plt.plot(x, vec_integral(x, 2, 1),'--', color='red', label='theoretical solution')
plt.legend()
plt.show()

enter image description here

Как видите, даже если результирующие параметры a и b образуют соответствие «не хорошо», сюжет очень похож.

1 голос
/ 14 января 2020

Это три переменные a, b, c, которые не являются независимыми. Один из них должен быть указан, если мы хотим вычислить два других благодаря регрессии. С учетом c решение для a, b простое:

enter image description here

Пример численного исчисления ниже сделан с небольшими данными (n = 10) чтобы облегчить проверку как для y (x), когда данные разбросаны (результат такой же, если нет разброса).

Если для y (x) абсолютно необходимо иметь регрессию, необходима нелинейная регрессия. Это включает в себя итеративный процесс, начинающийся с достаточно начального предположения для a, b. Вышеприведенное исчисление дает очень хорошие начальные значения.

ДОПОЛНИТЕЛЬНО:

Тем временем Андреа отправил соответствующий ответ. Конечно, согласование с его методом лучше, потому что это нелинейная регрессия, а не линейная, как уже указывалось в примечании выше.

Тем не менее, несмотря на различные значения (a = 1.881; b = 1.617) по сравнению с (a = 2.346, b = -0.361) соответствующие кривые, приведенные ниже, находятся недалеко друг от друга:

Синяя кривая: от линейной регрессии (метод выше)

Зеленая кривая: из нелинейной регрессии (Андреа)

enter image description here

Случай второго набора данных

https://mega.nz/#! echEjQyK ! tUEx0gpFND7gucvsTONiB_wn-ewBq-5k-pZlfLxmfvw

Регрессия не выполняется, поскольку предположение c = 3 неверно.

В случае c = 0 аналитическое c исчисление интеграла отличается от вышеуказанного:

enter image description here

...