Re: [Rd] optim: why is REPORT not used in SANN?

2008-04-06 Thread Thomas Petzoldt
Martin Maechler wrote:
>> "TP" == Thomas Petzoldt <[EMAIL PROTECTED]>
>> on Sun, 16 Mar 2008 13:50:55 +0100 writes:
> 
> TP> Hello, I wonder why the control parameter REPORT is not
> TP> supported by method SANN. Looking into optim.c I found
> TP> an internal constant:
> 
> TP> #define STEPS 100
> 
> TP> ... and decreasing this to 10 helped me fine-tuning the
> TP> annealing parameters in an actual problem.
> 
> TP> Is there any reason why not passing nREPORT to samin and
> TP> setting something like:
> 
> TP> STEPS = nREPORT / tmax
> 
> Sorry to reply late (but then, rather than never ..).
> 
> You ask for reasons... I see/guess :
> 
> - the  SANN  method also was contributed from "outside" 
>   (as ?optim mentions); and the original authors may not have
>   seen a use for such more flexible monitoring.
> 
> - the R core members are probably not using 'samin' very often
> 
> - If there is a need you can write the function you are
>   optimizing in a way that it prints info.
> 
> - Nobody has contributed a well-tested patch against R-devel to
>   both code and documentation
>   which would implement your proposal ___ BACK COMPATIBLY __
>   (i.e. the default for SANN should remain to print every 100th;
>and this differs from the default for BFGS where the default
>   'REPORT' leads to output every 10th eval).
> 
> Regards,
> Martin

Well, I see. The reasons are obviously more or less historical and not 
fundamental. I can, of course, contribute an idea for a modifications of 
samin and do_optim in optim.c. However, to convert my personal hack into 
a "well-tested patch" a few additional considerations have to be 
undertaken, especially about back-compatibility:

1) the patch requires to pass the additional argument nREPORT to 
function "samin".

- Is this still back-compatible? Is it likely that other functions (e.g. 
in packages call samin directly?

- if yes (i.e. direct call is likely), what is the preferred way to 
ensure compatibility? Rename samin to samin2 and add a new 
"compatibility function" function samin that calls samin2?

- the reporting interval of samin is STEPS * tmax where
   - tmax is as documented "number of function evaluations at each 
temperature" (10) and
   - STEPS is a hard-coded constant: #define STEPS 10
   - this means that reporting is every 10 per temperature.

2) if one starts patching SANN, then one migth also think about
"Nelder-Mead" and "CG" with totally different reporting, but smaller 
(!) reporting schedule.

In contrast to this, SANN that is used for difficult systems *and* that 
(see optim.Rd) "depends critically on the settings of the control 
parameters" has a rather long reporting interval.

Thomas P.

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


[Rd] manual 'Writing R Extensions': bug in example

2008-04-06 Thread Simon Anders
Hi,

I would like to report a small bug in one of the code examples given in 
the 'Writing R extensions' manual. Section 5.9.2 ("Calling .External") 
defines the function "showArgs" in order to demonstrate how to access in 
C the pair-listed arguments passed by .External. The function aims at 
printing echoing all passed arguments to the screen but only prints 
every second argument.

The loop to go through the pair list starts as follows

   args = CDR(args); /* skip 'name' */
   for(i = 0; args != R_NilValue; i++, args = CDR(args)) {
  args = CDR(args);/* <--- This line is wrong. */

As you can see, "args = CDR(args)" is called twice, resulting in every 
other argument skipped. The marked line should be deleted.

Furthermore, the function requires all arguments to be named, as it 
gives an error if TAG returns NIL in the line

name = CHAR(PRINTNAME(TAG(args)));

It might improve the example to introduce here a test "if( TAG(args) != 
R_NilValue )".

May I suggest to the maintainer of the manual to change the example to 
the following code? (Changes to original: lines 1, 2, 15, 17, 18 added; 
one line deleted after line 13.)

  1  #include 
  2  #include 
  3  #include 
  4
  5  SEXP showArgs(SEXP args)
  6  {
  7int i, nargs;
  8Rcomplex cpl;
  9const char *name;
10SEXP el;
11
12args = CDR(args); /* skip 'name' */
13for(i = 0; args != R_NilValue; i++, args = CDR(args)) {
14
15  if( TAG(args) != R_NilValue )
16 name = CHAR(PRINTNAME(TAG(args)));
17  else
18 name = "";
19  switch(TYPEOF(CAR(args))) {
20  case REALSXP:
21Rprintf("[%d] '%s' %f\n", i+1, name, REAL(CAR(args))[0]);
22break;
23  case LGLSXP:
24  case INTSXP:
25Rprintf("[%d] '%s' %d\n", i+1, name, INTEGER(CAR(args))[0]);
26break;
27  case CPLXSXP:
28cpl = COMPLEX(CAR(args))[0];
29Rprintf("[%d] '%s' %f + %fi\n", i+1, name, cpl.r, cpl.i);
30break;
31  case STRSXP:
32Rprintf("[%d] '%s' %s\n", i+1, name,
33   CHAR(STRING_ELT(CAR(args), 0)));
34break;
35  default:
36Rprintf("[%d] '%s' R type\n", i+1, name);
37  }
38}
39return(R_NilValue);
40  }
41

