score检验统计

				clear;
				clc;
				tic;
				randn('state',0);
				
				for n=25:25:100;
				    a=1
				   for gamma=-0.9:0.1:0.9;
				% delta=1.0e-003;
				
				      for j=1:1000
				          [x,y,m,epsilon,u,G]=data(n,gamma);
				   %求score检验统计量
				          C=eye(n)-(1/n)*ones(n);
				          G1=C*G';
				          H=inv(G1'*G1);
				          SC(j)=u*G1*H*G1'*u'/2;
				        if (abs(SC(j))>3.841)
				            sc(j)=1;
				        else
				           sc(j)=0;
				        end
				      
				      end
				    
				 psc(n/25,a)=sum(sc)/1000; 
				 a=a+1;
				
				%   psc=sum(sc)/1000;
				    end
				end