对于大型数据集，以下代码非常慢...任何人都可以告诉我如何更有效地编写以下内容？逻辑是：首先，符号化变量和参数。其次，使用符号化的变量和参数来编写模型。然后，使用雅可比矩阵来获得每个参数的一阶偏导数。最后，将真实数据代入雅可比矩阵，并将符号结果转换为数值结果。
非常感谢
% Standard errors and significant tests
% I am using firm-level panel data. I have reduced the sample size into
% 3500 firms with 18000 firm-year observations. But it is just too slow.
% Firstly, I define all variables and parameters into the symbolic style.
% sparms1 sparms2 ... sparms14 are parameters. Others are given variables
% contained in vectors.
syms sparms1 sparms2 sparms3 sparms4 sparms5 sparms6 sparms7 sparms8 sparms9 sparms10 sparms11 sparms12 sparms13 sparms14 sparms15 sparms16...
sRstt sRlstt sRvc sRlvc sRik sRlcfk sRlik sRik2 sRmdrate sRsrgrowth sR2007 sR2003 sR2004 sR2005 sR2006...
sRcfk sRl2ik sRlsrgrowth sRl2007 sRl2003 sRl2004 sRl2005 sRl2006...
sRik3 sRlik3 sRlik2 sRl2ik2 sRtltd sRltltd sRhind sRlhind sRpk100 sRfk100 sRcollk100 sRlpk100 sRlfk100 sRlcollk100
% Then, I write the model into the symbolic style. This is the
% first-differenced residual at the firm-year level.
obsres=(((sparms1+sparms2.*sRtltd+sparms3.*sRsrgrowth+sparms4.*sRcfk+sparms13.*sRhind+sparms14.*sRpk100+sparms15.*sRfk100+sparms16.*sRcollk100)...
.*(sRstt-sparms5.*sRvc+sparms6.*sRik2./2+sparms7.*sRik3./(2/3)+(1-sRmdrate).*(sparms6.*sRik+sparms7.*sRik2+1))...
-1-sparms6.*sRlik-sparms7.*sRlik2-sparms8.*sR2004-sparms9.*sR2005...
-sparms10.*sR2006-sparms11.*sR2007-sparms12.*sR2003)...
-((sparms1+sparms2.*sRltltd+sparms3.*sRlsrgrowth+sparms4.*sRlcfk+sparms13.*sRlhind+sparms14.*sRlpk100+sparms15.*sRlfk100+sparms16.*sRlcollk100)...
.*(sRlstt-sparms5.*sRlvc+sparms6.*sRlik2./2+sparms7.*sRlik3./(2/3)+(1-sRmdrate).*(sparms6.*sRlik+sparms7.*sRlik2+1))...
-1-sparms6.*sRl2ik-sparms7.*sRl2ik2-sparms8.*sRl2004-sparms9.*sRl2005...
-sparms10.*sRl2006-sparms11.*sRl2007-sparms12.*sRl2003));
% I use 'jacobian' to calculate the first partial derivative with respect
% to each parameter.
gradf=jacobian(obsres,[sparms1,sparms2,sparms3,sparms4,sparms5,sparms6,sparms7,sparms8,sparms9,sparms10,sparms11,sparms12,sparms13,sparms14,sparms15,sparms16]);
% I substitue all the firm-year level observations from dataset into the jacobian matrix.
gg=subs(gradf,{sparms1 sparms2 sparms3 sparms4 sparms5 sparms6 sparms7 sparms8 sparms9 sparms10 sparms11 sparms12 sparms13 sparms14 sparms15 sparms16...
sRstt sRlstt sRvc sRlvc sRik sRlcfk sRlik sRik2 sRmdrate sRsrgrowth sR2007 sR2003 sR2004 sR2005 sR2006...
sRcfk sRl2ik sRlsrgrowth sRl2007 sRl2003 sRl2004 sRl2005 sRl2006...
sRik3 sRlik3 sRlik2 sRl2ik2 sRtltd sRltltd sRhind sRlhind sRpk100 sRfk100 sRcollk100 sRlpk100 sRlfk100 sRlcollk100},...
{parms(1),parms(2),parms(3),parms(4),parms(5),parms(6),parms(7),parms(8),parms(9),parms(10),parms(11),parms(12),parms(13),parms(14),parms(15),parms(16),...
Rstt,Rlstt,Rvc,Rlvc,Rik,Rlcfk,Rlik,Rik2,Rmdrate,Rsrgrowth,R2007,R2003,R2004,R2005,R2006,...
Rcfk,Rl2ik,Rlsrgrowth,Rl2007,Rl2003,Rl2004,Rl2005,Rl2006,...
Rik3,Rlik3,Rlik2,Rl2ik2,Rtltd,Rltltd,Rhind,Rlhind,Rpk100,Rfk100,Rcollk100,Rlpk100,Rlfk100,Rlcollk100});
% then, I transfer the matrix 'gg' from symbolic style into the numerical style.
ggg=double(gg);
% finstru is the instrumental matrix. 'g' is the average moment condition.
% 'p' is the number of firms.
g=ggg'*finstru./p;
% this is the variance-covariance matrix of parameters
% Here, I use the weighting matrix calcuated from the parameters obtained
% from the first-step optimization. I am doing two-step GMM estimator.
% is it correct?
var=pinv(g*Wsum*g')./p;
SE=sqrt(diag(var)); % standard errors
results=[parms,SE];