【数学建模】MATLAB入门教程:插值与拟合(下)

前言

插值与拟合在数据处理和科学计算中扮演着非常重要的角色,它们用于估算未知数据点的值,帮助我们理解和预测数据趋势

一、一维插值

1、一维插值定义

已知n+1个节点($x_j$,$y_j$)(j=0,1,...,n,其中$x_j$互不相同,不妨设a=$x_0$<$x_1$<...<$x_n$=b),求任一插值点$x^*$(!=$x_j$)处的插值$y^*$

解决方法:构造一个相对简单的函数y=f(x),通过全部节点,即f($x_j$)=$y_j$(j=0,1,...,n)

再用f(x)计算插值,即$y^*$=f($x^*$)

2、常见方法

1、拉格朗日插值:已知函数f(x)在n+1个节点$x_0$$x_1$,···,$x_n$处的函数值为$y_0$$y_1$,···,$y_n$.求一n次多项式$p_n$(x),使其满足:$p_n$($x_i$)=$y_i$,i=0,1,...,n.解决此问题的拉格朗日插值多项式公式如下:

 特别地:两点一次(线性)插值多项式

例题        g(x)=1/(1+x^2),-5<=x<5.采用拉格朗日多项式插值:选取不同插值节点n+1个,其中n为插值多项式的次数,当n分别区2,4,6,8,10时,绘出插值结果图形

function lagrange_interpolation_example% 定义函数 g(x)g = @(x) 1 ./ (1 + x.^2);% 定义插值区间和次数x_interval = [-5, 5];n_values = [2, 4, 6, 8, 10];% 初始化图形figure;hold on;% 绘制原始函数x_fine = linspace(x_interval(1), x_interval(2), 1000);y_fine = g(x_fine);plot(x_fine, y_fine, 'b-', 'LineWidth', 2);% 对于每个n值,计算并绘制拉格朗日插值for n = n_values% 在区间上均匀选择n+1个插值节点x_nodes = linspace(x_interval(1), x_interval(2), n + 1);y_nodes = g(x_nodes);% 计算拉格朗日插值多项式并绘制结果[x_lagrange, y_lagrange] = lagrange_poly(x_fine, x_nodes, y_nodes);plot(x_lagrange, y_lagrange, 'DisplayName', sprintf('n=%d', n), 'LineWidth', 1.5);end% 添加图例、标签和标题legend('Original Function', sprintf('Lagrange Interpolation, n=%d', n_values));xlabel('x');ylabel('g(x)');title('Lagrange Interpolation of g(x) = 1 / (1 + x^2)');grid on;hold off;
endfunction [x_lagrange, y_lagrange] = lagrange_poly(x, x_nodes, y_nodes)% 拉格朗日插值多项式计算函数n = length(x_nodes) - 1;y_lagrange = zeros(size(x));for i = 1:n+1% 计算拉格朗日基函数L_i = 1;for j = 1:n+1if j ~= iL_i = L_i * (x - x_nodes(j)) / (x_nodes(i) - x_nodes(j));endend% 累加基函数与对应函数值的乘积y_lagrange = y_lagrange + L_i * y_nodes(i);end% 返回插值点x和对应的插值结果y_lagrangex_lagrange = x;
end% 调用函数以运行示例
lagrange_interpolation_example

2、分段线性插值:将每两个相邻的点用直线连接起来,如此形成的一条折线就是分段线性插值函数$I_n$(x),它满足$I_n$($x_i$)=$y_i$

计算量与n无关,n越大,误差越小

 例题        g(x)=1/(1+x^2),-6<=x<6,采用分段线性插值法求插值

x=linspace(-6,6,100);
y=1./(1+x.^2);
x1=linspace(-6,6,5);
y1=1./(1+x1.^2);
plot(x,y,x1,y1,x1,y1,'o','LineWidth',1.5),
gtext('n=4'),

3、样条插值:比分段插值更加光滑,在数学上,光滑程度的定量描述是:函数(曲线)的k阶导数存在且连续,则称该曲线具有k接光滑性/

光滑性的阶次越高,则越光滑。而三次样条插值就是一个较低次的分段多项式达到较高阶光滑性的典型例子!

