From: laila khalfan <byblus412_at_hotmail.com>

Date: Thu, 03 Jul 2008 19:36:40 +0000

R-help_at_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. Received on Thu 03 Jul 2008 - 19:39:18 GMT

Date: Thu, 03 Jul 2008 19:36:40 +0000

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_at_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. Received on Thu 03 Jul 2008 - 19:39:18 GMT

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.2.0, at Thu 03 Jul 2008 - 20:31:48 GMT.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*