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]]