我参考了网上曹老师的程序后，自己在细节的地方修改了一下，但是不知道哪里出现了错误，一直提示错误，希望懂的同学和老师可以帮帮忙，先谢谢大家了！
function [paraest,t_sta,V,it,Chi_sta,Pvalue]=gmmestimation

(moment,para0,Y,X,Z,number,K)

nlag=round(size(Y,1)^(1/3));

W(:,:,1)=eye(K);

[para(:,1),fv(:,1)]=fminsearch(moment,para(:,1),[],1,Y,X,Z,W(:,:,1));

% update weight matrix and find optimal solution

for i=2:number

mom=feval_r(moment,para(:,i-1),2,Y,X,Z,W(:,:,i-1));

W(:,:,i)=gmmweightmatrix(mom,nlag);

[para(:,i),fv(:,i)]=fminsearch(moment,para(:,i-1),[],1,Y,X,Z,W(:,:,i));

if abs(fv(:,i)-fv(:,i-1))/abs(fv(:,i-1))<1e-4|fv(i)<=1e-10

  break

end

end

it=i;

if it==number

    error('number of iteration exceeds defined maximum number')

else

    paraest=para(:,it);% optimal parameter

    f0=feval_r(moment,paraest,3,Y,X,Z,W(:,:,it));% optimal function value

% find covariance matrix of estimated parameter vector, using

numerical method

    for j=1:length(para0)

        a=zeros(length(para0),1);

        eps=max(paraest(j)*1e-6,1e-5);

        a(j)=eps;

        M(:,j)=(feval_r(moment,paraest+a,3,Y,X,Z,W(:,:,it))-f0)/eps;

    end

end

V=pinv(M'*W(:,:,it)*M)/size(Y,1);

stderror=sqrt(diag(V));

t_sta=paraest./(stderror);

Chi_sta=size(Y,1)*fv(it);

Pvalue=1-chi2cdf(Chi_sta,K-length(para0));


%find the weight matrix using Newey and West method

function W=gmmweightmatrix(mom,nlag)

q=size(mom,2);T=size(mom,1);

a2=zeros(q,q);a3=zeros(q,q);

for j=1:nlag

    a1=zeros(q,q);

    for i=1:(T-j)

       a1=mom(i+j,:)'*mom(i,:)+a1;

    end

    S(:,:,j)=1/T*a1;

    a2=(1-j/(nlag+1))*S(:,:,j)+a2;

    a3=(1-j/(nlag+1))*S(:,:,j)'+a3;

end

b1=zeros(q,q);

for i=1:T

    b1=mom(i,:)'*mom(i,:)+b1;

end

if nlag==0

    newS=b1*1/T;

else

    newS=a2+a3+b1*1/T;

end

W=pinv(newS);


%%%%%%%%%%%%%%%%%%%%%
function f=cklstest01(para,num,Y,X,Z,W)

[T,q]=size(Y);

alpha=para(1);beta=para(2);sigsq=para(3);

eta=[Y-(alpha+beta*X)/12,(Y-(alpha+beta*X)/12).^2-sigsq.^2/12];

for i=1:T

    m_t(i,:)=kron(eta(i,:),Z(i,:));

end

m=mean(m_t)';

obj=m'*W*m;

if num==1

     f=obj;

elseif num==2
     f=m_t;

elseif num==3

     f=m;

end
%%%%%%%%%%%%%
clear
data=xlsread('gmm.xlsx');
Y=data(:,1);
X=data(:,2);
para0=[0.05;-0.6;0.5];number=20;moment='cklstest01';K=3;
Z=[ones(length(Y),1),X];
[beta,stderror,V,it,Chi_sta,Pvalue]=gmmestimation

(moment,para0,Y,X,Z,number,K);