首页 文章

用python拟合Tanh曲线

提问于
浏览 次
0

我需要适合像这样的tanh曲线：

import numpy as np
import matplotlib.pyplot as plt
from lmfit import Model


def f(x, a1=0.00010, a2=0.00013, a3=0.00013, teta1=1, teta2=0.00555, teta3=0.00555, phi1=-50, phi2=600, phi3=-900,
      a=0.000000019, b=0):

    formule = a1 * np.tanh(teta1 * (x + phi1)) + a2 * np.tanh(teta2 * (x + phi2)) + a3 * np.tanh(
        teta3 * (x + phi3)) + a * x + b

    return formule

# generate points used to plot
x_plot = np.linspace(-10000, 10000, 1000)

gmodel = Model(f)

result = gmodel.fit(f(x_plot), x=x_plot, a1=1,a2=1,a3=1,teta1=1,teta2=1,teta3=1,phi1=0,phi2=0,phi3=0)

plt.plot(x_plot, f(x_plot), 'bo')
plt.plot(x_plot, result.best_fit, 'r-')
plt.show()

我尝试做那样的事，但我得到了这个结果：

enter image description here

还有另一种方法可以拟合这条曲线吗？我不知道我做错了什么？

python curve-fitting data-fitting

4 回答

  • 1

    你的功能有点混乱，你没有真正的功能值 . 你基本上想要适应你自己的功能 . 理想情况下，您希望用真实的实验数据替换 curve_fit() 中的 f(x_plot) .

    适合函数的好方法是使用 scipy.optimize.curve_fit

    from scipy.optimize import curve_fit

popt, pcov = curve_fit(f, x_plot, f(x_plot), p0=[0.00010, 0.00013, 0.00013, 1, 0.00555, .00555, -50, 600, -900,
  0.000000019, 0])

plt.plot(f(x_plot, *popt))

    结果看起来像这样
    enter image description here

    回复于
  • 1

    与真实数据：

    test_X = np.array(
    [-9.77073e+03, -9.29706e+03, -8.82339e+03, -8.34979e+03, -7.87614e+03, -7.40242e+03, -6.92874e+03, -6.45506e+03,
     -5.98143e+03, -5.50771e+03, -5.03404e+03, -4.56012e+03, -4.08674e+03, -3.61304e+03, -3.13937e+03, -2.66578e+03,
     -2.19210e+03, -1.71845e+03, -1.24478e+03, -9.78925e+02, -9.29077e+02, -8.79059e+02, -8.29082e+02, -7.79092e+02,
     -7.29080e+02, -6.79084e+02, -6.29061e+02, -5.79078e+02, -5.29103e+02, -4.79089e+02, -4.29094e+02, -3.79071e+02,
     -3.29074e+02, -2.79062e+02, -2.29079e+02, -1.92907e+02, -1.72931e+02, -1.52930e+02, -1.32937e+02, -1.12946e+02,
     -9.29511e+01, -7.29438e+01, -5.29292e+01, -3.29304e+01, -1.29330e+01, 7.04455e+00, 2.70676e+01, 4.70634e+01,
     6.70526e+01, 8.70340e+01, 1.07056e+02, 1.27037e+02, 1.47045e+02, 1.67033e+02, 1.87039e+02, 2.20765e+02,
     2.70680e+02, 3.20699e+02, 3.70693e+02, 4.20692e+02, 4.70696e+02, 5.20704e+02, 5.70685e+02, 6.20710e+02,
     6.70682e+02, 7.20705e+02, 7.70707e+02, 8.20704e+02, 8.70713e+02, 9.20691e+02, 9.70700e+02, 1.23926e+03,
     1.73932e+03, 2.23932e+03, 2.73926e+03, 3.23924e+03, 3.73926e+03, 4.23952e+03, 4.73926e+03, 5.23930e+03,
     5.71508e+03, 6.21417e+03, 6.71413e+03, 7.21412e+03, 7.71410e+03, 8.21405e+03, 8.71402e+03, 9.21423e+03])

