2022年 11月 16日

使用Python进行数据拟合

文章目录

使用Python进行数据拟合
多项式拟合
非多项式拟合

多项式拟合

任何一个函数都可以拆分成近似于这个函数的多项式表达。

多项式拟合需要用到的函数是numpy库当中的np.polyfit，它的使用方法为：

np.polyfit(x, y, deg, rcond=None, full=False, w=None, cov=False)

使用最小二乘法原理，根据已知的x与y对应值，拟合一个下面形式的多项式。

(

)

−

⋅

P(x)=P{_0}x^{deg}+P{_1}x^{deg-1}···+P_{deg}

$P (x) = P_{0} x^{d e g} + P_{1} x^{d e g - 1} \cdot \cdot \cdot + P_{d e g}$
返回一系列的系数

$P$ .

参数说明：
一般情况下，我们只需要用到前三个参数。

x	array类型，形状（M，），M 个样本点的 x 坐标`(x[i], y[i])`
y	array类型，形状 (M,) 或 (M, K)，样本点的 y 坐标。
deg	int型常量，拟合多项式的最高次项。

返回值：

	各个系数，

我们还可以使用polyval()来计算我们需要预测多项式的值.

import matplotlib.pyplot as plt
import numpy as np

if __name__ == "__main__":
  x = np.arange(1, 31, 1)
  y = np.array([20, 23, 26, 29, 32, 35, 38, 45, 53, 62, 73, 86, 101, 118, 138, 161, 188, 220, 257, 300, 350, 409, 478, 558, 651, 760, 887, 1035, 1208, 1410])

  z1 = np.polyfit(x, y, 3)              # 曲线拟合，返回值为多项式的各项系数
  p1 = np.poly1d(z1)                    # 返回值为多项式的表达式，也就是函数式子
  print(p1)
  y_pred = p1(x)                        # 根据函数的多项式表达式，求解 y
  #print(np.polyval(p1, 29))             #根据多项式求解特定 x 对应的 y 值
  #print(np.polyval(z1, 29))             #根据多项式求解特定 x 对应的 y 值
  plot1 = plt.plot(x, y, '*', label='original values')
  plot2 = plt.plot(x, y_pred, 'r', label='fit values')
  plt.title('')
  plt.xlabel('')
  plt.ylabel('')
  plt.legend(loc=3, borderaxespad=0., bbox_to_anchor=(0, 0))
  plt.show()

输出结果:

        3         2
0.1215 x - 3.045 x + 28.62 x - 34.47

则其拟合的函数为

0.1215

−

3.045

28.62

−

34.47

y=0.1215 x{^3} – 3.045 x {^2}+ 28.62 x – 34.47

$y = 0.1215 x^{3} - 3.045 x^{2} + 28.62 x - 34.47$
将拟合的曲线与原数据各点进行对比发现拟合效果良好。
在这里插入图片描述

非多项式拟合

如果需要进行多项式拟合，前提是必须大体上知道散点的大致曲线形式，也就是大致的函数的形式。
比如，例子中的散点看起来像是指数的函数分布，因此可以给出假设的函数

y=ba{^x}+c

$y = b a^{x} + c$
所以，只要给出具体的函数形式(可以是任意的，只要能写的出来皆可)，用最小二乘的方式去逼近和拟合，即求出函数的各项系数。

此时用到的是scipy.optimize包下的curve_fit函数了：

代码示例：

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

def func(x, a, b, c):
  return b * np.power(a, x) + c

if __name__ == "__main__":
  x = np.arange(1, 31, 1)
  y = np.array([20, 23, 26, 29, 32, 35, 38, 45, 53, 62, 73, 86, 101, 118, 138, 161, 188, 220, 257, 300, 350, 409, 478, 558, 651, 760, 887, 1035, 1208, 1410])

  popt, pcov = curve_fit(func, x, y)                # 曲线拟合，popt为函数的参数list
  y_pred = [func(i, popt[0], popt[1], popt[2]) for i in x]    # 直接用函数和函数参数list来进行y值的计算
  print(popt)

  plot1 = plt.plot(x, y, '*', label='original values')
  plot2 = plt.plot(x, y_pred, 'r', label='fit values')
  plt.title('')
  plt.xlabel('')
  plt.ylabel('')
  plt.legend(loc=3, borderaxespad=0., bbox_to_anchor=(0, 0))
  plt.show()

输出结果为：

[ 1.16791847 13.39168878  1.24633841]

那么就有

{

1.16791847

13.39168878

1.24633841

\left\{

\begin{aligned} a & = & 1.16791847 \\ b & = & 13.39168878 \\ c & = & 1.24633841 \end{aligned}

$\begin{aligned} a & = & \ 1.16791847 \\ b & =&\ 13.39168878 \\ c & = & 1.24633841 \end{aligned}$ \right.

$⎩ ⎪ ⎨ ⎪ ⎧ a b c = = = 1.16791847 13.39168878 1.24633841$
于是我们得到了拟合的函数为

13.39168878

∗

1.16791847

1.24633841

y=13.39168878*1.16791847{^x}+1.24633841

$y = 13.39168878 * 1.16791847^{x} + 1.24633841$
将拟合的曲线与原数据各点进行对比发现拟合效果良好。
在这里插入图片描述