Best regards
   Simon Anders


+---
| Dr. Simon Anders, Dipl. Phys.
| European Bioinformatics Institute, Hinxton, Cambridgeshire, UK
| office phone +44-1223-494478, mobile phone +44-7505-841692
| preferred (permanent) e-mail: [EMAIL PROTECTED]

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] Matrix multiplication in a C code

2008-04-06 Thread Prof Brian Ripley
On Sun, 6 Apr 2008, Mathieu Ribatet wrote:

>Dear list members,
>
> I've got a small question on matrix multiplications in a C code. Because
> of a really cpu demanding likelihood, I had to use a C code within an R
> function wrapper. I'm pretty sure that there is already one good code
> for matrix multiplication in C - maybe in the R source code itself - but
> I wasn't able to find it.
>
> Anyone as an idea on how to handle matrix multiplications?

Well, R does use fast C code where available: see matmult in 
src/main/array.c which calls the BLAS.  It works best if you have an 
optimized BLAS: otherwise it may not be as fast as the simple version used 
when there are NAs present.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


[Rd] Matrix multiplication in a C code

2008-04-06 Thread Mathieu Ribatet
Dear list members,

I've got a small question on matrix multiplications in a C code. Because 
of a really cpu demanding likelihood, I had to use a C code within an R 
function wrapper. I'm pretty sure that there is already one good code 
for matrix multiplication in C - maybe in the R source code itself - but 
I wasn't able to find it.

Anyone as an idea on how to handle matrix multiplications?

Best,
Mathieu

-- 
Institute of Mathematics
Ecole Polytechnique Fédérale de Lausanne
STAT-IMA-FSB-EPFL, Station 8
CH-1015 Lausanne   Switzerland
http://stat.epfl.ch/
Tel: + 41 (0)21 693 7907

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] [R] How to improve the "OPTIM" results

2008-04-06 Thread Spencer Graves
 

  What would people suggest someone do to improve the optimization 
capabilities in R?  A few ideas are mentioned below.  I just got 578 
hits from RSiteSearch("optimization", "fun"), and I wonder if there's a 
better way to get a reasonable overview of what's currently available in 
R -- and outside R?  My initial thought was to develop an enhanced 
version of 'optim' that would check for negative curvature (via 
indefinite Hessians), do auto-scaling (see 
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/125518.html), and output 
an object of a class like 'mle' or a new class like 'optimMLE' for which 
appropriate methods could be written.  I have to often encountered 
problems with overparameterization, which could be fairly easily 
diagnosed via a singular or indefinite hessian.  It would help me, and I 
believe other users, if functions like 'optim' and 'nls' would include 
simple diagnostics like this -- rather than returning a cryptic error or 
warning message. 

  However, the comment from Hans Borchers (below) raises a related 
question:  What might be the best way to build a collaboration with 
existing optimization initiatives, possibly making R a platform of 
choice for making it easy for users to access and compare alternative 
optimization methods? 

  Thanks,
  Spencer 