test_Y = np.array(
    [-3.17679e-04, -3.27541e-04, -3.51184e-04, -3.60672e-04, -3.75965e-04, -3.86888e-04, -4.03222e-04, -4.23262e-04,
     -4.38526e-04, -4.51187e-04, -4.61081e-04, -4.67121e-04, -4.96690e-04, -4.94811e-04, -5.10110e-04, -5.18985e-04,
     -5.11754e-04, -4.90964e-04, -4.36904e-04, -3.93638e-04, -3.83336e-04, -3.71110e-04, -3.57207e-04, -3.39643e-04,
     -3.24155e-04, -2.97296e-04, -2.74653e-04, -2.43700e-04, -1.95574e-04, -1.60716e-04, -1.43363e-04, -1.33610e-04,
     -1.30734e-04, -1.26332e-04, -1.26063e-04, -1.24228e-04, -1.23424e-04, -1.20276e-04, -1.16886e-04, -1.21865e-04,
     -1.16605e-04, -1.14148e-04, -1.14728e-04, -1.14660e-04, -1.16927e-04, -1.10380e-04, -1.09836e-04, 4.24232e-05,
     8.66095e-05, 8.43905e-05, 9.09867e-05, 8.95580e-05, 9.02585e-05, 8.87033e-05, 8.86536e-05, 8.92236e-05,
     9.24438e-05, 9.27929e-05, 9.24961e-05, 9.72166e-05, 1.00432e-04, 1.05457e-04, 1.11278e-04, 1.14716e-04,
     1.25818e-04, 1.40721e-04, 1.62968e-04, 1.91776e-04, 2.28125e-04, 2.57918e-04, 2.88941e-04, 3.85003e-04,
     4.91916e-04, 5.32483e-04, 5.50929e-04, 5.45350e-04, 5.38903e-04, 5.27765e-04, 5.15592e-04, 4.95717e-04,
     4.81722e-04, 4.69538e-04, 4.58643e-04, 4.41407e-04, 4.29820e-04, 4.07784e-04, 3.92236e-04, 3.81761e-04])

    我试试这个：

    import numpy,
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from scipy.optimize import differential_evolution
import warnings


def function(x, a1, a2, a3, teta1, teta2, teta3, phi1, phi2, phi3, a, b):
    import numpy as np
    formule = a1 * np.tanh(teta1 * (x + phi1)) + a2 * np.tanh(teta2 * (x + phi2)) + a3 * np.tanh(teta3 * (x + phi3)) + a * x  + b
    return formule

# function for genetic algorithm to minimize (sum of squared error)
def sumOfSquaredError(parameterTuple):
    warnings.filterwarnings("ignore")  # do not print warnings by genetic algorithm
    val = function(test_X, *parameterTuple)
    return numpy.sum((test_Y - val) ** 2.0)


def generate_Initial_Parameters():

    parameterBounds = []
    parameterBounds.append([1.4e-04, 1.4e-04])
    parameterBounds.append([2.00e-04,2.0e-04])
    parameterBounds.append([2.5e-04, 2.5e-04])
    parameterBounds.append([0, 2.0e+01])
    parameterBounds.append([0, 4.0e-03])
    parameterBounds.append([0, 4.0e-03])
    parameterBounds.append([-8.e+01, 0])
    parameterBounds.append([0, 9.0e+02])
    parameterBounds.append([-2.1e+03, 0])
    parameterBounds.append([-3.4e-08, -2.4e-08])
    parameterBounds.append([-2.2e-05*2, 4.2e-05])

    # "seed" the numpy random number generator for repeatable results
    result = differential_evolution(sumOfSquaredError, parameterBounds)
    return result.x


# generate initial parameter values
geneticParameters = generate_Initial_Parameters()

# curve fit the test data
fittedParameters, pcov = curve_fit(function, test_X, test_Y, geneticParameters)

print('Parameters', fittedParameters)

modelPredictions = function(test_X, *fittedParameters)

absError = modelPredictions - test_Y

SE = numpy.square(absError)  # squared errors
MSE = numpy.mean(SE)  # mean squared errors
RMSE = numpy.sqrt(MSE)  # Root Mean Squared Error, RMSE
Rsquared = 1.0 - (numpy.var(absError) / numpy.var(test_Y))
print('RMSE:', RMSE)
print('R-squared:', Rsquared)

ytry = ftry(test_X)

