原代码:
n=1000000;
ml=100:50:2000;
alpha=0.05;
srl=0.4:0.01:0.6;
t(length(srl),n,length(ml))=0;
q(ml(length(ml)))=0;
h=1;
for j=ml
for k=1:n
sq=[];
u=rand(j,1);
for sr=srl
for l=1:j
if u(l)<=sr
q(l)=0;
else
q(l)=1;
end
end
q=q(:);
sq=cat(2,sq,q);
sumsq=sum(sq);
end
for i=1:length(srl)
t(i,k,h)=(sumsq(i)/j-0.5)/((sumsq(i)/j)*(1-sumsq(i)/j)/j)^0.5;
end
end
h=h+1;
end
p=normpdf(t);
v(length(srl),length(ml))=0;
pi(length(srl),length(ml))=0;
for m=1:length(ml)
for r=1:length(srl)
for s=1:n
if p(r,s,m)<=alpha
v(r,m)=v(r,m)+1;
end
end
pi(r,m)=1-v(r,m)/n;
end
end
px(length(srl),length(ml))=0;
for w=1;length(ml)
px(:,:,w)=sort(p(:,:,w),2);
end
pxtq(length(srl),length(ml))=0;
for y=1:length(ml)
for z=1:length(srl)
pxtq(z,y)=px(z,0.95*n,y);
end
end
a=0.4:0.01:0.6;
b=100:20:2000;
[x,y]=meshgrid(a,b);
z=10;
mesh(x,y,z);
xlabel('pi');
ylabel('n');
zlabel('alpha');
能否把n改成500000或者100000