Блог

У меня есть набор данных, и я хочу наложить на него параболическую подгонку. Это уже работает с функцией polyfit от numpy следующим образом:

 fit = np.polyfit(X, y, 2)
formula = np.poly1d(fit)
 

Теперь я хочу, чтобы парабула имела свое пиковое значение при фиксированном значении x и чтобы подгонка по-прежнему выполнялась как можно лучше с этим фиксированным пиком. Есть ли способ добиться этого?

Из моих данных я знаю, что парабола всегда будет открыта вниз.

Я думаю, что это довольно сложная проблема, поскольку координата x пика полинома второго порядка (a x^2 b x c) всегда лежит в x = -b/2a.

Что вы могли бы сделать, так это отбросить член b и компенсировать его на желаемое пиковое значение x при подгонке полинома, как показано в приведенном ниже коде. Обратите внимание, что раньше я подходил scipy.optimize.curve_fit для пользовательской функции func .

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

# generating a parabola with noise
np.random.seed(42)
x = np.linspace(-10, 10, 100)
y = 10 -(x-2)**2   np.random.normal(0, 5, x.shape)

# function to fit
def func(x, a, c):
    return a*x**2   c

# desired x peak value
x_peak = 2

popt, pcov = curve_fit(func, x - x_peak, y)

y_fit = func(x - x_peak, *popt)

# plotting
plt.plot(x, y, 'k.')
plt.plot(x, y_fit)
plt.axvline(x_peak)
plt.show()
 

Выводит изображение:

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

 y = A(x - B)**2   C
 

Учитывая коэффициенты a b , c в вашем первоначальном неограниченном соответствии у вас есть соотношения

 a = A
b = -2AB
c = AB**2   C
 

Единственная разница в том, что, поскольку B это константа, и у вас нет x - B члена в уравнении, вам нужно самостоятельно решить задачу наименьших квадратов. Учитывая массивы x y и константу B , проблема выглядит следующим образом:

 m = np.stack((x - B, np.ones_like(x)), axis=-1)
(A, C), *_ = np.linalg.lstsq(m, y, rcond=None)
 

Затем вы можете извлечь нормальный коэффициент из приведенных выше формул для a , b , c .

Вот полный пример, точно такой же, как в другом ответе:

 B = 2

np.random.seed(42)
x = np.linspace(-10, 10, 100)
y = 10 -(x - B)**2   np.random.normal(0, 5, x.shape)

m = np.stack(((x - B)**2, np.ones_like(x)), axis=-1)
(A, C), *_ = np.linalg.lstsq(m, y, rcond=None)

a = A
b = -2 * A * B
c = A * B**2   C

y_fit = a * x**2   b * x   c
 

Вы можете a b полностью бросить c и сделать

 y_fit = A * (x - B)**2   C
 

Результат будет одинаковым.

 plt.plot(x, y, 'k.')
plt.plot(x, y_fit)
 

Без условия расположения пика функция, которая должна быть установлена, будет :

y = a x^2 b x c

С условием расположения пика при x=p , учитывая p :

-b/(2a)=p

b=-2 a p

y = a x^2 -2 a p x c

y = a (x^2 — 2 p x) c

Зная p , одно изменение переменной :

X = x^2 -2 p x

Итак, из данных (x,y) сначала вычисляются новые данные (X,y).

Затем a и c вычисляются с помощью линейной регрессии

y = a X c