0

I'm trying to estimate the parameters of my Log-Likelihood function given a set of constraints and using the Newton-Raphson method.

My actual target function is more complex than the one that I will be showing in this question and it has a fairly complex set of constraints to implement (namely, the matrix of the ratios of the estimated coefficients needs to have spectral radius smaller then one.) Because of this constraint, I cannot simply use "maxNR" or "nlm" or "constrOptim" or "optim" and I have decided to implement this constrain as a penalty in my Likelihood function. So far, for the simplified case (in which the spectral radius constraint translates into the constraint that $\alpha/\beta<1$ ) I have written this for my Log-likelihood function: 

neg.loglik.penalty <- function(params, data, opt=TRUE) {
 mu <- params[1]
 alpha <- params[2]
 beta <- params[3]
 t <- sort(data)
 r <- rep(0,length(t))
 for(i in 2:length(t)) {
     r[i] <- exp(-beta*(t[i]-t[i-1]))*(1+r[i-1])
 }
 loglik <- -tail(t,1)*mu
 loglik <- loglik+alpha/beta*sum(exp(-beta*(tail(t,1)-t))-1)
 loglik <- loglik+sum(log(mu+alpha*r))
 loglik <- loglik+if( alpha <0 || beta <0 || alpha/beta <1 ) {-loglik/2} else {0}
 if(!opt) {
     return(list(negloglik=-loglik, mu=mu, alpha=alpha, beta=beta, t=t,
                 r=r))
 }
 else {
     return(-loglik)
 }
}

I would like to know whether this is the right way of doing things and, in the specific, ask for your help with the following issues:

  • shouldn't the value of the penalty be changing with each iteration of the algorithm (instead of being fixed as I have put it)?
  • how should I choose the value for the penalty? In principle, I would like it to be as high as possible since finding a solution outside my constraints implies that my process is non stationary.
  • How does the choice of the value for my penalty affect the shape of the solution space? Won't the algorithm be fooled into going in the wrong direction if I add the wrong penalty value?
g_puffo
  • 613
  • 3
  • 11
  • 21
  • are you looking for statistical or coding answers – rawr Sep 19 '14 at 18:55
  • mostly coding, in the sense that I would rather have my penalty function increasing in value with each iteration of the algorithm rather than fixed but I have no idea of how to do it in R; also, I'm not sure whether this is the right way of implementing a penalty function in the first place. However, given that I'm not an expert in function optimization, if there is anything wrong with my approach I would love to hear it. – g_puffo Sep 19 '14 at 19:00

0 Answers0