3、用MATLAB作插值计算(重点)

一维插值函数:

yi=interp1(x,y,xi,'method')

yi表示xi处的插值结果

x、y表示插值节点

xi表示被插值点

'method'表示插值方法:

'nearest'      最近临近插值  

'linear'        线性插值

'spilne'        三次样条插值

'cubic'        立方插值

缺省时默认采用分段线性插值

例题    从1点到12点的11个小时内,每个1小时测量一次温度,测得的温度数值依次为5,8,9,15,25,29,31,30,22,25,27,24,试估计每隔1/10消失的温度值

x=1:12;
y=[5 8 9 15 25 29 31 30 22 25 27 24];
xi=1:0.1:12;
yi=interp1(x,y,xi,'spline');
plot(x,y,'+',xi,yi,x,y,'r:')
xlabel('小时'),ylabel('温度')

例  在车辆行驶中,从驾驶员看到障碍物开始,到作出判断而采取制动措施停车所需的最短距离叫停车视距,停车视距由三部分组成:一是驾驶员反应时间内行驶的距离(即反应距离);三是开始制动到车辆完全停止所行驶的距离(即制动距离);三是车辆停止时与障碍物应该保持的安全距离。其中,制动距离主要与行驶速度和路面类型有关。根据测试,某型年辆在潮湿天气于沥青路面行驶时,其行车速度(单位:km/h) 与制动距离(单位:m)的关系如下表所示. 

速度2030405060708090100110120130140150
制动距离3.157.0812.5919.6828.3438.5750.463.7578.7195.22113.29132.93154.12176.87

假设驾驶员的反应时间为10s,安全距离为10m。请问(1)根据某驾驶员的实际视力和视觉习惯,其驾驶时的有效视距为120m,则其在该路面行车时,时速最高不能超过多少?(结果取整)(2)若以表中数据参考,设计一条最高时速为125km/h的高速公路,则设计人员应该保证驾驶者在公路上任一点的可视距离为多少米?

我们设速度为v,停车视距为d,反应距离为d1,制动距离为d2,安全距离为d3,已知反应时间为10s,安全距离为10m,则:10v+d2+10=120

v=20:10:150;
vs=v.*(1000/3600);
d1=10.*vs;
d2=[3.15,7.08,12.59,19.68,28.34,38.57,50.4,63.75,78.71,95.22,113.29,132.93,154.12,176.87];
d3=10;
d=d1+d2+d3;
vi=20:1:150;
di=interp1(v,d,vi,'spline');
x=abs(di-120);
[y,i]=sort(x);
vi(i(1))
plot(vi,di,vi(i(1)),di(i(1)),'rp')
j=find(vi==125);
di(j)

二、二位插值 

 二位插值我们就不再细讲了,讲一下如何用matlab来计算即可

用MATLAB作网格节点数据的插值

 z=interp2(x0,y0,z0,x,y,'method')

其中x0、y0、z0表示插值节点

x、y表示被插值点

z表示被插值点的函数值

'method'表示插值方法

例题 测得在平板表面3*5网格点处的温度分别为:

82        81        80        82        84

79        63        61        65        81

84        84        82        85        86 

 试作出平板表面的温度分布曲面z=f(x,y)的图形

x=1:5;
y=1:3;
temps=[82 81 80 82 84;79 63 61 65 81;84 84 82 85 86];
mesh(x,y,temps)
pause
xi=1:0.2:5;
yi=1:0.2:3;
zi=interp2(x,y,temps,xi,yi,'cubic');
mesh(xi,yi,zi)

数据导入 

 xlsread('filename',n,'al:bl')

数据残缺、异常数据

新建实时脚本 - 导入数据 - 任务 - 清理缺失(清理离群) 

三、拟合(最小二乘法)

拟合问题引例 已知热敏电阻数据:

温度t(摄氏度)20.532.751.073.095.7
电阻R(\Omega)7658268739421032

求60摄氏度时的电阻

设R=at+b,a,b为待定系数

曲线拟合问题的提法:已知一组(二维数据),即平面上n个点($x_i$,$y_i$)i,···,你,寻求一个函数(曲线)y=f(x),使f(x)在某种准则下与所有数据点最为接近,即曲线拟合的最好。

