matlab 解决数据预测作业

1、1某一个农场自1990-2000年的各年份的粮食的总产量分别为8.2,8.7,9.4,9.1,9.8,9.8,10.5,10.6,10.8,10.6,12.4.试用加权平滑法，预测该农场2005年的产量。年份19901991199219931994199519961997199819992000产量8.28.79.49.19.89.810.510.610.810.612.4首先我们画出折线图和散点图，可以比较直观的看出年份和产量之间的一些走势。如图：>> x=1990;1991;1992;1993;1994;1995;1996;1997;1998;1999;2000;y=8.2;8

2、.7;9.4;9.1;9.8;9.8;10.5;10.6;10.8;10.6;12.4;>> plot(x,y) 折线图 散点图从上面的图形可以看出我们得到的年份和产量间的关系符合直线关系，故我们采用3次加权平滑法来处理相关数据。从图知道近期数据的变化趋势应在预测模型中得到从分的反应，加权系数应该取得大一点，这里我们取0.3.做3次加权平滑程序（加权系数去0.3）：close allclear,clc% 统计数据 实际值arr = 8.2 8.7 9.4 9.1 9.8 9.8 10.5 10.6 10.8 10.6 12.4'm,nn=size(arr); alpha =

3、 0.3; % 平滑常数的范围为0,1% 1次指数平滑s1 = zeros(m,1);s1(1,1) = arr(1,1);for i=2:m s1(i) = alpha*arr(i,1)+(1-alpha)*s1(i-1);endsx1 = s1% 2次指数平滑s2 = zeros(m,1);s2(1,1) = arr(1,1);for i=2:m s2(i) = alpha*s1(i,1)+(1-alpha)*s2(i-1);endsx2 = s2% 3次指数平滑s3 = zeros(m,1);s3(1,1) = arr(1,1);for i=2:m s3(i) = alpha*s2(i,

4、1)+(1-alpha)*s3(i-1);endsx3 = s3% 计算二次曲线中的参数a = zeros(m,1);b = zeros(m,1);c = zeros(m,1);beta=alpha/(2*(1-alpha)*(1-alpha);a = 3*sx1-3*sx2+sx3;b = beta*(6-5*alpha)*sx1-2*(5-4*alpha)*sx2+(4-3*alpha)*sx3);c = beta*alpha*(sx1-2*sx2+sx3);% 二次曲线模型 a+b*t+c*t*tt = :11; % sf = zeros(m,1);% sf(1,1) = arr(1,1

5、);sf = a+b*t+c*t*t % 预测得到的结果如下：sx1 = 8.2000 8.3500 8.6650 8.7955 9.0968 9.3078 9.6655 9.9458 10.2021 10.3215 10.9450sx2 = 8.2000 8.2450 8.3710 8.4983 8.6779 8.8669 9.1064 9.3583 9.6114 9.8244 10.1606sx3 = 8.2000 8.2135 8.2607 8.3320 8.4358 8.5651 8.7275 8.9167 9.1251 9.3349 9.5826sf = 8.2000 8.6500

6、9.4600 9.4870 10.0563 10.2217 10.8316 11.0832 11.2590 11.0469 12.4869带入相关数据到公式中得：a=0.019b=0.160c=11.94y=0.019T2+0.16T+11.94把T=5带入公式得y=13.215，即是2005年的产量预计在13.215左右。注：相关的数据综合见手写部分！2.已知大豆亩产量与鼓粒期的土壤湿度的关系如图表所示，试求他们的相关变化规律，并预测湿度为22%时的大豆产量。湿度%26.929.019.524.224.323.024.819.424.523.924.117.6产量14713871110989

7、199831009810659方法1 我们严格的按照一元线性回归的方法来做。一元线性回归的所需的各数据程序和结果如下：x=(26.9+29+19.5+24.2+24.3+23.0+24.8+19.4+24.5+23.9+24.1+17.6)/12x = 23.4333 （x的平均值） y=(147+138+71+110+98+91+99+83+100+98+106+59)/12y = 100 （y的平均值）x=(26.92+292+19.52+24.22+24.32+23.02+24.82+19.42+24.52+23.92+24.12+17.62)/12x = 558.6183 （x2的平均

8、值）y=(1472+1382+712+1102+982+912+992+832+1002+982+1062+592)/12y = 1.0558e+004 （y2的平均值）xy=(26.9*147+29*138+19.5*71+24.2*110+24.3*98+23.0*91+24.8*99+19.4*83+24.5*100+23.9*98+24.1*106+17.6*59)/12xy = 2.4107e+003 （x*y的平均值）b=(2410.7-23.4333*100)/(558.6183-23.43332)b = 7.0925a=100-23.4333*7.0925a = -66.200

