1
$\begingroup$

I need to know how can we estimate the tuning parameter in penalized likelihood?

I write my own code but there is a mistake I could not find it.!! Could you please help, I do my best to make the code clear..

Given the likelihood that

$l(\beta)= \sum_i y_i log (\pi_i)+(1-y_i)log (1-\pi) -1/2 \lambda \sum_j {\beta_j}^2$

where $\pi_i= logit^{-1} (x\beta$) and $y_i$ belong to (0,1)

I need to chose the lambda using cross validation ; my question is that

  • I need to include the Code B inside Code A. then I can run the code.
  • is the code correct or no. Is there any way to make it faster.\

Code A : The Cross Validation to choose $\lambda$

   PE=c()
lambda=c(0.005,0.001,0.1,0.5,0.2)

for(j in 1:length(lambda)){
print(lambda[j])
tem.PE=c()
for(i in 1:2014){
cv.y=y[-i]
cv.x=x[-i,]
tem.lambda=lambda[j]
## I need to put the function /WNRfit/  here to calculate the Coefficients -Betas-
f.hat=fit$fitted
#   print(f.hat)
tem.PE[i]=mean(cv.y-f.hat)
rm(cv.y,cv.x,fit,f.hat,tem.lambda)
}
PE[j]=mean(tem.PE)
rm(tem.PE)
}

Code B : The Newton raphson for estimate beta $\hat{\beta}$

   WNRfit <- function(y,X,m=rep(1,length(y)),start=0,n.iter=10,tol=1e-4,verbose=TRUE) {
   W=matrix(0,nrow=length(y),ncol=length(y))
   diag(W)<-m
  n.pred <- ncol(X)
  b <-  matrix(NA,ncol=n.iter,nrow=n.pred)
  b[,1] <- start
  for (i in 2:n.iter) {
  if (verbose) cat(i,"\n")
  p <- plogis(X %*% b[,i-1])
  v.2 <- diag(c(p*(1-p)))
  score.2 <- t(X) %*% W %*% (y - p) # score function
  increm <- solve(t(X)%*%W%*%v.2%*%W%*%X) 
  b[,i] <- b[,i-1]+increm%*%score.2
  if (all(abs(b[,i]-b[,i-1]) < tol))

   return(b)
   fitted=exp(X%*%b[,10])/(1+exp(X%*%b[,10]))

   }
  print(b)
  }
$\endgroup$
  • 1
    $\begingroup$ Please see SO help on how to write a minimal, complete, verifiable example. In particular you need to show why you think you've made a mistake. In addition, making code clear entails explaining what it's doing in comments - especially important when not all readers are familiar with the programming language you're using. $\endgroup$ – Scortchi Oct 13 '15 at 10:16
3
$\begingroup$

It is better to choose a method before writing the code. To answer your overall question, I prefer to choose $\lambda$ to optimize the "effective AIC" when can be stated as the penalized model likelihood ratio $\chi^2$ minus twice the effective degrees of freedom. This is discussed in my book Regression Modeling Strategies and details, with R code, may be found at http://biostat.mc.vanderbilt.edu/wiki/pub/Main/FHHandouts/iscb98.pdf

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.