HiI have a specific sample coming from a gamma(alpha,theta1) distribution and 
then divided into two parts first part follows a gamma(alpha,theta1) the second 
is gamma(alpha,theta2) then I would like to find the mle`s for theta1 and 
theta2 which I found. Now I would like to simulate those estimates 500 or 1000 
times.I tried for loop but it did not work It wont do the loop the problem is 
that I need to evaluate n1 which is the number of units in the first part. n1 
could be different each time. here is the code 
r<-100n<-100shape<-2theta1<-exp(1)theta2<-exp(.5)m0<- function(XX)  #a function 
that generates the estimates{    loglik<-function(xx,alpha,theta1,theta2)    
-1*( -r*lgamma(alpha)-alpha*n1*log(theta1)-alpha*(r-n1)*log(theta2)+(alpha-1)   
 
*sum(log(Ti))+(alpha-1)*sum(log(Tj-tau+(theta2/theta1)*tau))-(1/theta1)*sum(Ti)-
    (1/theta2)*sum(Tj-tau+(theta2/theta1)*tau)+(n-r)*log(1-pgamma((max(Tj)-tau+ 
   (theta2/theta1)*tau)/theta2,alpha,1)))    V<-mle2(minuslogl = loglik, star!
 t = list(alpha= 2, theta1= 3, theta2= 2), data = list(size = 100))    
Est<-coef(V)}estimates<-matrix(,)for (k in 1:2){      
X<-rgamma(n,shape,scale=theta1)      Xs<-sort(X)      tau<-5      for (i in 
1:n) {        if (tau-Xs[i]>0)        n1=i        }      n1      X1<-Xs[1:n1]   
   Ti<-X1      u=n1+1      X2<-Xs[u:n]      X3<-X2*theta2/theta1      
Tj<-X3+tau-tau*theta2/theta1      c1<-matrix(Ti,ncol=1)      
c2<-matrix(Tj,ncol=1)      cc<-data.frame(rbind(c1,c2))[,1]      cc             
                            # the special sample that I need to find the mle`s 
for      estimates<- as.data.frame(t(m0(cc)))      }estimates Thanks in 
advanceLaila

Have fun while connecting on Messenger! Click here to learn more. 

Have fun while connecting on Messenger! Click here to learn more. 
_________________________________________________________________


        [[alternative HTML version deleted]]

______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Reply via email to