possion.m

本程序是有关泊松过程的模拟、检验及参数估计

function possion(n,r)
%r为指数分布的参数,n为产生随机数个数
a=random('unif',0,1,n,1);
b=zeros(n,1); 
for i=1:n;
    b(i)=-log(a(i))/r;
end
s=zeros(n+1,1);
for i=2:n+1;
    s(i)=s(i-1)+b(i-1);
end
t=0:0.0001:s(n+1)+0.0001;
p=length(t);
nt=zeros(1,p);
j=1;
for i=1:p;
    if t(i)<s(j+1)&t(i)>s(j)
        nt(i)=j-1;
    end
    if t(i)>=s(j+1);
        nt(i)=j;
        j=j+1;
    end
end
save time t ;
save ntime nt ;
axis([0 t(p) 0 n+1]);
plot(t,nt) 
grid        
load time  		
load ntime 	
m=length(t);
N=nt(m);
T=t(m);
s=zeros(1,N+1);
k=1;
s(1)=0;
for i=2:m;
    if nt(i)==k;
        s(k+1)=t(i);
        k=k+1;
    end
end
deta=zeros(1,N);
deta(1)=s(2);
for i=2:N;
    deta(i)=deta(i-1)+s(i+1);
end
check=deta(N)/T;
if check<=(N+1.96*sqrt(N/3))/2&check>=(N-1.96*sqrt(N/3))/2;
    disp('检验成功，该序列是泊松分布');
end
if check>(N+1.96*sqrt(N/3))/2||check<(N-1.96*sqrt(N/3))/2;
    disp('检验成功，该序列不是泊松分布!');
end
check
checkdown=(N-1.96*sqrt(N/3))/2
checkup=(N+1.96*sqrt(N/3))/2
%估计参数lamda
lamda=N/T;
sprintf('极大似然参数估计为%f ',lamda)