[R] Help with faster optimization for large parameter problem

2007-03-02 Thread James Fowler
Hello all,

I have a large parameter problem with the following very simple likelihood 
function:

fn-function(param) {
  x1-param[1:n]
  g1-param[(n+1):(2*n)]
  beta-param[(2*n+1):(2*n+k)]
  sigma2-param[2*n+k+1]^2
  meang1sp-mean(g1[sp])
  mu-beta%*%matrix(x1,1,n)-(g1[sp]-meang1sp)%*%matrix(g1,1,n)
  return(sum((ydc-mu)^2)/(2*sigma2) + n*k*log(sqrt(sigma2)) +
   mean(x1)^2 + mean(g1)^2 + 1000*(x1[1]x1[n]))
}

There's nothing special here -- it's just plain old OLS, except all the 
variables on the right hand side are parameters (only the variable ydc is 
observed).  I have no problems getting this to recover parameter estimates 
from data I myself have generated for smaller problems (e.g. where ydc is 
a k=500 by n=50 matrix and there are 601 parameters).  But the target 
problem with real data will be k=6000 by n=400 with 6801 parameters.

I am currently using optim(method=BFGS) but it is slow.  I can get good 
starting values for x1 and g1 with some svd techniques, and these help me 
generate the starting values for the betas via lm().

I then use optim() on a modified likelihood function to find g1,x1,sigma2 
while holding beta fixed and then use optim() again to find beta while 
holding the other variables fixed.  But eventually, I have to run optim on 
the unmodified likelihood function above and it is very slow, taking 
several days for large problems.

I have also tried metrop() in mcmc, but I find this needs to be very 
close to the mode of the likelihood to be efficient (in fact, MCMCpack's 
metropolis function calls optim first and even requires it to invert the 
hessian before even starting the metropolis algorithm, unless we can 
provide our own covariance matrix).  I will probably use metrop() to 
generate standard errors once I find a mode

In the mean time, I can't help thinking that there is some easy way to 
make this much faster than I am currently doing, especially since the 
likelihood is normal.  I am sure I have missed something obvious so I'd 
very much appreciate any advice you could give on packages in R or code 
that might help.

Thanks!
james

--
James H. Fowler, Associate Professor   web:   http://jhfowler.ucsd.edu
Department of Political Scienceemail: [EMAIL PROTECTED]
University of California, San Diegophone: (858) 534-6807
Social Sciences Building 383, #0521fax:   (858) 534-7130
9500 Gilman Drive, La Jolla, CA 92093

__
R-help@stat.math.ethz.ch 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.


[R] parameter constraints with sem

2006-03-16 Thread James Fowler
I would like to constrain three of the parameters I am estimating using
sem (strucutral equation model package) but do not know how.
Specifically, I would like to use the following constraint:

beta1 + beta2 + beta3 = 1

can I assign variable names as the values for fixed parameters?  If so, I
could define the third parameter to be equal to 1 - beta1 - beta2.

thanks,
james


-
James H. Fowler, Assistant Professor   web:   http://jhfowler.ucdavis.edu
Department of Political Scienceemail: [EMAIL PROTECTED]
University of California, Davisphone: (530) 752-1649
One Shields Avenue fax:   (530) 752-8666
Davis, CA  95616

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html