插值和拟合的关系:

问题:给定一批数据点,需要确定满足特定要求的曲线或曲面

解决方案:

·若要求所求曲线(面)通过所给特定数据点,就是插值问题;

·若不要求曲线(面)通过所有数据点,而是要求它反映对象整体的变化趋势,这就是数据拟合,又称曲线拟合或曲面拟合。

二者区别:

·插值试图去通过已知点了解未知点处的函数值;

·拟合则在于整体上用某种已知函数去拟合数据点列所在未知函数的状态。

·关键区别在于插值要求必须经过已知点列,拟合只求尽量靠近不必经过!

用MATLAB解决线性拟合问题

线性最小二乘法

 polyfit(x,y,m)

x、y表示输入同长度的数组

m表示拟合多项式次数

a表示输出拟合多项式系数

a=[$a_1$,...,$a_m$,$a_{m+1}$]

z=polyval(a,x)来算多项式在x处的值

非线性最小二乘拟合 

MATLAB提供了两个求非线性最小二乘拟合的函数:lsqcurvefit和lsqnonlin。两个命令都要先建立M文件fun.m,在其中定义函数f(x),但两者定义f(x)的方式是不同的。

1、lsqcurvefit('fit',x0,xdata,ydata,options);

fun是一个实现建立的定义函数F(x,xdata)的M文件,自变量为x和xdata

x0表示迭代初值

options表示选项见无约束优化

2、lsqnonlin('fun',x0,options);

同上 

MATLAB解应用问题实例

1、电阻问题

例  由数据

温度t(摄氏度)20.532.751.073.095.7
电阻R(\Omega)7658268739421032

来拟合y=ax+b;用命令:polyfit(x,y,m)

t=[20.5 32.5 51 73 95.7];
r=[765 826 873 942 1032];
n=polyfit(t,r,1);
a=n(1)
b=n(2)
y=polyval(n,t);
plot(t,r,'k+',t,y,'r')

2、水塔流量估计问题 

