非参数计量经济学附录5程序选.docx

1、非参数计量经济学附录5程序选附录5 MATLAB和GAUSS程序精选 本附录为本书各章例子计算的MATLAB和GAUSS程序。利用MATLAB和GAUSS软件进行计算必须编程，详细介绍MATLAB和GAUSS软件的编程语言不是本书的目的，读者想完全理解本附录的程序并举一反三，可仔细阅读附录3和附录4，并参考MATLAB和GAUSS软件的帮助文件。精通GAUSS编程的捷径之一是学习林光平的著作计算计量经济学计量经济学家和金融分析师GAUSS编程与应用(清华大学出版社，2003)。1.MATLAB程序。第2章例2.2中随机变量密度函数的核估计及其置信水平为95的置信区间。h为窗宽，f为密度函数的核

2、估计，y=sortrows(x f)是对矩阵x f按照第一列由小到大排序。x=0;-0.04;0.18;-0.91;-0.95;0.16;-0.13;0.08;0.7;3.01;-1.15;0.21;-0.11;-0.31;-0.67;-1.22;-1.49;1.45;-0.88;0.72;0.46;-2.03;-0.41;0.99;-0.4;-0.93;1.1;0.66;0.61;2.32;-0.48;-1.42;-0.89;0.93;0.72;0.6;1.06;-0.19;-0.3;-0.38;-0.86;0.99;-1.04;-0.11;1.05;-1.37;0.28;-0.08;-0.

3、38;0.73;-1.52;-0.3;0.96;2.98;0.55;-0.03;0.96;-0.34;0.12;0.2;0.23;0.37;-0.1;1.35;-0.98;0.29;-1.6;0.76;1.17;0.49;-1.64;0.69;-0.75;1.97;-0.61;0.23;-1.07;0.23;0.95;0.56;0.8;-0.62;-0.19;0.95;-0.6;-2.07;0.47;0.11;-1.3;1.25;-0.77;0.96;-1.18;0.28;0.34;-0.48;-0.38;0.26;-0.79;1.82;n=100;%the number of the obs

