[R] Simulated annealing using optim()
Ben Bolker
bolker at zoo.ufl.edu
Fri Apr 13 23:07:55 CEST 2007
<otalora <at> ugr.es> writes:
> Finally, I tried to write my "gr" function to compute new random test
> configurations. Here I found a problem: In the simulated annealing
> algorithm, changes between succesive iterations during the optimization
> must be progressively smaller. The "distance" from the "current"
> configuration to the next one to be tested must be somehow proportional to
> the "temperature" but I don't see how to acces the value of the
> "temperature" from my gr function. Unfortunately, this is not described in
> the help page and the examples shown don't need a gr function (in the case
> of the "wild" function example) or uses an "temperature" independent
> function (in the case of the Traveling salesman problem example).
>
> I had a look at the C source code for optim(), and I confirmed that the
> shrinking "temperature" is used by the default Gaussian Markov kernel
> (under the name "scale") but I didn't found any clue on how to access this
> value from the user defined gr function (probably because of my lack of
> knowledge of the R API).
Unfortunately, it looks at the moment optim assumes that "gr" only
has a single argument -- the parameter vector. It would take some
hacking, which I'm (a) not quite sure how to do and (b) really not
sure how to handle without breaking backward compatibility
(i.e. how do you check in C code to see whether an R function
has one or two arguments?)
>From optim.c:
static void genptry(int n, double *p, double *ptry, double scale, void *ex)
{
[SNIP]
/* if user-defined call */
/* set x equal to the scaled parameter vector */
SETCADR(OS->R_gcall, x); /* set argument of "gr" to x */
PROTECT_WITH_INDEX(s = eval(OS->R_gcall, OS->R_env), &ipx);
/* call the "gr" function */
/* etc. */
}
else { /* default Gaussian Markov kernel */
for (i = 0; i < n; i++)
ptry[i] = p[i] + scale * norm_rand(); /* new candidate point */
}
}
Possibly useful:
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/21978.html
further discussion should probably go to R-devel ...
More information about the R-help
mailing list