9、7得出关系式 y=7.0925*x-66.2007 把X=22带入方程得X=89.8N=12 查表得R0.05=0.576，R0.01=0.708由于R>R0.01，故线性关系极为明显。故在湿度为22%时的产量为89.8，查t分布表得t(a/2,n-2)=t(0.025,10)=2.23由于N<30故属于小样本。残差：y1=147-(7.0925*26.9-66.2007)y2=138-(7.0925*29-66.2007)y3=71-(7.0925*19.5-66.2007)y4=110-(7.0925*24.2-66.2007)y5=98-(7.0925*24.3-66.200

10、7)y6=91-(7.0925*23-66.2007)y7=99-(7.0925*24.8-66.2007)y8=83-(7.0925*19.4-66.2007)y9=100-(7.0925*24.5-66.2007)y10=98-(7.0925*23.9-66.2007)y11=106-(7.0925*24.1-66.2007)y12=59-(7.0925*17.6-66.2007)y1 = 22.4124y2 =-1.4818y3 = -1.1031y4 = 4.5622y5 = -8.1471y6 = -5.9268y7 = -10.6933y8 = 11.6062y9 = -7.565

11、6y10 = -5.3101y11 = 1.2714y12 = 0.3727s=(22.41242+1.48182+1.10312+4.56222+8.14712+5.92682+10.69332+11.60622+7.56562+5.31012+1.27142+0.37272)/10)0.5s = 9.8198c=(1+1/12+(22-23.433)/(26.9-22)2+(29-22)2+(19.5-22)2+(24.2-22)2+(24.3-22)2+(23-22)2+(24.8-22)2+(19.4-22)2+(24.5-22)2+(23.9-22)2+(24.1-22)2+(17.

12、6-22)2)0.5c = 1.0359r=(2410.7-2343.33)/(558.6183-23.43332)*(10558-10000)0.5r = 0.9254所以把数值带引入的得22%湿度时产量为Y=89.8±22.6843即当湿度为22%时的产量由95%的可能落在Y=89.8±22.6843区间。方法2.我们先画出湿度和大豆产量的的散点图，来观察他们之间的关系走势，可以大致判断他们有一定的线性关系存在，但是我们先用曲线的方式来进行拟合，看能否得出一条比较全面的反应各个数据的关系的曲线图。我们做3、4次高阶函数的图像来进行比较看看。： 散点图首先是3次程序：x=

13、26.9 29.0 19.5 24.2 24.3 23.0 24.8 19.4 24.5 23.9 24.1 17.6;y=147 138 71 110 98 91 99 83 100 98 106 59;plot(x,y,'ro');p=polyfit(x,y,3);%于是拟合的曲线就是入下面的公式所示的那样p(1)x4+p(2)x3+p(3)x2+p(4)x+p(5)，想拟合成其它次数的多项式只需将3改为相应的次数即可f=poly2sym(p);xinterp=16 18 20 22 24 26 28 30;yinterp=subs(f,xinterp);hold on;p

14、lot(xinterp,yinterp,'o');ezplot(f,15,30)3次项拟合曲线4次：x=26.9 29.0 19.5 24.2 24.3 23.0 24.8 19.4 24.5 23.9 24.1 17.6;y=147 138 71 110 98 91 99 83 100 98 106 59;plot(x,y,'ro');p=polyfit(x,y,4);%于是拟合出的曲线就是p(1)x4+p(2)x3+p(3)x2+p(4)x+p(5)，想拟合成其它次数的多项式只需将4改为相应的次数即可f=poly2sym(p);xinterp=16 18 2

15、0 22 24 26 28 30;yinterp=subs(f,xinterp);hold on;plot(xinterp,yinterp,'o');ezplot(f,15,30) 4次项拟合曲线可见，4此项拟合已经比较符合实际图形的真是面目，故采用4次项拟合曲线结果。 拟合出的相关系数为：a = Columns 1 through 3 -6.505434592838499e-002 5.971906071853749e+000 -2.031637491078744e+002 Columns 4 through 5 3.042760977261934e+003 -1.68775

16、1264312053e+004故方程为： Y=-0.0651x4+5.9719x3-203.1637x2+3042.7609x-16877.5126把X=22带入得到： Y=70.722（数据误差较大，慎用！误差大的原因是系数的精确位数（因为一个数的4次方后结果变化很大，对应的系数对结果的影响就比较大）方法3由于他们存在一定的曲线关系（其实线性不算太明显），我们用y=a*xb这样的幂函数来进行拟合，得到下面的图形： （ 基本思路是将方程两边取对数，构成性关系（其实就是一元线性回归思路，这里不再赘述）。相关结果和程序如下所示：程序：x=26.9 29.0 19.5 24.2 24.3 23.0

17、24.8 19.4 24.5 23.9 24.1 17.6;y=147 138 71 110 98 91 99 83 100 98 106 59;beta0=0.1,3'beta,r,J=nlinfit(x',y','yut',beta0);betaYY,delta=nlpredci('yut',x',beta,r,J);plot(x,y,'k+',x,YY,'r');beta = 0.4219 1.7301即 y=0.4219x1.7301把x=22带入公式得y=88.6614.（置信区间等验算过