Hans W Borchers wrote:
>> MORE GENERAL OPTIM ISSUES
>>
>>   I'm considering creating a package 'optimMLE' that would automate 
>> some of this and package it with common 'methods' that would assume that 
>> sum(fn(...)) was either a log(likelihood) or the negative of a 
>> log(likelihood), etc.  However, before I do, I need to make more 
>> progress on some of my other commitments, review RSiteSearch("optim", 
>> "fun") to make sure I'm not duplicating something that already exists, 
>> etc.  If anyone is interested in collaborating on this, please contact 
>> me off-line. 
>>
>>   Hope this helps. 
>>   Spencer
>> 
>
> Thanks for your tips on using `optim()'. I believe `optim' is a reasonable 
> good
> implementation of numerical optimization techniques such as quasi-Newton BFGS 
> or
> conjugate gradients. Maybe methods using modern line searches or trust regions
> can follow someday.
>
> Everybody applying `optim' should be aware that it is a *Local Optimization*
> (LO) approach. What you describe appears to be a step towards *Global
> Optimization* (GO) in R. And indeed more and more requests in to the R-help 
> list
> are about `optim' as a tool for global optimization, not always being fully
> aware of the differences.
>
> I am wondering whether it would be more useful to provide one or two global
> optimization approaches to the R community. And as this is an active research
> area, there are many and with different advantages and drawbacks.
>
> I would like to propose IPOPT as one of he newer and more advanced methods for
> global optimization. This powerful software package is open source and 
> available
> through the COIN-OR initiative and its Web pages:
>
> http://www.coin-or.org/Ipopt/documentation/
>
> ``Ipopt (Interior Point OPTimizer, pronounced I-P-Opt) is a software 
> package 
> for large-scale nonlinear optimization. Ipopt is written in C++ and is 
> released as open source code under the Common Public License (CPL). It is
> available from the COIN-OR initiative. The code has been written by Carl
> Laird (Carnegie Mellon University) and Andreas Wachter (IBM's T.J. Watson 
> Research Center), who is the COIN project leader for Ipopt.''
>
> PERHAPS the COIN project would agree for IPOPT to be integrated into the open
> source project R as a package. For testing it right now there is an AMPL-based
> interface to IPOPT at the NEOS solver:
>
> http://neos.mcs.anl.gov/neos/solvers/nco:Ipopt/AMPL.html
>
> There may be other rewarding projects in the COIN-OR initiative, such as
> `SYMPHONY' for solving mixed-integer linear programs (MILP, stronger than 
> `glpk'
> and `lp-solve'), or the BONMIN open source *non-linear* mixed integer
> programming (MINLP) solver. I could imagine R to be a good platform for
> integrating some of these algorithms.
>
> 
> Hans W. Borchers
> Control and Optimization Group
> ABB Corporate Research Germany
> [EMAIL PROTECTED]
>
> __
> [EMAIL PROTECTED] 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-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] [R] How to improve the "OPTIM" results

2008-04-06 Thread Dirk Eddelbuettel

On 6 April 2008 at 08:37, Spencer Graves wrote:
|  
[...]
|   However, the comment from Hans Borchers (below) raises a related 
| question:  What might be the best way to build a collaboration with 
| existing optimization initiatives, possibly making R a platform of 
| choice for making it easy for users to access and compare alternative 
| optimization methods? 

Are you (used in the plural) aware of Kurt's more recent advances into OR?
There are new-ish packages RSymphony, Rglpk and other things at CRAN.

Dirk

--
Three out of two people have difficulties with fractions.

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] [R] How to improve the "OPTIM" results

2008-04-06 Thread Gabor Grothendieck
On Sun, Apr 6, 2008 at 12:59 PM, Kurt Hornik <[EMAIL PROTECTED]> wrote:
> > Dirk Eddelbuettel writes:
>
> > On 6 April 2008 at 08:37, Spencer Graves wrote:
> > | 
> > [...]
> > |   However, the comment from Hans Borchers (below) raises a related
> > | question:  What might be the best way to build a collaboration with
> > | existing optimization initiatives, possibly making R a platform of
> > | choice for making it easy for users to access and compare alternative
> > | optimization methods?
>
> > Are you (used in the plural) aware of Kurt's more recent advances into
> > OR?  There are new-ish packages RSymphony, Rglpk and other things at
> > CRAN.
>
> Not only by me: there's also Rcplex and TSP and there are at least three
> more COIN-OR package interfaces in the pipeline.
>

It would be nice if there were a CRAN Task View to pull all this together.

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] [R] How to improve the "OPTIM" results

2008-04-06 Thread Kurt Hornik
> Dirk Eddelbuettel writes:

> On 6 April 2008 at 08:37, Spencer Graves wrote:
> |  
> [...]
> |   However, the comment from Hans Borchers (below) raises a related 
> | question:  What might be the best way to build a collaboration with 
> | existing optimization initiatives, possibly making R a platform of 
> | choice for making it easy for users to access and compare alternative 
> | optimization methods? 

> Are you (used in the plural) aware of Kurt's more recent advances into
> OR?  There are new-ish packages RSymphony, Rglpk and other things at
> CRAN.

Not only by me: there's also Rcplex and TSP and there are at least three
more COIN-OR package interfaces in the pipeline.

-k

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


Re: [Rd] [R] How to improve the "OPTIM" results

2008-04-06 Thread Spencer Graves


Gabor Grothendieck wrote:
> On Sun, Apr 6, 2008 at 12:59 PM, Kurt Hornik <[EMAIL PROTECTED]> wrote:
>   
>>> Dirk Eddelbuettel writes:
>>>   
>>> On 6 April 2008 at 08:37, Spencer Graves wrote:
>>> | 
>>> [...]
>>> |   However, the comment from Hans Borchers (below) raises a related
>>> | question:  What might be the best way to build a collaboration with
>>> | existing optimization initiatives, possibly making R a platform of
>>> | choice for making it easy for users to access and compare alternative
>>> | optimization methods?
>>>   
>>> Are you (used in the plural) aware of Kurt's more recent advances into
>>> OR?  There are new-ish packages RSymphony, Rglpk and other things at
>>> CRAN.
>>>   
>> Not only by me: there's also Rcplex and TSP and there are at least three
>> more COIN-OR package interfaces in the pipeline.
>>
>> 
>
> It would be nice if there were a CRAN Task View to pull all this together.
>   
How about R Wiki entries also?  The enormous advantage of a Wiki is that 
anyone can make contributions, subject to moderation when things get 
extreme.  I'd like to see copies of all the CRAN Task Views ported to 
the R Wiki, with emails sent to the maintainers of all the packages 
mentioned, inviting updates, comparisons, clarifications, 
amplifications, etc.  Similarly, I'd like to see people posting 
questions to R Wiki accompanied by a post to R-Help.  Then someone else 
would edit the question with an answer.  Officially designated R Wiki 
editors could notify people whose posts seem inappropriate for whatever 
reason that their entries may be removed unless clarified, etc.  With 
the right selection of "editors" and a little time, I think the R Wiki 
could become the most valuable documentation available on R. 
> __
> R-devel@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>

__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel