[R] Estimating parameters of 2 phase Coxian using optim

[Gad Abraham]

Andy Fugard wrote:
> Hi There,
> 
> Perhaps the problem is the line
> 
>          loglik<-log(p %*% expm(Q * y(i)) %*% q)
> 
> You mention that y is a vector but here you're treating it as a  
> function.  Maybe try
> 
>          loglik<-log(p %*% expm(Q * y[i]) %*% q)
> 
> ?
> 
> Don't have a clue about the correctness of the contents of cox2.lik...
> 
> Andy
> 
> 
> On 6 Mar 2007, at 08:54, Laura Hill wrote:
> 
>> Hi,
>>
>> My name is Laura. I'm a PhD student at Queen's University Belfast  
>> and have
>> just started learning R. I was wondering if somebody could help me  
>> to see
>> where I am going wrong in my code for estimating the parameters  
>> [mu1, mu2,
>> lambda1] of a 2-phase Coxian Distribution.
>>
>> cox2.lik<-function(theta, y){
>>     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)
>>
>>     for (i in 1:length(y)){
>>         loglik<-log(p %*% expm(Q * y(i)) %*% q)
>>     return(loglik)}
>>
>>     sumloglik<-sum(loglik)
>>
>>     return(-sumloglik)
>>     }

Just to add my 2 AU cents regarding the for loop:

You're trying to create a vector of log likelihoods to sum up later, but 
that's not what's happening there. Instead, assign an empty vector of 
same length as y, then assign the loglik from each iteration to a 
different cell.

Lastly, there's no need to return anything from a for loop, it's not a 
function.

HTH,
Gad

-- 
Gad Abraham
Department of Mathematics and Statistics
The University of Melbourne
Parkville 3010, Victoria, Australia
email: g.abraham at ms.unimelb.edu.au
web: http://www.ms.unimelb.edu.au/~gabraham