MATLAB 数据拟合 (使用 polyfit 多项式曲线拟合、polyval)

本文主要是介绍MATLAB 数据拟合 (使用 polyfit 多项式曲线拟合、polyval),希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

解决数据拟合问题最重要方法是最小二乘法和回归分析。如,我们需要从一组测定的数据(例如N个点(xi,yi)(i=0,1,…,m))去求得自变量 x 和因变量 y 的一个近似解表达式 y=f(x),这就是由给定的 N 个点(xi,yi)(i=0,1,…,m)求数据拟合的问题。(注意数据拟合和数据插值是不同的,举个例子:因为测量数据往往不可避免地带有测试误差,而插值多项式又通过所有的点(xi,yi),这样就使插值多项式保留了这些误差,从而影响逼近精度,使得插值效果不理想)
所以使用最小二乘法曲线拟合法:即寻求已知函数的一个逼近函数y=f(x),使得逼近函数从总体与已知函数的偏差按某种方法度量能达到最小,而又不一定通过全部的点(xi,yi),这个时候就需要使用最小二乘法曲线拟合法。
数据拟合的具体做法是:对给定的数据(xi,yi)(i=0,1,…,m),在取定的函数类 ϕ \phi ϕ中使误差 r i = p ( x i ) − y i ( i = 0 , 1 , … , m ) r_{i}=p\left(x_{i}\right)-y_{i}(i=0,1, \ldots, m) ri=p(xi)yi(i=0,1,,m)的平方和最小,即
[ ∑ i = 0 m r i 2 = ∑ i = 0 m [ p ( x i − y i ) ] 2 ] min ⁡ \left[\sum_{i=0}^{m} r_{i}^{2}=\sum_{i=0}^{m}\left[p\left(x_{i}-y_{i}\right)\right]^{2}\right]_{\min } [i=0mri2=i=0m[p(xiyi)]2]min
从几何意义讲,即寻求与给定点 x i − y i ( i = 0 , 1 , … , m ) x_{i}-y_{i}(i=0,1, \ldots, m) xiyi(i=0,1,,m) 的距离平方和为最小的曲线y=p(x)。函数p(x)称为拟合函数或最小二乘解,求拟合函数p(x)的方法称为曲线拟合的最小二乘法。
在曲线拟合中,函数类 ϕ \phi ϕ 可有不同的选取方法。
MATLAB工具箱中提供了最小二乘拟合函数 polyfit() -->多项式曲线拟合
具体调用格式有三种:

  1. P = polyfit(X,Y,N)
  2. [P,S] = polyfit(X,Y,N)
  3. [P,S,MU] = polyfit(X,Y,N)

(1)P = polyfit(X,Y,N) 返回次数为 n 的多项式 p(x) 的系数,该阶数是 y 中数据的最佳拟合(在最小二乘方式中)。p 中的系数按降幂排列,p 的长度为 n+1.
其中 X 为输入的向量x,Y 为得到的函数值,N 表示拟合的最高次数,返回的P值为拟合的多项式:

P ( 1 ) × X n + P ( 2 ) × X n − 1 + … + P ( N ) × X + P ( N + 1 ) P(1) \times X^{n}+P(2) \times X^{n-1}+\ldots+P(N) \times X+P(N+1) P(1)×Xn+P(2)×Xn1++P(N)×X+P(N+1)

(2) [P,S] = polyfit(X,Y,N)还返回一个结构体 S,S可用作 polyval 的输入来获取误差估计值(S为由范德蒙矩阵的QR分解的R分量)。
其中 X 为输入的向量x,Y 为得到的函数值,N 表示拟合的最高次数,返回的P值为拟合的多项式:

P ( 1 ) × X n + P ( 2 ) × X n − 1 + … + P ( N ) × X + P ( N + 1 ) P(1) \times X^{n}+P(2) \times X^{n-1}+\ldots+P(N) \times X+P(N+1) P(1)×Xn+P(2)×Xn1++P(N)×X+P(N+1)

