几种常见的插值、多元线性回归

 一、利用matlab命令进行以下数据的不同插值（如最近邻插值、分段线性插值、分段三次样条插值、分段三次多项式插值）并作图比较。离散数据如下：

x	0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1
y	0.3	0.5	1	1.4	1.6	1.9	0.6	0.4	0.8	1.5	2

 

x = [0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1];
y = [0.3 05 1 1.4 1.6 1.9 0.6 0.4 0.8 1.5 2];

% 最近邻插值
xq_nn = 0:0.01:1;  % 插值间隔
yq_nn = interp1(x, y, xq_nn, 'nearest');

% 分段线性插值
xq_linear = 0:0.01:1;  % 插值间隔
yq_linear = interp1(x, y, xq_linear, 'linear');

% 分段三次样条插值
xq_spline = 0:0.01:1;  % 插值间隔
yq_spline = interp1(x, y, xq_spline, 'spline');

% 分段三次多项式插值
xq_pchip = 0:0.01:1;  % 插值间隔
yq_pchip = interp1(x, y, xq_pchip, 'pchip');

% 绘制原始数据
figure;
plot(x, y, 'ro', 'DisplayName', '原始数据');
hold on;

% 绘制最近邻插值曲线
plot(xq_nn, yq_nn, 'g-', 'DisplayName', '最近邻插值');

% 绘制分段线性插值曲线
plot(xq_linear, yq_linear, 'b-', 'DisplayName', '分段线性插值');

% 绘制分段三次样条插值曲线
plot(xq_spline, yq_spline, 'm-', 'DisplayName', '分段三次样条插值');

% 绘制分段三次多项式插值曲线
plot(xq_pchip, yq_pchip, 'c-', 'DisplayName', '分段三次多项式插值');

% 添加图例和标签
legend('Location','best');
xlabel('x');
ylabel('y');
title('不同插值方法的比较');

% 显示图形
hold off;

二、某人做一种材料的伸缩实验,t为温度(℃),L为长度(mm),实验数据见下表

t	L
20	81
25	82.3
30	84
35	86.8
40	89

用二次多项式拟合求L与t的表达式.要求:1、写出L与t的关系式；2、估计t=38时L的值.

 

%原始数据
t = [20 25 30 35 40];
L = [81 82.3 84 86.8 89];

%进行二次多项式拟合
coeff = polyfit(t, L, 2);

%关系式系数
L_equation = sprintf('L = %.4f*t^2 + %.4f*t + %.4f', coeff(1), coeff(2), coeff(3));
disp(L_equation);

%估计t=48时的值
L_estimate = polyval(coeff, 38);
disp(L_estimate);

输出结果：

三、财政收入预测问题：财政收入与国民收入、工业总产值、农业总产值、总人口、就业人口、固定资产投资等因素有关。下表列出了1952-1981年的原始数据，试构造多元线性回归模型。

年份	国民收入（亿元）	工业总产值(亿元)	农业总产值（亿元）	总人口（万人）	就业人口（万人）	固定资产投资（亿元）	财政收入(亿元)
1952	598	349	461	57482	20729	44	184
1953	586	455	475	58796	21364	89	216
1954	707	520	491	60266	21832	97	248
1955	737	558	529	61465	22328	98	254
1956	825	715	556	62828	23018	150	268
1957	837	798	575	64653	23711	139	286
1958	1028	1235	598	65994	26600	256	357
1959	1114	1681	509	67207	26173	338	444
1960	1079	1870	444	66207	25880	380	506
1961	757	1156	434	65859	25590	138	271
1962	677	964	461	67295	25110	66	230
1963	779	1046	514	69172	26640	85	266
1964	943	1250	584	70499	27736	129	323
1965	1152	1581	632	72538	28670	175	393
1966	1322	1911	687	74542	29805	212	466
1967	1249	1647	697	76368	30814	156	352
1968	1187	1565	680	78534	31915	127	303
1969	1372	2101	688	80671	33225	207	447
1970	1638	2747	767	82992	34432	312	564
1971	1780	3156	790	85229	35620	355	638
1972	1833	3365	789	87177	35854	354	658
1973	1978	3684	855	89211	36652	374	691
1974	1993	3696	891	90859	37369	393	655
1975	2121	4254	932	92421	38168	462	692
1976	2052	4309	955	93717	38834	443	657
1977	2189	4925	971	94974	39377	454	723
1978	2475	5590	1058	96259	39856	550	922
1979	2702	6065	1150	97542	40581	564	890
1980	2791	6592	1194	98705	41896	568	826
1981	2927	6862	1273	100072	73280	496	810

 （）

%原始输入test;
x_test = test(2:7,:);
y_test = test(8,:); %原始输出
y_use = y_test';
temp = ones(30,1);
x_use = x_test';
x_end = [temp,x_use];
%多元线性回归
[b, bint,r,rint,stats]=regress(y_use,x_end);
rcoplot(r,rint)    %第一个图像

y_temp = x_end*b;
year = test(1,:);
y_out = y_temp';

% 添加图例和标签

plot(year,y_test,'ob-',year,y_temp,'*r-' )
legend('原始数据','回归预测');
xlabel('year');
ylabel('y');
title('多元线性回归模型');
hold on

（创建变量test）

输出结果：

 四、题目

 

clear
clc

%初始数据
a=100;
m=1400;
g=9.78;
w=m*g/a;
N=100; 
tol=10^(-6);

%定义函数
f=inline('a*(1+8/3*(x/a)^2-32/5*(x/a)^4+256/7*(x/a)^6)-L','a','x','L');
df=inline('(8/3*2*(x/a)-32/5*4*(x/a)^3+256/7*6*(x/a)^5)','a','x');

i=1;
x0=a/2;
for j=1:10
    L(j)=a+5*j;
    while(i<=100)
        x=iteration(f,df,a,x0,L(j));
        if abs(x-x0)<=tol
            break;
        end;
        x0=x;
        i=i+1;
    end
    T(j)=w*a/2*sqrt(1+(a/x)^2); 
end
plot(L,T)

function [x]=iteration(f,df,a,x0,L)
x=x0-f(a,x0,L)/df(a,x0);
end

 

 

 

 利用dsolve求得真值，并绘制图像

fplot('exp(2*t)',[0,1])

 

 

dt=0.05;
y(1)=1;
t=0:dt:1;
     
for i=1:20  
    y(i+1)=y(i)+2*y(i)*dt;
end

plot(t,y,'m*')

欧拉法 

 

 优化后的欧拉法

dt=0.05;y(1)=1;
 t=0:dt:1;

 for i=1:20  
     y1(i+1)=y(i)+2*y(i)*dt;
     y(i+1)=y(i)+1/2*(2*y1(i+1)+2*y(i))*dt;
 end

 plot(t,y,'b+')