mardi 28 mars 2017

I am getting error msg as Error in if (abs(p11_p10) < 1e-04 & abs(alpha1_alpha0) < 1e-04) { : missing value where TRUE/FALSE needed

I am getting error msg as mentioned above. I am trying to do simulation and estimation and for that I have wrote the above program. Kindly help me to solve this problem.

p1=0.2;p2=0.3;p3=1-p1-p2;alpha=0.4

n=10;k=1;count=0

p11_p10=0.01

alpha1_alpha0=0.01

alpha0=0.5;p10=0.3;p11_p10=0.1;alpha1_alpha0=0.1

repeat

{

repeat{




    x1<-rgeom(n-k,p3/(1-p2))




    y1<-rnbinom(n-k,x1+1,1-p2)




    p3star<-1-alpha*p1-alpha*p2




    xstar<-rgeom(k,p3star/(1-alpha*p2))




    ystar<-rnbinom(k,xstar+1,1-alpha*p2)





    xsamp<-c(x1,xstar)




    ysamp<-c(y1,ystar)




    x=rep(0,n)




    y=rep(0,n)




    z=rep(0,n)




for(i in 1:n)




        {




          x[i]=(xsamp[i]+ysamp[i])*(alpha0^(xsamp[i]+ysamp[i]-1))




          y[i]=alpha0^(xsamp[i]+ysamp[i])




          z[i]=((xsamp[i]+ysamp[i])*(xsamp[i]+ysamp[i]-1)*
                (alpha0^(xsamp[i]+ysamp[i]-2)))




        }




          x;y;z




        num=sum(x)




        den=sum(y)




        num1=sum(z)




        term=num1/den




        if(term>0)




        {




          T1<-sum(xsamp); T2<-sum(ysamp)




        if (T1!=0 || T2!=0)




            {




            break




            }




        }




    }

print(T1)

print(T2)

count=0

repeat




    {




    count=count+1




    b1=(k/n)




    A=alpha0*(((T1+T2)^2)+n*(T1+T2))




    B=(-T1*((alpha0+1)*(T1+T2)+(n+alpha0-1)))




    C=(T1^2)





    t1_p1_alpha<-A*(p10^2)+(B*p10)+C




    t2_p1_alpha= -((p10*(T1+T2))/(T1-alpha0*(T1+T2)*p10))+(num/den)






    t1dash_wrt_p1=2*A*p10+B




    t1dash_wrt_alpha=(p10^2)*((T1+T2)^2)+n*(p10^2)*(T1+T2)-p10*T1*(T1+T2+1)




    t2dash_wrt_p1=(((-T1*(T1+T2))/(T1-alpha0*(T1+T2)*p10))^2)




    t2dash_wrt_alpha= (-((p10*(T1+T2))^2)/((T1-alpha0*(T1+T2)*p10)^2))+(num1/den)-((num^2)/(den^2))




    a<-t1_p1_alpha;b<-t1dash_wrt_p1;c<-t1dash_wrt_alpha




    d<-t2_p1_alpha;e<-t2dash_wrt_p1;f<-t2dash_wrt_alpha






    alpha1=alpha0+(1/c)*(-a-b*((a*f-d*c)/(e*c-b*f)))
    p11=p10+((a*f-d*c)/(e*c-b*f))





    p11_p10=p11-p10




    alpha1_alpha0<-alpha1-alpha0





    if(abs(p11_p10)<0.0001 & abs(alpha1_alpha0)<0.0001)




            {




             break




            print(p11)




            print(alpha1)




            }




            else




            {




            p10<-p11




            alpha0<-alpha1




            }




    if (count>50){break}




    }

}

Aucun commentaire:

Enregistrer un commentaire