(3)[P,S,MU] = polyfit(X,Y,N)还返回 mu,mu是一个二元素向量,包含中心化值和缩放值。mu(1) 是 mean(x),mu(2) 是 std(x)。使用这些值时,polyfit 将 x 的中心置于零值处并缩放为具有单位标准差:

x ^ = x − x ˉ σ x \hat{x}=\frac{x-\bar{x}}{\sigma_{x}} x^=σxxxˉ

这种中心化和缩放变换可同时改善多项式和拟合算法的数值属性。

下面来看一些具体的例子:(来源于帮助文档,改编)

  • 将多项式与三角函数拟合:
    在区间 [0,4*pi] 中沿余弦曲线生成 10 个等间距的点。
>> x = linspace(0,4*pi,10);
>> y = cos(x);

使用 polyfit 将一个 7 次多项式与这些点拟合。

>> p = polyfit(x,y,7);

在更精细的网格上计算多项式并绘制结果图。

>> x1 = linspace(0,4*pi);
>> y1 = polyval(p,x1);
>> figure
>> plot(x,y,'o')
>> hold on
>> plot(x1,y1,'m')
>> hold off

运行结果如下:
在这里插入图片描述

  • 将多项式与点集拟合
    创建一个由区间 [0,1] 中的 5 个等间距点组成的向量,并计算这些点处的 y(x)= (2+x)^-1 。
>> x=0:0.2:1;   //或者写成  >> x=0:0.2:1;  也可以实现相同的作用
>> y = 1./(2+x);

将 4 次多项式与 5 个点拟合。通常,对于 n 个点,可以拟合 n-1 次多项式以便完全通过这些点。

>> p = polyfit(x,y,4);

在由 0 和 2 之间的点组成的更精细网格上计算原始函数和多项式拟合。

>> x1 = linspace(0,2);
>> y1 = 1./(2+x1);
>> f1 = polyval(p,x1);

在更大的区间 [0,2] 中绘制函数值和多项式拟合,其中包含用于获取以圆形突出显示的多项式拟合的点。多项式拟合在原始 [0,1] 区间中的效果较好,但在该区间外部很快与拟合函数出现差异。

>> figure
>> plot(x,y,'o')
>> hold on
>> plot(x1,y1)
>> plot(x1,f1,'b--')
>> legend('y','y1','f1')

运行结果如下:
在这里插入图片描述

  • 对误差函数进行多项式拟合:
    首先生成 x 点的向量,在区间 [0,1.5*pi] 内等间距分布;然后计算这些点处的 erf(x)。
    注:erf 是误差函数,也称高斯误差函数。
x = (0:0.10:1.5*pi)';y=erf(x);

确定6次逼近多项式的系数。

p=polyfit(x,y,6)

p =

0.0019   -0.0259    0.1242   -0.1785   -0.3442    1.2684   -0.0095

为了查看拟合情况如何,在各数据点处计算多项式,并生成说明数据、拟合和误差的一个表。

f=polyval(p,x);
T=table(x,y,f,y-f,'VariableNames',{'X','Y','Fit','FitError'})

T =

48×4 table

 X        Y          Fit         FitError  
___    _______    __________    ___________0          0    -0.0094811      0.0094811
0.1    0.11246       0.11375     -0.0012876
0.2     0.2227       0.22919     -0.0064911
0.3    0.32863       0.33619     -0.0075599
0.4    0.42839       0.43431      -0.005914
0.5     0.5205       0.52334     -0.0028408
0.6    0.60386       0.60326     0.00059392
0.7     0.6778        0.6742      0.0035981
0.8     0.7421       0.73643      0.0056695
0.9    0.79691       0.79033      0.00657941     0.8427       0.83637      0.0063316
1.1    0.88021        0.8751      0.0051058
1.2    0.91031       0.90712       0.003194
1.3    0.93401       0.93307     0.00093826
1.4    0.95229       0.95361      -0.001323
1.5    0.96611        0.9694     -0.0032971
...

