Q=[4889.5 5330.5 5985.6 7243.8 9040.7 10274.4 12050.6 15036.8 17000.9 18718.3 21826.2 26937.3 35260.0...
48108.5 59810.5 70142.5 78060.8 83024.3 88479.2 98000.5 108068.2 119095.7 135174.0];
IQ=Q/18718.3;
K=[961 1230.4 1430.1 1832.9 2543.2 3120.6 3791.7 4653.8 4410.4 4517.5 5594.5 8080.1 13072.3 17827.1 20524.9...
23358.6 25259.7 28716.9 29754.6 33110.4 37987.0 45046.9 58616.3];
IK=K/4517.5;
L=[4.5126 4.6358 4.7286 4.8179 4.9873 5.1282 5.2783 5.4334 5.5329 6.4749 6.5491 6.6152 6.6808 6.7455...
6.8065 6.8951 6.9821 7.0637 7.1394 7.2085 7.3025 7.3741 7.4432];
IL=L/6.4749;
Et=zeros(23);
Et=Et(1,1:23);
for t=1:1:23;
Et(t)=log(IL(t)/IK(t));
end
Et;
Wt=Et;
for t=1:1:23;
Wt(t)=log(IQ(t)/IK(t));
end
Wt;
x=Et;
y=Wt;
plot(x,y,'*');
xlabel('E');
ylabel('W');
a=[0.2612 0.2848 0.3198 0.3870 0.4830 0.5489 0.6438 0.8033...
0.9083 1.0000 1.1660 1.4391 1.8837 2.5701 3.1953 3.7473...
4.1703 4.4355 4.7269 5.2355 5.7734 6.3625 7.2215];
y=[0.2127 0.2724 0.3166 0.4057 0.5630 0.6908 0.8393 1.0302...
0.9763 1.0000 1.2384 1.7886 2.8937 3.9462 4.5434 5.1707...
5.5915 6.3568 6.5865 7.3294 8.4089 9.9716 12.9754...
0.6969 0.7160 0.7303 0.7441 0.7703 0.7920 0.8152 0.8391...
0.8545 1.0000 1.0115 1.0217 1.0318 1.0418 1.0512 1.0649...
1.0783 1.0909 1.1026 1.1133 1.1278 1.1389 1.1495];
curvefun=inline('x(1)*(y(1,:).^x(2)).*(y(2,:).^x(3))','x','y')
x0=[0.1,0.1,0.2];
x=lsqcurvefit(curvefun,x0,y,a)
a=x(1),alpha=x(2),beta=x(3)
Q=[4889.5 5330.5 5985.6 7243.8 9040.7 10274.4 12050.6 15036.8 17000.9 18718.3 21826.2 26937.3 35260.0...
48108.5 59810.5 70142.5 78060.8 83024.3 88479.2 98000.5 108068.2 119095.7 135174.0];
IQ=Q/18718.3;
K=[961 1230.4 1430.1 1832.9 2543.2 3120.6 3791.7 4653.8 4410.4 4517.5 5594.5 8080.1 13072.3 17827.1 20524.9...
23358.6 25259.7 28716.9 29754.6 33110.4 37987.0 45046.9 58616.3];
IK=K/4517.5;
L=[4.5126 4.6358 4.7286 4.8179 4.9873 5.1282 5.2783 5.4334 5.5329 6.4749 6.5491 6.6152 6.6808 6.7455...
6.8065 6.8951 6.9821 7.0637 7.1394 7.2085 7.3025 7.3741 7.4432];
IL=L/6.4749;
Qt=zeros(23);
Qt=Qt(1,:);
for t=1:1:23;
Qt(t)=0.9889*IL(t)^0.2167*IK(t)^0.7738;
end
Qt
y1=[0.2612 0.2848 0.3198 0.3870 0.4830 0.5489 0.6438...
0.8033 0.9083 1.0000 1.1660 1.4391 1.8837 2.5701...
3.1953 3.7473 4.1703 4.4355 4.7269 5.2355 5.7734...
6.3625 7.2215];
y2=[0.2761 0.3362 0.3793 0.4615 0.5991 0.7061 0.8262...
0.9742 0.9382 0.9889 1.1697 1.5580 2.2655 2 8862...
3.2250 3.5745 3.8079 4.2159 4.3433 4.7276 5.2727...
6.0289 7.4063];
x=1981:1:2003;
p1=polyfit(x,y1,2);
p2=polyfit(x,y2,2);
xi=1981:0.01:2003;
y3=polyval(p1,xi);
y4=polyval(p2,xi);
plot(x,y1,'*r',xi,y4,'-b')
legend('实际值','预测曲线')
xlabel('年份');
title('预测值与实际值比较图');