利用scipy求解方程组、拟合直线、圆、椭圆、抛物线

scipy.optimize作为优化模块可以实现任意曲线拟合，方程求根、非线性方程组求解、自定义代价函数求解等功能，下面给出了optimize中常用的几个子模块：
minimize：需要自己构建代价函数（有时也称损失函数，目标函数等），理论上可以求解任意最优化问题
curve_fit:可以拟合任意的显式函数曲线，对于隐式函数曲线不能拟合
fsolve：方程根、求解适定方程组，需要满足未知数数量等于方程数量
least_squares、leastsq：这两个可以用于求解最小二乘问题
root：求解方程根
minimize模块灵活度最多，功能最为强大，理论上上述子模块都可以看作minimize的特例：下面分别给出使用案例：
使用minimize拟合直线、圆、椭圆，使用fsolve、least_squares求解方程组，使用curve_fit拟合抛物线

from scipy.optimize import minimize
from scipy.optimize import curve_fit
from scipy.optimize import fsolve
from scipy.optimize import least_squares
from scipy.optimize import leastsq
from matplotlib import pyplot as plt
import numpy as np
import math
def line_fun(params,x,y):
    a,b=params
    return np.sum(np.square(x*a+b-y))

def circle_fun(params,x,y):
    #圆的一般方程：x^2+y^2+ax+by+c=0
    a,b,c=params
    return np.sum(np.square(x*x+y*y+a*x+b*y+c))

def ellipse_model2d(params, x, y):
    #椭圆一般方程：A * x ** 2 + B * y ** 2 + C * x * y + D * x + E * y + F=0
    A, B, C, D, E, F= params
    error=A * x ** 2 + B * y ** 2 + C * x * y + D * x + E * y + F
    squares = np.square(error)  # 将数组a每个元素平方
    sum_of_squares = np.sum(squares)  # 对平方后的数组求和
    return sum_of_squares


#超定方程组
def non_linear_equations(var):
    x = var[0]
    y = var[1]
    Func= np.empty((3))
    Func[0] = x**2+y-5
    Func[1] = x+y-3
    Func[2] = 4*x+y**2-9
    return Func

#适定方程组
def non_linear_equations2(var):
    x = var[0]
    y = var[1]
    Func= np.empty((2))
    Func[0] = x**2+y-5
    Func[1] = x+y-3
    return Func


def curve_line(x, a, b):
    return a*x**2+b*x+10


# 按间距中的绿色按钮以运行脚本。
if __name__ == '__main__':

    #生成直线数据
    x = np.linspace(-20, 20, 100)

    y = x * 2 + 1
    # plt.scatter(x,y)
    # plt.show()
    x0 = np.array([1.0, 0.8])
    #求解直线方程，当然这里也可以用curve_fit模块
    result = minimize(line_fun, x0, args=(x, y))
    print(r"line ab:", result.x)

    # 生成圆数据
    m = np.linspace(0, 2 * math.pi, 100)
    a = 3
    b = -6
    r = 5
    x = a + r * np.cos(m)
    y = b + r * np.sin(m)
    figure, axes = plt.subplots(1)
    axes.plot(x, y)
    axes.set_aspect(1)
    plt.show()
    #拟合圆
    result = minimize(circle_fun, np.ones(3), args=(x, y))
    print(r"circle abc:", result.x)

    #ellipse
    #data
    theta_samples = np.linspace(0, 20, 100)
    # 椭圆方位角
    alpha_samples = -45.0 / 180.0 * np.pi
    # 长轴长度
    a_samples = 1.0
    # 短轴长度
    b_samples = 2.0
    # 样本x 序列，并叠加正态分布的随机值
    x_samples = a_samples * np.cos(theta_samples) * np.cos(alpha_samples) \
                - b_samples * np.sin(theta_samples) * np.sin(alpha_samples) +1\
               # + np.random.randn(100) * 0.05 * a_samples
    # 样本y 序列 ，并叠加正态分布的随机值
    y_samples = b_samples * np.sin(theta_samples) * np.cos(alpha_samples) \
                + a_samples * np.cos(theta_samples) * np.sin(alpha_samples) +2\
               # + np.random.randn(100) * 0.05 * b_samples
    z_samples=np.zeros(100)
    plt.axes([0.16, 0.15, 0.75, 0.75])
    plt.scatter(x_samples, y_samples, color="magenta", marker="+",
                zorder=1, s=80, label="samples")
    plt.show()
    #fit ellipse
    result=minimize(ellipse_model2d,np.ones(6),args=(x_samples,y_samples))
    print("ellipse abcdef:",result.x)

    #solve non_linear_equations
    a = np.array([0.0,0])
    b = least_squares(non_linear_equations, a)
    print("non_linear_equations x0,x1:", b.x)
    a = np.array([0.0,0])
    b = fsolve(non_linear_equations2, a)#方程数量和未知数数量要保持一致
    print("non_linear_equations x0,x1:",b)

    #fit curve_line
    # 这部分生成样本点，对函数值加上高斯噪声作为样本点
    xdata = np.linspace(-5, 5, 50)
    #
    # a=2.5, b=1.3
    y = curve_line(xdata, 2.5, 1.3)

    np.random.seed(1)
    err_stdev = 0.2
    # 生成均值为0，标准差为err_stdev为0.2的高斯噪声
    y_noise = err_stdev * np.random.normal(size=len(xdata))
    ydata = y + y_noise
    plt.scatter(xdata, ydata, label='data')

    # 利用curve_fit作简单的拟合，popt为拟合得到的参数,pcov是参数的协方差矩阵
    popt_1, pcov = curve_fit(curve_line, xdata, ydata)
    print("curve_line_result:",popt_1)
    plt.plot(xdata, curve_line(xdata, *popt_1), 'r-', label='fit_1')
    plt.show()

结果：

line ab: [1.99999999 0.99999999]
circle abc: [-6.00000001 11.99999997 19.99999983]
ellipse abcdef: [-0.04833563 -0.04833556  0.05800266 -0.01933401  0.13533947 -0.04833538]
non_linear_equations x0,x1: [2. 1.]
non_linear_equations x0,x1: [-1.  4.]
curve_line_result: [2.50137612 1.30813455]

椭圆拟合结果：

抛物线拟合结果：

参考：
1
2
3
4
5
6