在该区间中,插值与实际值非常符合。创建一个绘图,以显示在该区间以外,外插值与实际数据值如何快速偏离。

 x1=(0:0.10:2*pi)';y1=erf(x1);f1=polyval(p,x1);plot(x,y,'o')plot(x1,y1,'-')hold onplot(x1,y1,'-')plot(x1,f1,'g-')axis([0  5.5  0  3])

在这里插入图片描述

  • 使用中心化和缩放改善数值属性
    创建一个由 1750 - 2000 年的人口数据组成的表,并绘制数据点。
 year = (1750:25:2000)';pop = 1e6*[791 856 978 1050 1262 1544 1650 2532 6122 8170 11560]';T = table(year, pop)

T =

11×2 table

year       pop   
____    _________1750     7.91e+08
1775     8.56e+08
1800     9.78e+08
1825     1.05e+09
1850    1.262e+09
1875    1.544e+09
1900     1.65e+09
1925    2.532e+09
1950    6.122e+09
1975     8.17e+09
2000    1.156e+10
 plot(year,pop,'square')

在这里插入图片描述
使用带三个输入的 polyfit 拟合一个使用中心化和缩放的 5 次多项式,这将改善问题的数值属性。polyfit 将 year 中的数据以 0 为进行中心化,并缩放为具有标准差 1,这可避免在拟合计算中出现病态的 Vandermonde (范德蒙)矩阵。

[p,~,mu] = polyfit(T.year, T.pop, 5);

使用带四个输入的 polyval,根据缩放后的年份 (year-mu(1))/mu(2) 计算 p。绘制结果对原始年份的图。

f = polyval(p,year,[],mu);
hold on
plot(year,f)

在这里插入图片描述

  • 简单线性回归
    将一个简单线性回归模型与一组离散二维数据点拟合。
    创建几个由样本数据点 (x,y) 组成的向量。对数据进行一次多项式拟合。
>> x=1:40;
>> y=0.4*x-1.5*randn(1,60);
注意:这里 矩阵维度必须一致。 所以第二行代码的输入是错误的y=0.4*x-1.5*randn(1,40);p=polyfit(x,y,1);

计算在 x 中的点处拟合的多项式 p。用这些数据绘制得到的线性回归模型。

 f = polyval(p,x); plot(x,y,'o',x,f,'-')legend('data','linear fit')

在这里插入图片描述

  • 具有误差估计值的线性回归
    将一个线性模型拟合到一组数据点并绘制结果,其中包含预测区间为 95% 的估计值。
    创建几个由样本数据点 (x,y) 组成的向量。使用 polyfit 对数据进行一次多项式拟合。指定两个输出以返回线性拟合的系数以及误差估计结构体。
x = 1:100; 
y = -0.3*x + 2*randn(1,100); 
[p,S] = polyfit(x,y,1);

计算以 p 为系数的一次多项式在 x 中各点处的拟合值。将误差估计结构体指定为第三个输入,以便 polyval 计算标准误差的估计值。标准误差估计值在 delta 中返回。

[y_fit,delta] = polyval(p,x,S);

绘制原始数据、线性拟合和 95% 预测区间 y ± 2 Δ y \pm 2 \Delta y±2Δ

plot(x,y,'go')
hold on
plot(x,y_fit,'r-')
plot(x,y_fit+2*delta,'m--',x,y_fit-2*delta,'m--')
title('Linear Fit of Data with 95% Prediction Interval')
legend('Data','Linear Fit','95% Prediction Interval')

在这里插入图片描述

  • 输入参数解释:
  • x – 查询点 (向量)
    查询点,指定为一个向量。x 中的点对应于 y 中包含的拟合函数值。如果 x 不是向量,则 polyfit 将其转换为列向量 x( : )。😃
    x 具有重复(或接近重复)的点或者如果 x 可能需要中心化和缩放时的警告消息结果。
    **数据类型:**single | double
    复数支持:
  • y – 查询点位置的拟合值 (向量)
    查询点位置的拟合值,指定为向量。y 中的值对应于 x 中包含的查询点。如果 y 不是向量,则 polyfit 将其转换为列向量 y ( : ) 。😃
    **数据类型:**single | double
    复数支持:
  • n --多项式拟合的次数 (正整数标量)
    多项式拟合的次数,指定为正整数标量。n 指定 p 中最左侧系数的多项式幂
  • 输出参数解释:
  • p – 最小二乘拟合多项式系数(向量)

