假设原函数由一个三角函数和一个线性项组成

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline def f(x):
return np.sin(x) + 0.5 * x x = np.linspace(-2 * np.pi, 2 * np.pi, 50) plt.plot(x, f(x), 'b')
plt.grid(True)
plt.xlabel('x')
plt.ylabel('f(x)')

一、用回归方式逼近

1. 作为基函数的单项式

最简单的情况是以单项式为基函数——也就是说,b1=1,b2=x,b3=x2,b4=x3,... 在这种情况下,Numpy有确定最优参数(polyfit)和以一组输入值求取近似值(ployval)的内建函数。ployfit函数参数如下:

参数 描述

x

x坐标(自变量值)
y y坐标(因变量值)
deg 多项式拟合度
full 如果有真,返回额外的诊断信息
w 应用到y坐标的权重
cov 如果为真,返回协方差矩阵

典型向量化风格的polyfit和polyval线性回归(deg=7)应用方式如下:

reg = np.polyfit(x, f(x), deg=7)
ry = np.polyval(reg, x)
plt.plot(x,f(x),'b',label='f(x)')
plt.plot(x,ry,'r.',label='regression')
plt.grid(True)
plt.xlabel('x')
plt.ylabel('f(x)')

2. 单独的基函数

当选择更好的基函数组时,可以得到更好的回归结果。单独的基函数必须能通过一个矩阵方法定义(使用Numpy ndarray对象)。

numpy.linalg子库提供lstsq函数,以解决最小二乘化做问题
matrix = np.zeros((3+1,len(x)))
matrix[3,:] = np.sin(x)
matrix[2,:] = x **2
matrix[1,:] = x
matrix[0,:] = 1 reg = np.linalg.lstsq(matrix.T,f(x))[0]
ry = np.dot(reg,matrix) plt.plot(x,f(x),'b',label='f(x)')
plt.plot(x,ry,'r.',label='regression')
plt.legend(loc = 0)
plt.grid(True)
plt.xlabel('x')
plt.ylabel('f(x)')

对于有噪声的数据,未排序的数据,回归法都可处理。

3. 多维

以fm函数为例

def fm(x,y):
return np.sin(x) + 0.25 *x +np.sqrt(y) + 0.05 * y **2 x = np.linspace(0,10,20)
y = np.linspace(0,10,20)
X, Y = np.meshgrid(x,y)
Z = fm(X,Y)
x = X.flatten()
y = Y.flatten() from mpl_toolkits.mplot3d import Axes3D
import matplotlib as mpl fig = plt.figure(figsize=(9, 6))
ax = fig.gca(projection='3d')
surf = ax.plot_surface(X,Y,Z,rstride=2,cstride=2,cmap=mpl.cm.coolwarm,linewidth=0.5,antialiased=True)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('f(x,y)')
fig.colorbar(surf,shrink=0.5,aspect=5)

为了获得好的回归结果,编译一组基函数,包括一个sin函数和sqrt函数。statsmodels库提供相当通过和有益的函数OLS,可以用于一维和多维最小二乘回归。

matrix = np.zeros((len(x),6+1))
matrix[:,6] = np.sqrt(y)
matrix[:,5] = np.sin(x)
matrix[:,4] = y ** 2
matrix[:,3] = x **2
matrix[:,2] = y
matrix[:,1] = x
matrix[:,0] = 1 import statsmodels.api as sm
model = sm.OLS(fm(x,y),matrix).fit()
OSL函数的好处之一是提供关于回归及其而非结果的大量几百万来看信息。调用model.summary可以访问结果的一个摘要。单独统计数字(如确定系数)通常可以直接访问model.rsquared。最优回归系数,保存在model对象的params属性中。
model.rsquared
a = model.params def reg_func(a,x,y):
f6 = a[6] * np.sqrt(y)
f5 = a[5] * np.sin(x)
f4 = a[4] * y ** 2
f3 = a[3] * x ** 2
f2 = a[2] * y
f1 = a[1] * x
f0 = a[0] * 1
return (f6+f5+f4+f3+f2+f1+f0)
# reg_func返回给定最优回归参数和自变量数据点的回归函数值
RZ = reg_func(a,X,Y)
fig = plt.figure(figsize=(9, 6))
ax = fig.gca(projection='3d')
surf1 = ax.plot_surface(X,Y,Z,rstride=2,cstride=2,cmap=mpl.cm.coolwarm,linewidth=0.5,antialiased=True)
surf2 = ax.plot_wireframe(X,Y,RZ,rstride=2,cstride=2,label = 'regression')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('f(x,y)')
fig.colorbar(surf,shrink=0.5,aspect=5)

最新文章

  1. vs发布的程序不依赖运行时库msvcp100.dll
  2. 【C语言入门教程】3.4 循环控制语句
  3. 李洪强漫谈iOS开发[C语言]-045-循环结构
  4. 问题:如何在固定大小的DIV层插入N多个图片
  5. CUDA基本概念
  6. Enter password: ERROR 2002 (HY000): Can't connect to local MySQL server through socket '/var/lib/mysql/mysql.sock' (2)
  7. java转义xml中的多余尖括号
  8. [HTML5实现人工智能]小游戏《井字棋》发布,据说IQ上200才能赢
  9. STM32学习笔记(一)——点亮一个LED
  10. java常见面试题(二)
  11. Leetcode解题-链表(2.2.2)ReverseLinkedList
  12. 华为Java机试题
  13. Fragment调用startActivityForResult导致的回调Activity无法获取正确的requestId的问题
  14. @RequestParam 和@RequestBody 的区别?
  15. P2661 信息传递 二分图的最小环
  16. 深入JVM之类的加载过程
  17. 微信公众号 access_token 没有过期 却失效
  18. Linux内存压力测试stressapptest
  19. YAML 语言教程(转载)
  20. Redis:ERR operation not permitted

热门文章

  1. 使用OpenLDAP部署目录服务
  2. day-11函数的形参与实参
  3. 树莓派3 Raspberry系统安装samba
  4. 谷歌浏览器可以google了
  5. 深度系统 deepin 15.9 关闭桌面
  6. Java 异常: SimpleDateFormat java.lang.NumberFormatException: multiple points
  7. Java中级开发工程师知识点归纳
  8. Java中所涉及到的设计模式小记
  9. VMware 打开虚拟机的时候提示 internal error 内部错误 遇到这个问题时我的解决方法
  10. Kong(V1.0.2) Clustering Reference