[R] Estimating parameters of 2 phase Coxian using optim

[Laura Hill]

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)
    }

I have installed the Matrix package. y is a vector of 240 survival times. In
the R console I typed

> source("/private/var/automount/users/lhill07/Desktop/cox2.lik.R")

> optim(c(0.5, 0.5, 0.5), cox2.lik, y=y, method="BFGS")

Error: could not find function "y"

Error in expm(Q * y(i)) : error in evaluating the argument 'x' in selecting
a method for function 'expm'

Error in log(p %*% expm(Q * y(i)) %*% q) :
    error in evaluating the argument 'x' in selecting a method for function
'log'


I'm sorry if I have missed something really obvious and I would appreciate
any help that is offered. Hopefully I have given enough information.

Many thanks in advance,

Laura