最小二乘拟合多项式系数,以向量的形式返回。p 的长度为 n+1,包含按降幂排列的多项式系数,最高幂为 n。如果 x 或 y 包含 NaN 值且 n < length(x),则 p 的所有元素均为 NaN。

使用 polyval 计算 p 在查询点处的解。

  • S – 误差估计结构体 (结构体)
    误差估计结构体。此可选输出结构体主要用作 polyval 函数的输入,以获取误差估计值。S 包含以下字段:
字段说明
RVandermonde 矩阵 x 的 QR 分解的三角因子
df自由度
normr残差的范数

如果 y 中的数据是随机的,则 p 的估计协方差矩阵是 (Rinv*Rinv’)*normr^2/df,其中 Rinv 是 R 的逆矩阵。

如果 y 中数据的误差呈独立正态分布,并具有常量方差,则 [y,delta] = polyval(…) 可生成至少包含 50% 的预测值的误差边界。即 y ± delta 至少包含 50% 对 x 处的未来观测值的预测值。

  • mu – 中心化值和缩放值

中心化值和缩放值,以二元素向量形式返回。mu(1) 为 mean(x),mu(2) 为 std(x)。这些值以单位标准差将 x 中的查询点的中心置于零值处。

使用 mu 作为 polyval 的第四个输入以计算 p 在缩放点 (x - mu(1))/mu(2) 处的解。

局限性

  • 在涉及很多点的问题中,使用 polyfit 增加多项式拟合的次数并不总能得到较好的拟合。高次多项式可以在数据点之间振动,导致与数据之间的拟合较差。在这些情况下,可使用低次多项式拟合(点之间倾向于更平滑)或不同的方法,具体取决于该问题。
  • 多项式在本质上是无边界的振荡函数。所以,它们并不非常适合外插有界的数据或单调(递增或递减)的数据。

算法:
polyfit 使用 x 构造具有 n+1 列和 m = length(x) 行的 Vandermonde 矩阵 V 并生成线性方程组
( x 1 n x 1 n − 1 … 1 n 2 x 2 − 1 … 1 ⋮ ⋮ ⋱ ⋮ n m x m − 1 … 1 ) ( p 1 p 2 ⋮ p n + 1 ) = ( y 1 y 2 ⋮ y m ) \left(\begin{array}{cccc} x_{1}^{n} & x_{1}^{n-1} & \ldots & 1 \\ n_{2} & x_{2}-1 & \ldots & 1 \\ \vdots & \vdots & \ddots & \vdots \\ n_{m} & x_{m}-1 & \ldots & 1 \end{array}\right)\left(\begin{array}{c} p_{1} \\ p_{2} \\ \vdots \\ p_{n+1} \end{array}\right)=\left(\begin{array}{c} y_{1} \\ y_{2} \\ \vdots \\ y_{m} \end{array}\right) x1nn2nmx1n1x21xm1111p1p2pn+1=y1y2ym

其中 polyfit 使用 p = V\y 求解。由于 Vandermonde 矩阵中的列是向量 x 的幂,因此条件数 V 对于高阶拟合来说通常较大,生成一个奇异系数矩阵。在这些情况下,中心化和缩放可改善系统的数值属性以产生更可靠的拟合。

这篇关于MATLAB 数据拟合 (使用 polyfit 多项式曲线拟合、polyval)的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/948027

相关文章

解决Maven项目idea找不到本地仓库jar包问题以及使用mvn install:install-file