4、ervationsh=0.5;f=zeros(n,1);for i=1:n; f(i)=mean(1/sqrt(2*pi)*exp(-0.5*(x-x(i)/h).2)/h;end;y=sortrows(x f);m=y(:,2);b=1/sqrt(4*pi);m0=(1.96/sqrt(n*h)*sqrt(m*b);plot(y(:,1),m,-,y(:,1),m+m0,y(:,1),m-m0);2. MATLAB程序。第2章例2.3中选择3阶核函数的减少核估计偏的方法对随机变量密度函数的核估计及其置信水平为95的置信区间。h为窗宽，f为密度函数的减少偏方法的核估计，y=sortrows(x

5、 f)是对矩阵x f按照第一列由小到大排序。x=0;-0.04;0.18;-0.91;-0.95;0.16;-0.13;0.08;0.7;3.01;-1.15;0.21;-0.11;-0.31;-0.67;-1.22;-1.49;1.45;-0.88;0.72;0.46;-2.03;-0.41;0.99;-0.4;-0.93;1.1;0.66;0.61;2.32;-0.48;-1.42;-0.89;0.93;0.72;0.6;1.06;-0.19;-0.3;-0.38;-0.86;0.99;-1.04;-0.11;1.05;-1.37;0.28;-0.08;-0.38;0.73;-1.52;-

6、0.3;0.96;2.98;0.55;-0.03;0.96;-0.34;0.12;0.2;0.23;0.37;-0.1;1.35;-0.98;0.29;-1.6;0.76;1.17;0.49;-1.64;0.69;-0.75;1.97;-0.61;0.23;-1.07;0.23;0.95;0.56;0.8;-0.62;-0.19;0.95;-0.6;-2.07;0.47;0.11;-1.3;1.25;-0.77;0.96;-1.18;0.28;0.34;-0.48;-0.38;0.26;-0.79;1.82;n=100;%the number of the observationsh=0.5;

7、f=zeros(n,1);for i=1:n; f(i)=mean(0.5/sqrt(2*pi)*(3-(x-x(i)/h).2).*exp(-0.5*(x-x(i)/h).2)/h;end;y=sortrows(x f);m=y(:,2);b=0.476;m0=(1.96/sqrt(n*h)*sqrt(m*b);plot(y(:,1),m,-,y(:,1),m+m0,y(:,1),m-m0);3.MATLAB程序。第2章例2.5的多元密度函数的核估计。x= mvnrnd(mu,sigma,n)为生成均值为mu,协方差阵为sigma的n个二元正态随机向量序列，h为窗宽，f0为实际密度函数值，f

8、为密度函数的核估计。n=120;%the number of the observationsmu = 0 0;sigma = 1 1.5; 1.5 3;x= mvnrnd(mu,sigma,n);h=0.4;f0=zeros(n,1);f=zeros(n,1);p=ones(n,1);for i=1:n; w=zeros(n); for j=1:n; e=(x(j,1)-x(i,1)/h)2+(x(j,2)-x(i,2)/h)2; w(j,j)=(2/pi)*(1-e)*(e1); end; f(i)=(p*w*p)/(n*h2);%It is the estimator of f(x) f

9、0(i)=(1/sqrt(3*pi2)*exp(-2*(x(i,1)2-x(i,1)*x(i,2)+x(i,2)2/3);end;x(:,1) x(:,2) f0 f4.MATLAB程序。计算第3章例3.11的局部线性估计及其置信水平为95%的置信区间，其中核函数选择，窗宽h=0.08。d为1203的数据阵，mp为回归函数的局部线性估计，vp是mp的标准差的估计。n=120;%the number of the observationsd= ; %input data in itx=d(:,3);y=d(:,2);h=0.08;mp=zeros(n,1);p=ones(n,1);for i=1

10、:n; x0=p x-x(i); w=zeros(n); for j=1:n; e=(x(j)-x(i)/h)2; w(j,j)=0.75*(1-e)*(e1); end; a=inv(x0*w*x0)*(x0*w*y); mp(i)=a(1);%It is the estimator of m(x) end;u=y-mp;vp=zeros(n,1);p=ones(n,1);for i=1:n; x0=p x-x(i); w=zeros(n); v=zeros(2); for j=1:n; e=(x(j)-x(i)/h)2; w(j,j)=0.75*(1-e)*(e1); v=v+w(j,j)

11、2*u(j)2*x0(j,:)*x0(j,:); end; o=inv(x0*w*x0); vv=o*v*o;%It is the estimator of m(x) vp(i)=sqrt(vv(1);end;plot(x,mp,:,x,mp-1.96*vp,x,mp+1.96*vp,x,y,.)5.MATLAB程序。计算第3章例3.15的变窗宽局部线性估计，变窗宽为解释变量密度函数的核估计，不变窗宽为0.05, 核函数为。d为1203的数据阵，fp为解释变量密度函数的核估计，mp为回归函数的局部线性估计，hv为变窗宽，vp是mp的标准差的估计。n=120;%the number of the

12、 observationsd= ;%input data in itx=d(:,3);y=d(:,2);h=0.6;fp=zeros(n,1);p=ones(n,1);for i=1:n; w=zeros(n); for j=1:n; e=(x(j)-x(i)/h)2; w(j,j)=0.75*(1-e)*(e1); end; fp(i)=(n*h)(-1)*p*w*p; end;%plot(x,fp,.)mp=zeros(n,1);p=ones(n,1);for i=1:n; hv=0.054/fp(i); x0=p x-x(i); w=zeros(n); for j=1:n; e=(x(j

13、)-x(i)/hv)2; w(j,j)=0.75*(1-e)*(e1); end; a=inv(x0*w*x0)*(x0*w*y); mp(i)=a(1);%It is the estimator of m(x) end;u=y-mp;vp=zeros(n,1);p=ones(n,1);for i=1:n; x0=p x-x(i); w=zeros(n); v=zeros(2); for j=1:n; e=(x(j)-x(i)/h)2; w(j,j)=0.75*(1-e)*(e0);end;q=p ones(n,1) x x.2 x.3;a=inv(q*q);b=a*q*y;f=q*b;se=

14、zeros(m,1);z=zeros(m,1);for i=1:m; se(i)=sqrt(y-f)*(y-f)/n)*sqrt(a(i,i); z(i)=b(i)/se(i);end;t zplot(x,f,.);程序2：固定节点的多项式样条估计。d为1203的数据阵，t为节点，q为样条函数的数据阵，f为回归函数的样条估计，z为检验节点显著与否的统计量。d= ; %input data in itn=120;x=d(:,3);y=d(:,2);t= ; %input ts in its=size(t);m=s(1);p=zeros(n,m);for i=1:m; p(:,i)=(x-t(i)

15、.3.*(x-t(i)0);end;q=p ones(n,1) x x.2 x.3;a=inv(q*q);b=a*q*y;f=q*b;se=zeros(m,1);z=zeros(m,1);for i=1:m; se(i)=sqrt(y-f)*(y-f)/n)*sqrt(a(i,i); z(i)=b(i)/se(i);end;t zplot(x,f,.);10.MATLAB程序。计算第4章例4.1的不变窗宽核估计。d为1213的数据阵，m0为回归函数值，m为回归函数的核估计。n=121;d= ; %input data in ity=d(:,1);x=d(:,2:3);m0=10*(x(:,1)

16、+1).*(x(:,2)+1);h=0.24;m=zeros(n,1);p=ones(n,1);for i=1:n; w=zeros(n); for j=1:n; e=(x(j,1)-x(i,1)/h)2+(x(j,2)-x(i,2)/h)2; w(j,j)=(2/pi)*(1-e)*(e1); end; m(i)=(p*w*y)/(p*w*p);%It is the estimator of m(x)end;y x(:,1) x(:,2) m0 msqrt(mean(m0-m).*(m0-m)11.MATLAB程序。计算第4章例4.1的不变窗宽局部线性估计。d为1213的数据阵，m0为回归函

17、数值，m为回归函数的局部线性估计。n=121;d= ; %input data in ity=d(:,1);x=d(:,2:3);m0=10*(x(:,1)+1).*(x(:,2)+1);h=0.6;m=zeros(n,1);p=ones(n,1);for i=1:n; x0=p x(:,1)-x(i,1) x(:,2)-x(i,2); w=zeros(n); for j=1:n; e=(x(j,1)-x(i,1)/h)2+(x(j,2)-x(i,2)/h)2; w(j,j)=(2/pi)*(1-e)*(e1); end; a=inv(x0*w*x0)*(x0*w*y); m(i)=a(1);

18、%It is the estimator of m(x)end;y x(:,1) x(:,2) m0 msqrt(mean(m0-m).*(m0-m) 12.MATLAB程序。计算第4章例4.1的近邻估计。采用一元核的k近邻估计，选择三角形核函数。d为1213的数据阵，m0为回归函数值，m为回归函数的近邻估计。n=121;d= ; %input data in ity=d(:,1);x=d(:,2:3);m0=10*(x(:,1)+1).*(x(:,2)+1);k=10;m=zeros(n,1);p=ones(k,1);for i=1:n; q=(x(:,1)-x(i,1).2+(x(:,2)

19、-x(i,2).2; d0=q x(:,1) x(:,2) y; d=sortrows(d0); w=zeros(k); for j=1:k; e=d(j,1)/d(k,1); w(j,j)=(1-e)*(en;p=(xi,1*ones(n,1)(xi,2*ones(n,1)(xi,3*ones(n,1); x0=ones(n,1)(x-p);w=zeros(n,n); j = 1; do until jn; h=-500+40000*j/n; e=(xj,1-xi,1)/h)2+(xj,2-xi,2)/h)2+(xj,3-xi,3)/h)2; wj,j=1*(1-e)*(e1); j=j+1; endo;q=x0*w*x0;b=inv(q);f=x0*w*y;ai,.=(b*f);ui=yi-ai,1;i=i+1; endo;emc=(y-a.,1)2;print sqrt(meanc(emc);print a;end;16.GAUSS程序。计算第