某居民区有一供居民用水的圆柱形水塔,一般可以通过测量其水位来估计水的流量,但面临的困难是,当水塔水位下降到设定的最低水位时,水泵自动启动向水塔供水,到设定的最高水位时停止供水,这段时间无法测量水塔的水位和水泵的供水量.通常水泵每天供水一两次,每次约两小时.
水塔是一个高12.2m,直径17.4m的正圆柱.按照设计,水塔水位降至约8.2m时,水泵自动启动,水位升到约10.8m时水泵停止工作.
表1是某一天的水位测量记录,试估计任何时刻(包括水泵正供水时)从水塔流出的水流量,及一天的总用水量.(符号//表示水泵启动)

时刻(h)00.921.842.953.874.985.907.017.938.97

水位(cm)

968948931913898881869852839822
时刻(h)9.9810.9210.9512.0312.9513.8814.9815.9016.8317.93
水位(cm)////108210501021994965941918892
时刻(h)19.0419.9620.8422.0122.9623.8824.9925.91
水位(cm)866843822////105910351018

从测量记录看,一天有两个供水时段(以下称第1供水时段和第2供水时段),和3个水泵不工作时段(以下称第1时段 =0到 =8.97,第2次时段 =10.95到 =20.84和第3时段/23以后),对第1、2时段的测量数据直接分别作多项式拟合,得到水位函数.为使拟合曲线比较光滑,多项式次数不要太高,一般在3~6. 由于第3时段只有3个测量记录,无法对这一时段的水位作出较好的拟合

确定流量:对于第1、2时段只需将水位函数求导数即可,对于两个供水时段的流量,则用供水时段前后水泵(不工作时段)的流量拟合得到,并且将拟合得到的第2供水时段流量外推,将第3时段流量包含在第2供水时段内. 

一天总用水量的估计:总用水量等于两个水泵不工作时段和两个时段用水量之和,他们都可以由流量对时间积分得到

这个问题是相对比较复杂的,我们这里只给出结果,计算就交给大家

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mzph.cn/diannao/26767.shtml

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈email:809451989@qq.com,一经查实,立即删除!

相关文章

鸿蒙开发:任务(Mission)与启动模式

任务&#xff08;Mission&#xff09;与启动模式 如前文所述&#xff0c;一个UIAbility实例对应一个任务。UIAbility实例个数与UIAbility配置的启动模式有关。在FA模型下&#xff0c;通过config.json配置文件中的“launchType”属性配置&#xff1b;在Stage模型下&#xff0c;…

C#聊天室②

客户端 桌面 MyClient client;public Form1(){InitializeComponent();}// 进入聊天室按钮方法private void button1_Click(object sender, EventArgs e){if (!string.IsNullOrEmpty(textBox1.Text)){// 开始连接服务器 封装一个自定义客户端类client new MyClient(); // 给cl…

基于jeecgboot-vue3的Flowable流程-流程处理(一)

因为这个项目license问题无法开源&#xff0c;更多技术支持与服务请加入我的知识星球。 这部分修正一些流程处理中VForm3线上的一些bug问题 1、初始化流程提交与现实的前端页面代码 <!--初始化流程加载默认VForm3表单信息--><el-col :span"16" :offset&qu…

nlp学习笔记

目录 很多入门例子 bert chinese 很多入门例子 https://github.com/lansinuote/Huggingface_Toturials bert chinese import torch import torch.nn as nn from transformers import AutoTokenizer, AutoModel, BertModel, TFBertModel, BertTokenizer# youpath = D:/bert-…

【Mybatis】动态SQL的绑定和公共sql语句片段

Mybatis还有三个标签&#xff0c;分别是bind&#xff0c;sql和include ①bind&#xff1a;这个标签作用就是将OGNL标签里的值&#xff0c;进行二次加工&#xff0c;在绑定到另一个变量里&#xff0c;供其他标签使用&#xff0c;举个例子 调用getUsers方法的时候&#xff0c;我…

网络编程(二)TCP

一、TCP网络编程 网络编程模型&#xff1a; C/S模型&#xff1a;客户端服务器模型 优点&#xff1a; 客户端可以缓存一些数据&#xff0c;使用时直接在本地读取&#xff0c;无需每次重新下载&#xff1b; 由于客户端和服务器都是自己开发的&#xff0c;可以自定义协议 缺点&a…

文案提取小帮手轻松将视频为转文字!而且不限时长

作为一个自媒体的资深用户总在一个一个的敲字真的太慢了&#xff0c;而且很多创作者都知道追热点是和时间赛跑。如果你嫌弃自己手抄效率太低&#xff0c;看视频又嫌时间太长。 今天叫教你一个可以将视频转文字的工具&#xff0c; 这个工具就叫文案提取小帮手&#xff0c;而且…

Vue前端ffmpeg压缩视频再上传(全网唯一公开真正实现)

1.Vue项目中安装插件ffmpeg 1.1 插件版本依赖配置 两个插件的版本 "ffmpeg/core": "^0.10.0", "ffmpeg/ffmpeg": "^0.10.1"package.json 和 package-lock.json 都加入如下ffmpeg的版本配置&#xff1a; 1.2 把ffmpeg安装到项目依…

Java老人护理上门服务类型系统小程序APP源码

&#x1f338; 老人上门护理服务系统&#xff1a;温暖与专业并存 &#x1f338; 一、&#x1f3e0; 走进老人上门护理服务系统 随着社会的快速发展&#xff0c;我们越来越关注老年人的生活质量。老人上门护理服务系统应运而生&#xff0c;它结合了现代科技与人性化服务&#…

最小生成树prim算法详解

prim算法解决的是最小生成树问题&#xff0c;即在一个给定的无向图G中求一棵生成树T&#xff0c;使得这棵树拥有图G中的所有顶点&#xff0c;且所有边都是来自图G中的边&#xff0c;并且满足整棵树的边权之和最小。 prim算法的基本思想是对图G设置集合S来存放已被访问的顶点&a…

【源码】【Spring+SpringMVC+MyBatis】电子商城网上购物平台的设计与开发

学生成绩管理系统 系统功能开发环境开发技术前端技术后端技术 系统展示登录界面注册界面系统首页商品详情页下单界面付款界面购物车界面 源码获取↓↓↓↓&#xff1a; 源码可在后台私信联系博主或文末添加博主微信获取帮助 系统功能 登录、注册模块&#xff1a;如果用户第一次…

基于Java+Swing+mysql幼儿园信息管理系统V2

博主介绍&#xff1a; 大家好&#xff0c;本人精通Java、Python、C#、C、C编程语言&#xff0c;同时也熟练掌握微信小程序、Php和Android等技术&#xff0c;能够为大家提供全方位的技术支持和交流。 我有丰富的成品Java、Python、C#毕设项目经验&#xff0c;能够为学生提供各类…

【网络编程】基于UDP的服务器端/客户端

UDP可看作是信件邮寄&#xff0c;邮寄过程可能会信件丢失&#xff0c;是一种不可靠的数据传输服务。 但UDP性能更高&#xff0c;实现更加简洁。流控制是区分UDP和TCP的最重要标志。 IP的作用就是让离开主机B的UDP数据包传递给主机B&#xff0c;UDP根据端口号将传到主机的数据包…

Python - OS模块+sys模块

一、OS模块基本用法&#xff1a; import osprint(os.getcwd()) # 获取当前工作目录os.chdir(data) # 改变当前脚本工作目录&#xff1b;相当于终端里面的cd命令 print(os.getcwd()) # 获取当前工作目录 运行结果&#xff1a; D:\__TC22008_VBF\FOTA-vFlash-AutoTest D:\__TC22…

conda虚拟环境报错总结

创建conda虚拟环境 文章前景&#xff08;小白篇&#xff09;为什么要安装Anaconda&#xff1f;&#xff1f;&#xff1f; Conda创建虚拟环境遇到的错误总结错误1&#xff1a;jupyter 里面没有显示我的虚拟环境怎么办&#xff1f;错误2&#xff1a;配置pycharm的时候conda虚拟环…

破解发展难题 台山这家合作社以农业社会化服务助推乡村振兴

风吹稻田千层浪&#xff0c;眼下&#xff0c;台山四九镇的早稻长势喜人&#xff0c;沉甸甸的稻穗迎风而动&#xff0c;已进入破口抽穗的关键期&#xff0c;即将在6月底陆续迎来丰收。在台山市明华汇种养专业合作社管理的稻田里&#xff0c;合作社负责人梁明喜正仔细观察着稻苗的…

山东大学软件学院项目实训-创新实训-基于大模型的旅游平台(三十二)- 微服务(12)

目录 12.8 RestClient查询文档 12.8.1 快速入门 12.8.2 match&#xff0c; term&#xff0c;bool&#xff0c;range查询 12.8.3 排序和分页 12.8.4 高亮 12.8 RestClient查询文档 12.8.1 快速入门 Testvoid testMatchALL() throws IOException {// 1. 准备requestSearchReq…

以bert为例,了解Lora是如何添加到模型中的

以bert为例,了解Lora是如何添加到模型中的 一.效果图1.torch.fx可视化A.添加前B.添加后 2.onnx可视化A.添加前B.添加后 3.tensorboard可视化A.添加前B.添加后 二.复现步骤1.生成配置文件(num_hidden_layers1)2.运行测试脚本 本文以bert为例,对比了添加Lora模块前后的网络结构图…

Linux下的GPIO编程

目录 一、前言 二、sysfs方式 1、sysfs简介 2、基本目录结构 3、编号计算 4、sysfs方式控制GPIO 三、libgpiod库 1、libgpiod库简介 2、API函数 四、LED灯编程 一、前言 在Linux下&#xff0c;我们通常使用 sysfs 和 libgpiod库 两种方式进行控制GPIO&#xff0c;目前…

DDei在线设计器-属性编辑器

DDei-Core-属性编辑器 DDei-Core-属性编辑器插件包含了文本、大文本、数值、下拉、单选、勾选以及颜色等属性编辑。 图形和属性共同构成一个完整的定义&#xff0c;属性编辑器就是编辑属性值的控件。当选中图形实例时&#xff0c;属性面板就会展现当前实例的所有属性以及属性编…