[R] R computing speed

From: Carlo Fezzi <C.fezzi_at_uea.ac.uk>
Date: Tue, 11 Dec 2007 14:55:51 -0000


Dear helpers,

I am using R version 2.5.1 to estimate a multinomial logit model using my own maximum likelihood function (I work with share data and the default function of R cannot deal with that).

However, the computer (I have an Athlon XP 3200+ with 512 GB ram) takes quite a while to estimate the model.

With 3 categories, 5 explanatory variables and roughly 5000 observations it takes 2-3 min. For 10 categories and 10 explanatory variables (still 5000 obs) more than 1 hour.

Is there any way I can speed up this process? (Modifying the code or modifying some R options maybe?)

I would be really grateful if anybody could help me with this issue, I attach my code below.

Many thanks,

Carlo



Carlo Fezzi

Centre for Social and Economic Research on the Global Environment (CSERGE),
School of Environmental Sciences,
University of East Anglia,
Norwich, NR4 7TJ
United Kingdom.


# MULTILOGIT

# This function computes the estimates of a multinomial logit model

# inputs: a matrix vector of 1 and 0 (y) or of shares
# a matrix of regressors (x) - MUST HAVE COLUMN NAMES! -
# names of the variables, default = colnames(x)
# optimization methods, default = 'BFGS'

#		base category, default = 1
#		restrictions, default = NULL

# weights, default all equal to 1

# outputs: an object of class "multilogit.c"

# McFadden D. (1974) "Conditional logit analysis of qualitative choice
behavior", in Zarembka P. (ed.), Frontiers in Econometrics, Academic Press.

multilogit.c <- function(y, xi, xi.names = colnames(xi), c.base=1, rest=NULL, w = rep(1,nrow(y)), method='BFGS') {         

	n.obs <- sum(w)
	xi<-cbind(1,xi)
	colnames(xi)[1]<-"Intercept"

	nx<-ncol(xi)
	ny<-ncol(y)
	
	beta<-numeric(nx*ny)
	
	negll<- function(beta,y,xi)
	{
		beta[rest]<-0
		beta[(((c.base-1)*nx)+1):(c.base*nx)]<-0
		lli <- y  * (xi%*%matrix(beta,nx,ny) - log ( apply(exp(
xi%*%matrix(beta,nx,ny)) ,1,sum ) ) 	)
		lli<-lli*w
		-sum(lli)
	}

	pi<- apply((y*w),2,mean)/mean(w)
	
	ll0 <- (t(pi)%*%log(pi))*sum(w)
	
	result<-c(	optim(par = rep(0,nx*(ny)), fn = negll, y=y, xi=xi,
hessian=T, method=method),
			list(varnames=xi.names, rest=rest, nx=nx, ny=ny,
npar=nx*(ny-1)-length(rest), ll0=ll0, 	pi=pi, xi=xi,
n.obs=n.obs,c.base=c.base,w=w))         
	result$par <- result$par[-(((c.base-1)*nx)+1):-(c.base*nx)]
	result$hessian <-

result$hessian[-(((c.base-1)*nx)+1):-(c.base*nx),-(((c.base-1)*nx)+1):-(c.ba se*nx)]
	class(result)<-"multilogit.c"
	return(result)

}

R-help_at_r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Tue 11 Dec 2007 - 14:59:01 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Tue 11 Dec 2007 - 15:30:20 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.