python手写多项式拟合、曲线拟合

python 同时被 2 个专栏收录
41 篇文章 3 订阅
10 篇文章 0 订阅

上篇博客写完之后,终于发现自己线性回归入门!

然后洗澡的时候就在想一个问题,线性回归会了,写线性拟合是完全没问题的,但是np库的多项式拟合到底是怎么做出来的呢?突然灵光一闪多项式拟合?多变量的线性回归?好像发现了什么?重新理清一下思路。什么是多项式拟合?对的,这个问题以前没有好好思考过,以前的直观感觉是,不就是给一些x,y的点,让一个多项式去拟合,使得这个多项式的曲线看起来大致符合那些点。等等,好像忽略了什么,什么叫大致符合?对的,什么叫大致符合?答案已在心中了!就是让这些点距离这个拟合曲线的和最小!慢着,这其实就是线性回归的损失函数的定义嘛!最小化,损失函数!多项式多个参数怎么确定?多变量线性回归!每个一个次方项当成独立的一个变量!比如a*x^2+b*x+c这个标准的二次项,x^2当成x1,x当成x2,1当成x3!好了,整个思路有了,尝试一下实现吧!

import numpy as np

def my_fit(x, y, power):
    size = len(x)
    c = np.ones((1, size))
    x_arr = np.array(x)
    x_a=x_arr.copy()
    x_a.resize((1,size))
    x_mat = np.append(x_a, c, axis=0)
    y_arr = np.array(y)
    y_arr.resize((1, size))
    y_mat = np.mat(y_arr)
    for i in range(2, power+1):
        temp_x = x_arr**i
        temp_x.resize((1, size))
        x_mat = np.append(temp_x, x_mat, axis=0)
    x_mat=np.mat(x_mat)
    w = (x_mat * x_mat.T).I * x_mat * y_mat.T
    w0=w.T
    w0.resize(w0.size)
    return w0

def f(x):
    return 2*x ** 2 +x+3

x = np.linspace(-5, 5)
y = f(x) + 0.5*np.random.randn(len(x))  # 加入噪音
y_fit = np.polyfit(x, y, 2)  # 二次多项式拟合
print(y_fit)
w=my_fit(x,y,2)
print(w)
激动人心的时候到了,看下输出:
[1.99205043 1.00302882 3.08335708]
[1.99205043 1.00302882 3.08335708]

竟然和np库的结果完全一样!惊了!虽说实现代码不算优美,但是基本思路没啥区别!懒得调参用梯度下降了,直接解出来。

很好,还有一个疑问,scipy的曲线拟合怎么做的呢?一样嘛!当成多变量线性回归就好!

code:

import numpy as np
from scipy.optimize import curve_fit


def f_fit(x, a, b, c):
    return a*np.sin(x)+b*x+c


def f_test(x):
    return 2*np.sin(x)+3*x+1


def my_curve_fit1(fit_fun, x, y, p_num):  #p_num参数的个数
    size = len(x)
    if p_num <= 0 or p_num > size:
        print('no parameter to fit')
        return
    test_list = [0]*p_num
    test_list[0] = 1
    x_arr = np.array(fit_fun(x, *test_list))
    x_arr.resize((1, size))
    for i in range(1, p_num):
        test_list[i-1] = 0
        test_list[i] = 1
        temp_x = np.array(fit_fun(x, *test_list))
        temp_x.resize((1, size))
        x_arr = np.append(x_arr, temp_x, axis=0)
    x_mat = np.mat(x_arr)
    y_arr = np.array(y)
    y_arr.resize((1, size))
    y_mat = np.mat(y_arr)
    w = (x_mat * x_mat.T).I * x_mat * y_mat.T
    w0 = w.T
    w0.resize(w0.size)
    return w0

#为了强行接近库函数的接口,无需参数的个数,但实现有点捞
#很容易知道参数个数多于样本x的个数无法拟合
def my_curve_fit2(fit_fun, x, y):
    test_list = [0]
    size = len(x)
    p_num = 1  #参数个数默认为1
    for i in range(0, size):
        try:
            fit_fun(x, *test_list)
            break
        except:
            p_num += 1
            test_list.append(0)
    if p_num > size:
        print('can not fit')
        return
    test_list[0] = 1
    x_arr = np.array(fit_fun(x, *test_list))
    x_arr.resize((1, size))
    for i in range(1, p_num):
        test_list[i - 1] = 0
        test_list[i] = 1
        temp_x = np.array(fit_fun(x, *test_list))
        temp_x.resize((1, size))
        x_arr = np.append(x_arr, temp_x, axis=0)
    x_mat = np.mat(x_arr)
    y_arr = np.array(y)
    y_arr.resize((1, size))
    y_mat = np.mat(y_arr)
    w = (x_mat * x_mat.T).I * x_mat * y_mat.T
    w0 = w.T
    w0.resize(w0.size)
    return w0





x = np.linspace(-2*np.pi, 2*np.pi)
y = f_test(x)+0.3*np.random.randn(len(x))   #加入噪音
p_fit, prov = curve_fit(f_fit, x, y)  #曲线拟合
my_fit1 = my_curve_fit1(f_fit, x, y, 3)
my_fit2 = my_curve_fit2(f_fit, x, y)
print('sicpy库的曲线拟合')
print('a,b,c', p_fit)
print('手写曲线拟合1')
print('a,b,c', my_fit1)
print('手写曲线拟合2')
print('a,b,c', my_fit2)
输出:
sicpy库的曲线拟合
a,b,c [1.99131057 3.00473361 1.05173153]
手写曲线拟合1
a,b,c [1.99131057 3.00473361 1.05173153]
手写曲线拟合2
a,b,c [1.99131057 3.00473361 1.05173153]

结果竟然也就是一样!厉害了,竟然一样。好了,实现的小窍门用到了列表参数传递,想单独提出来的那个项,就令该项的参数为1,其他参数为0,然后就把每个项都单独提出来。一用多变量线性回归,搞掂!然后为了强行接近scipy库的曲线拟合的接口,不知道有几个参数拟合,怎么办呢?好像没有什么好的方法,只好用try来试了,反正传的参数个数不对就会报错的,就一个一个参数往上加咯。

以上。

  • 1
    点赞
  • 0
    评论
  • 14
    收藏
  • 一键三连
    一键三连
  • 扫一扫,分享海报

相关推荐
©️2020 CSDN 皮肤主题: 大白 设计师:CSDN官方博客 返回首页
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、C币套餐、付费专栏及课程。

余额充值