Can you please help me to estimate the parameters of a 2phase coxian distribution using R. I am Nelson Sibanda currently a student doing Msc Operations Research and Statistics  at National University of Science and Technology, Bulawayo, Zimbabwe. The following code were tried  and the data that l am trying fit is the length of stay of patients in a hospital from 1day to 15 days . Here the frequency table
Length of  stay 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 
Number Of  Patients 1827 835 257 200 195 130 101 76 51 44 31 22 18 8 22 



On Friday, January 31, 2014 4:01 PM, nelson Sibanda <nelsonsibbs78@yahoo.com> wrote:
 



On Saturday, January 18, 2014 9:20 PM, nelsonsibbs78 <nelsonsibbs78@yahoo.com> wrote:


Can you pliz  send me the R code for estimating parameters of a 2 phase distribution given some data to fit. Im nelson sibanda , a final year student at Nust

  

 x<-c(1827,835,257,200,195,130,101,76,51,44,31,22,18,8,22)
cox2.lik<-function(theta){
stopifnot(is.numeric(theta))
stopifnot(length(theta) == 3)
mu1<-theta[1]
mu2<-theta[2]
lambda1<-theta[3]
p<-matrix(c(1, 0), nrow=1, ncol=2)
Q<-matrix(c(-(lambda1 + mu1), 0, lambda1, -mu2), nrow=2, ncol=2)
q<-matrix(c(mu1, mu2), nrow=2, ncol=1)
loglik<-rep(0,length(x))
for(i in 1:length(x)){
loglik[i]<-log(p%*%expm(Q*x[i])%*%q)
sumloglik<-sum(loglik[i])
return(-sumloglik)
}}
i <- length(x)
mu1.start <- mean(sort(x)[seq(along = x) <= i / 2])
mu2.start <- mean(sort(x)[seq(along = x) <= i / 2])
lambda1.start <- var(sort(x)[seq(along = x) <= i / 2])
theta.start <- c(mu1.start, mu2.start, lambda1.start)
out <- nlm(cox2.lik, theta.start, print.level = 2,
fscale = length(x))
out <- nlm(cox2.lik, out$estimate, print.level = 2,
fscale = length(x), hessian = TRUE)

	[[alternative HTML version deleted]]

