pcancer <- read.table(file="cancerp.dat", header=T)

MixExp <- function(parametros, tipo, Ti){
J <- max(tipo)
nt <- length(tipo)
C.ij <- matrix(1:(nt * J), ncol = J) * 0
        for(i in 1:J) {
                C.ij[, i][tipo == i] <- 1
        }
uno <- rep(1, J)
ci <- C.ij %*% uno
Tij <- matrix(1:(nt * J), ncol = J) * 0 + 1
T.ij <- Tij * Ti
deltas <- parametros[1:(J-1)]
betas <- parametros[J:(2*J - 1)]
dj.u <- matrix(1:(J * nt), ncol = J) * 0
Uij <- rep(1, nt)
for(i in 1:(J - 1)) {
                                dj.u[, i] <- Uij * deltas[i]
                        }
bj.v <- matrix(1:(J * nt), ncol = J)
Vij <- rep(1, nt)
for(i in 1:J) {
                                bj.v[, i] <- Vij * betas[i]
                        }
p.ij <- matrix(1:(J * nt), ncol = J) * 0
exp.dj.u <- exp(dj.u)
sum.exps.i <- exp.dj.u %*% uno
for(i in 1:J) {
                        p.ij[, i] <- exp.dj.u[, i]/(sum.exps.i)
                }
l.ij <- exp(bj.v)
A <- sum(C.ij * (log(p.ij) + log(l.ij) - (l.ij * T.ij)))
F0.ij <- 1 - exp( - l.ij * T.ij)
p.F0.ij <- p.ij * F0.ij
sumJ.p.F0.ij <- p.F0.ij %*% uno
B <- sum((1 - ci) * log(1 - sumJ.p.F0.ij))
log.lik <- A + B
-1*log.lik
}


mle.null <- nlm(MixExp, c(0,rep(-3, 2)), tipo=pcancer$tipo,
                Ti=pcancer$tiempo, hessian=TRUE)