R help - Nov 2011 - How can I pick a matrix from a function? (Out Product of Gradient)

Hi all,

I would like to use optim() to estimate the equation by the log-likelihood
function and gradient function which I had written. I try to use OPG(Out
Product of Gradient) to calculate the Hessian matrix since sometime Hessian
matrix is difficult to calculate. Thus I want to pick the Gradient matrix
from the gradient function.

Moreover, could R show the process of calculation on gradient function I
written by using "optim( )" or other commands?

Thanks,

Yanghao

===========================================================
X <- cbind(rep(1,n),sex,age,yrmarry,children,rating)
dy <- (mydata$y>0)*1

# *********************************************
# Probit model: log-likelihood
# *********************************************
fprobit <- *function*(beta,y,X) {
n <- length(y)
k <- ncol(X)
b <- beta[1:k]
kk <- 2*y-1
 z <- X %*% b
L <- log(pnorm(kk*z))
L <- sum(L)
return(-L)
}

# *********************************************
# Probit model: analytic gradient
# *********************************************
gprobit <- *function*(b,y,X) {
kappa=2*y-1
z <- kappa*(X %*% b)
imr <- kappa*dnorm(z)/pnorm(z)
G<-matrix(ncol=ncol(X),nrow=nrow(X))
*for*(i *in* 1:nrow(imr)){
G[i,]<-imr[i,]*X[i,]
}
g <- apply(G,2,sum)
return(-g)
}

#**********************************************************
# How can I pick the G matrix from this function?
#**********************************************************

############For initial value#####################
xx <- cbind(sex,age,yrmarry,children,rating)
reg1 <- lm(dy~xx)
(b0 <- reg1$coef)
##################################################

(mle <- optim(b0,fprobit,gr=gprobit,
method="BFGS",hessian=T,y=dy,X=X,control=list(trace=T)))

	[[alternative HTML version deleted]]

Hi all,

I would like to use optim() to estimate the equation by the log-likelihood
function and gradient function which I had written. I try to use OPG(Out
Product of Gradient) to calculate the Hessian matrix since sometime Hessian
matrix is difficult to calculate. Thus I want to pick the Gradient matrix
from the gradient function.

Moreover, could R show the process of calculation on gradient function I
written by using "optim( )" or other commands?

Thanks,

Yanghao

===========================================================

X <- cbind(rep(1,n),sex,age,yrmarry,children,rating)
dy <- (mydata$y>0)*1

# *********************************************
# Probit model: log-likelihood
# *********************************************
fprobit <- function(beta,y,X) {
  n  <- length(y)
  k  <- ncol(X)
  b  <- beta[1:k]
  kk <- 2*y-1
  z  <- X %*% b
  L  <- log(pnorm(kk*z))
  L  <- sum(L)
  return(-L)
}

# *********************************************
# Probit model: analytic gradient
# *********************************************
gprobit <- function(b,y,X) {
  kappa=2*y-1
  z <- kappa*(X %*% b)
  imr <- kappa*dnorm(z)/pnorm(z)
  G<-matrix(ncol=ncol(X),nrow=nrow(X))
  for(i in 1:nrow(imr)){
    G[i,]<-imr[i,]*X[i,]
  }
  g <- apply(G,2,sum)
  return(-g)
}

############For initial value#####################
xx   <- cbind(sex,age,yrmarry,children,rating)
reg1 <- lm(dy~xx)
(b0  <- reg1$coef)
##################################################

(mle  <- optim(b0,fprobit,gr=gprobit,
              
method="BFGS",hessian=T,y=dy,X=X,control=list(trace=T)))

--
View this message in context:
http://r.789695.n4.nabble.com/How-can-I-pick-a-matrix-from-a-function-Out-Product-of-Gradient-tp4123965p4123965.html
Sent from the R help mailing list archive at Nabble.com.