##########################################################
# graphics output section
def ModelAndScatterPlot(graphWidth, graphHeight):
    f = plt.figure(figsize=(graphWidth / 100.0, graphHeight / 100.0), dpi=100)
    axes = f.add_subplot(111)

    # first the raw data as a scatter plot
    axes.plot(test_X, test_Y, 'D')

    # create data for the fitted equation plot

    yModel = function(test_X, *fittedParameters)

    # now the model as a line plot
    axes.plot(test_X, yModel)

    axes.set_xlabel('X Data')  # X axis data label
    axes.set_ylabel('Y Data')  # Y axis data label
    axes.plot(test_X, ytry)

    plt.show()
    plt.close('all')  # clean up after using pyplot


graphWidth = 800
graphHeight = 600
ModelAndScatterPlot(graphWidth, graphHeight)

    R平方：0.9978，不完美，但不是那么糟糕

    enter image description here

    回复于
  • 2

    基本上你的拟合很好（虽然从编码的角度来看不是很好） . 像往常一样，非线性拟合强烈依赖于初始参数 . 你的选择很糟糕 . 您可以考虑如何手动确定它们，或者使用 differential_evolution 中的 differential_evolution 预制包 . 我没有使用这个包，但你可以在SE上找到一个例子here

    回复于
  • 1

    我同意mikuszefski和F. Win的答案，但我想补充一点 .

    您的模型包括3行tanh函数 . 数据支持许多不同的tanh功能并不完全清楚 . 如果是这样（并回应mikuszefki），你需要告诉它们这些不相同 . 你的例子开始它们是相同的，这将使适合找到一个好的解决方案非常困难 . 无论哪种方式，能够轻松测试是否真的有1,2,3或更多tanh函数可能会有所帮助 .

    您可能还希望不仅为参数提供初始值，还要为它们提供真实的边界，以便tanh函数清晰分离，并且不会偏离它们应该的位置太远 .

    要清理代码并更好地允许您更改使用的tanh函数的数量并放置边界约束，我建议制作单独的模型并将其添加为：

    from lmfit import Model

def f_tanh(x, eta=1, phi=0):
    "tanh function"
    return np.tanh(eta * (x + phi))

def f_line(x, slope=0, intercept=0):
    "line function"
    return slope*x + intercept

# create model as line + 2 tanh functions
gmodel = Model(f_line) + Model(f_tanh, prefix='t1_') + Model(f_tanh, prefix='t2_')

    现在您可以轻松创建参数了

    params = gmodel.make_params(slope=0.003, intercept=0.001,
                            t1_eta=0.021, t1_phi=-2000,
                            t2_eta=0.013, t2_phi=600)

    通过定义fit参数，您可以使用以下边界放置边界：

    params['t1_eta'].min = 0
params['t2_eta'].min = 0

params['t1_phi'].min = -3000
params['t1_phi'].max = -1000

params['t2_phi'].min = 0
params['t2_phi'].max = 1000

    我认为所有这些都将帮助您更好地探索数据及其适应性 . 把这一切放在一起，你可能会：

    import numpy as np
import matplotlib.pyplot as plt
from lmfit import Model

def f_tanh(x, eta=1, phi=0):
    "tanh function"
    return np.tanh(eta * (x + phi))

def f_line(x, slope=0, intercept=0):
    "line function"
    return slope*x + intercept

# line + 2 tanh functions
gmodel = Model(f_line) + Model(f_tanh, prefix='t1_') + Model(f_tanh, prefix='t2_')

# generate "data"
x = np.linspace(-10000, 10000, 1000)
y = gmodel.eval(x=x, slope=0.0001,
                t1_eta=0.010, t1_phi=-2100,
                t2_eta=0.004, t2_phi=740)
y = y + np.random.normal(size=len(x), scale=0.02)

# make parameters with initial values
params = gmodel.make_params(slope=0.003, intercept=0.001,
                            t1_eta=0.021, t1_phi=-2000,
                            t2_eta=0.013, t2_phi=600)

# place realistic but generous constraints to keep tanhs separate
params['t1_eta'].min = 0
params['t2_eta'].min = 0

params['t1_phi'].min = -3000
params['t1_phi'].max = -1000

params['t2_phi'].min = 0
params['t2_phi'].max = 1000

result = gmodel.fit(y, params, x=x)

print(result.fit_report())

plt.plot(x, y, 'bo')
plt.plot(x, result.best_fit, 'r-')
plt.show()

    这将给出良好的拟合和绘图，并在噪声水平内找到预期值 . 希望有助于您指出正确的方向 .

    回复于

相关问题