18、程省略，与方法1相似）注：相关的数据综合见手写部分！3：某县1995-2000年的工业产值统计表如下：年份199519961997199819992000年产值140814211544178721602678根据此6年的工业产值进行线性回归分析，用来预测2001-2015年的工业产值。理论上讲这是不够的，但得出的条曲线比较符合实际情况，所以加以选用，考虑到2015年的工业产值要达到3.0亿元，才用回归分析，选用三次曲线比较合乎实际情况，即： Y=a+bx+cx2+dx3用最小二乘法求出各数，并对给工业发展进行预测。 首先我们来求各个系数。利用最小2乘法求系数的方法和程序如下：在MATLAB工作

19、窗口输入程序x=1 2 3 4 5 6;y=1408 1421 1544 1787 2160 2675;plot(x,y,'r*'),legend('实验数据(xi,yi)')xlabel('x'), ylabel('y'),title('数据点(xi,yi)的散点图')运行后屏幕显示数据的散点图编写下列MATLAB程序计算在处的函数值，即输入程序syms a1 a2 a3 a4x=1 2 3 4 5 6; fi=a1.*x.3+ a2.*x.2+ a3.*x+ a4运行后屏幕显示关于a1,a2, a3和a4的线性

20、方程组fi=a1+a2+a3+a4, 8*a1+4*a2+2*a3+a4, 27*a1+9*a2+3*a3+a4, 64*a1+16*a2+4*a3+a4, 125*a1+25*a2+5*a3+a4, 216*a1+36*a2+6*a3+a4写构造误差平方和的MATLAB程序y=1408 1421 1544 1787 2160 2675;fi=a1+a2+a3+a4,8*a1+4*a2+2*a3+a4,27*a1+9*a2+3*a3+a4,64*a1+16*a2+4*a3+a4,125*a1+25*a2+5*a3+a4,216*a1+36*a2+6*a3+a4fy=fi-y; fy2=fy.2

21、; J=sum(fy.2)fi =a1+a2+a3+a4,8*a1+4*a2+2*a3+a4,27*a1+9*a2+3*a3+a4,64*a1+16*a2+4*a3+a4, 125*a1+25*a2+5*a3+a4, 216*a1+36*a2+6*a3+a4运行后屏幕显示误差平方和如下J=(a1+a2+a3+a4-1408)2+(8*a1+4*a2+2*a3+a4-1421)2+(27*a1+9*a2+3*a3+a4-1544)2+(64*a1+16*a2+4*a3+a4-1787)2+(125*a1+25*a2+5*a3+a4-2160)2+(216*a1+36*a2+6*a3+a4-267

22、5)2为求使达到最小，只需利用极值的必要条件 ，程序如下：syms a1 a2 a3 a4;J=(a1+a2+a3+a4-1408)2+(8*a1+4*a2+2*a3+a4-1421)2+(27*a1+9*a2+3*a3+a4-1544)2+(64*a1+16*a2+4*a3+a4-1787)2+(125*a1+25*a2+5*a3+a4-2160)2+(216*a1+36*a2+6*a3+a4-2675)2Ja1=diff(J,a1); Ja2=diff(J,a2); Ja3=diff(J,a3); Ja4=diff(J,a4);Ja11=simple(Ja1), Ja21=simple(J

23、a2), Ja31=simple(Ja3), Ja41=simple(Ja4), syms a1 a2 a3 a4;J=(a1+a2+a3+a4-1408)2+(8*a1+4*a2+2*a3+a4-1421)2+(27*a1+9*a2+3*a3+a4-1544)2+(64*a1+16*a2+4*a3+a4-1787)2+(125*a1+25*a2+5*a3+a4-2160)2+(216*a1+36*a2+6*a3+a4-2675)2Ja1=diff(J,a1); Ja2=diff(J,a2); Ja3=diff(J,a3); Ja4=diff(J,a4);Ja11=simple(Ja1), J

24、a21=simple(Ja2), Ja31=simple(Ja3), Ja41=simple(Ja4), J =(a1+a2+a3+a4-1408)2+(8*a1+4*a2+2*a3+a4-1421)2+(27*a1+9*a2+3*a3+a4-1544)2+(64*a1+16*a2+4*a3+a4-1787)2+(125*a1+25*a2+5*a3+a4-2160)2+(216*a1+36*a2+6*a3+a4-2675)2运行后屏幕显示J分别对a1, a2 ,a3 ,a4的偏导数如下：Ja11 =-2033264+134342*a1+24402*a2+4550*a3+882*a4Ja21 =

25、-399760+24402*a1+4550*a2+882*a3+182*a4Ja31 =-85760+4550*a1+882*a2+182*a3+42*a4Ja41 =-21990+882*a1+182*a2+42*a3+12*a4解线性方程组Ja11 =0，Ja21 =0，Ja31 =0，Ja41 =0，输入下列程序：A=134342, 24402, 4550, 882; 24402, 4550, 882, 182; 4550, 882, 182,42; 882, 182, 42, 12;B=2033264, 399760, 85760, 21990; format long eC=B/A, f=poly2sym(C)运行后屏幕显示拟合函数f及其系数C如下C = Columns 1