《解决Maven项目idea找不到本地仓库jar包问题以及使用mvninstall:install-file》:本文主要介绍解决Maven项目idea找不到本地仓库jar包问题以及使用mvnin... 目录Maven项目idea找不到本地仓库jar包以及使用mvn install:install-file基

Spring 请求之传递 JSON 数据的操作方法

《Spring请求之传递JSON数据的操作方法》JSON就是一种数据格式,有自己的格式和语法,使用文本表示一个对象或数组的信息,因此JSON本质是字符串,主要负责在不同的语言中数据传递和交换,这... 目录jsON 概念JSON 语法JSON 的语法JSON 的两种结构JSON 字符串和 Java 对象互转

Python使用getopt处理命令行参数示例解析(最佳实践)

《Python使用getopt处理命令行参数示例解析(最佳实践)》getopt模块是Python标准库中一个简单但强大的命令行参数处理工具,它特别适合那些需要快速实现基本命令行参数解析的场景,或者需要... 目录为什么需要处理命令行参数?getopt模块基础实际应用示例与其他参数处理方式的比较常见问http

C 语言中enum枚举的定义和使用小结

《C语言中enum枚举的定义和使用小结》在C语言里,enum(枚举)是一种用户自定义的数据类型,它能够让你创建一组具名的整数常量,下面我会从定义、使用、特性等方面详细介绍enum,感兴趣的朋友一起看... 目录1、引言2、基本定义3、定义枚举变量4、自定义枚举常量的值5、枚举与switch语句结合使用6、枚

使用Python从PPT文档中提取图片和图片信息(如坐标、宽度和高度等)

《使用Python从PPT文档中提取图片和图片信息(如坐标、宽度和高度等)》PPT是一种高效的信息展示工具,广泛应用于教育、商务和设计等多个领域,PPT文档中常常包含丰富的图片内容,这些图片不仅提升了... 目录一、引言二、环境与工具三、python 提取PPT背景图片3.1 提取幻灯片背景图片3.2 提取

C++如何通过Qt反射机制实现数据类序列化

《C++如何通过Qt反射机制实现数据类序列化》在C++工程中经常需要使用数据类,并对数据类进行存储、打印、调试等操作,所以本文就来聊聊C++如何通过Qt反射机制实现数据类序列化吧... 目录设计预期设计思路代码实现使用方法在 C++ 工程中经常需要使用数据类,并对数据类进行存储、打印、调试等操作。由于数据类

使用Python实现图像LBP特征提取的操作方法

《使用Python实现图像LBP特征提取的操作方法》LBP特征叫做局部二值模式,常用于纹理特征提取,并在纹理分类中具有较强的区分能力,本文给大家介绍了如何使用Python实现图像LBP特征提取的操作方... 目录一、LBP特征介绍二、LBP特征描述三、一些改进版本的LBP1.圆形LBP算子2.旋转不变的LB

Maven的使用和配置国内源的保姆级教程

《Maven的使用和配置国内源的保姆级教程》Maven是⼀个项目管理工具,基于POM(ProjectObjectModel,项目对象模型)的概念,Maven可以通过一小段描述信息来管理项目的构建,报告... 目录1. 什么是Maven?2.创建⼀个Maven项目3.Maven 核心功能4.使用Maven H

Python中__init__方法使用的深度解析

《Python中__init__方法使用的深度解析》在Python的面向对象编程(OOP)体系中,__init__方法如同建造房屋时的奠基仪式——它定义了对象诞生时的初始状态,下面我们就来深入了解下_... 目录一、__init__的基因图谱二、初始化过程的魔法时刻继承链中的初始化顺序self参数的奥秘默认

SpringBoot使用GZIP压缩反回数据问题

《SpringBoot使用GZIP压缩反回数据问题》:本文主要介绍SpringBoot使用GZIP压缩反回数据问题,具有很好的参考价值,希望对大家有所帮助,如有错误或未考虑完全的地方,望不吝赐教... 目录SpringBoot使用GZIP压缩反回数据1、初识gzip2、gzip是什么,可以干什么?3、Spr