Topics in Model-Fitting for Omegahat

This page is a worksheet for the various interests within the Omega project in developing software for fitting and using statistical models, describing some of the current work with the main emphasis on what we would like to do and questions that need discussion. Some of the computational material here relates actually to the software in the Numerics/Optimization module.

Software for Fitting Nonlinear Models: Basics

As an initial focus for models, we have started by looking at general models corresponding to statistical criteria such as maximum likelihood, nonlinear least squares, or Bayesian modeling.

We have the basic concept of an ``optimizer'', an iterative fitting method, roughly broken down into 4 methods for the specific algorithm: initialize, refine, test convergence, and wrapup. The following is a possible eval method for an optimizer class, and expecting to be given an object, fitter having implementations for these methods. The continueIteration call will provide a hook for listeners to be notified on each iteration; see below.

 public ParameterSet eval(ParameterSet initialValues) {
   OptimizerAlgorithm fitter =  getFitter();
   ParameterSet currentValue = 
     fitter.initialize(initialValues, this);
   while( !fitter.hasConverged(currentValue, this)) {
       currentValue = fitter.refine(currentValue, this);
       if(!continueIteration(fitter, currentValue))
          break;
     } 
  return(fitter.wrapup(currentValue, this));
 }

Some questions and topics for discussion:

1. The structure for the algorithm object, fitter, seems general enough, at least for standard optimization. Is it? (And for other computational criteria?)

2. We can fill in actual algorithms by two strategies: Java code or native methods. At the moment it looks as if there is little or no Java optimization software that would be suitable, in terms of its structure or its solid, tested nature.

   So, are we happy enough about doing some implementation? If so, what methods are attractive starting points? Criteria would include relative simplicity to get decent numeric performance, and perhaps availability of a non-Java description or implementation that could be used as a good guide.

   The decision to write new algorithms for optimization raises some questions. Details count for these algorithms as much as for any, and we are not a priori interested in developing new computational methods. On the other hand, the structure above does provide opportunities for more flexible, data-oriented methods.

3. On the other hand, what would be the best approach to create native methods? In other words, to ``take apart'' a C, C++, or Fortran procedure sufficiently to allow the individual steps to be done via Java native calls. Are there good candidates here? What are the plus & minus considerations to this approach?

4. The fitting algorithm has to have access to the model definition, plus perhaps derivatives. The Omegahat Models module includes a Formula sub-module that supports a formula language, similar in spirit to that in S. This seems a natural starting point. Are there other ways to describe the model that would be convenient (right now or with some changes)?

5. In particular, what about automatic derivatives? (There is some elementary symbolic differentiation in S that would be easy to copy as a starting point.) Can we benefit from the work in Ampl? Or from other existing software?

6. Nonlinear least-squares seems to work into this general framework without too much difficulty, so we can ask questions 2 and 3 about that methodology as well. Or can we? Slightly farther afield, what about constrained optimization? It's interesting that statisticians have paid so little attention to this technique. Admittedly, the computational methods are somewhat harder to choose and use. Constrained fitting would probably also imply some extensions of the facilities for model definition.

7. Then, Bayesian model fitting. This is a much bigger question, so it seems better to deal with it separately.

Using the Java Event Structure

The continueIteration method provides for ``iteration listeners''. After each step in the iteration, the optimizer calls each action registered in a list of listener actions. The action is a Java method or a suitable Omegahat function, which can do essentially arbitrary computations. A mechanism is provided for running the listener actions in a background thread, in case one wants the fitting not to be blocked. Two types of listeners are likely to be useful:

1. Plots and other summaries that give feedback to the user about the progress of the fitting.

2. Specializations of the fitting method itself, such as tests that are appropriate to the particular class of models, or even to the specific example.

Bayesian Approaches to Models

Advances in the techniques (computational among other) for dealing with the Bayesian approach to statistical models are an important part of the recent history of the field. The same period has shown increased use of Bayesian analysis to challenging practical situations. In the Omega discussions of computations for models, we seem to have consensus that bringing Bayesian techniques into our scope is essential. Since Omega is a joint project, this also means we need to bring in active participants who are strongly interested in Bayesian techniques.

Here are some points that have come up for discussion (there are sure to be others we haven't yet uncovered).

1. At one extreme, the Bayesian module could be completely separate from the other components for models. The assumption so far is that instead we should try to integrate Bayesian and ``classical'' methods as far as possible, so long as the user's view does not suffer. Integration will make tools available for all approaches, providing a more uniform view of modeling than software.

2. A specific question that immediately follows: Does the basic optimizer's fitting loop make sense for Bayesian computation, in particular for Markov-Chain Monte-Carlo(MCMC)? An interpretation for MCMC would be that the iteration continues until a test of stationarity succeeds. The wrapup method might do some more sampling. Is this a sensible structure?
3. One major challenge is to integrate the various forms of user input defining models, again to the extent it makes sense. What do we need to add, if anything, to the existing formula package? The use of graphical representations for models, both as a data structure and as an interface (e.g., Stephan Lauer's work) is popular for Bayesian analysis, though not limited to that. We would like to integrate this with the Omegahat software.

Other Statistical Questions

1. What about other classes of models? Which would benefit from being included in this general structure? It seems that iteratively reweighted least squares in general would fit in. With all such models, we again need to consider the question of how best to incorporate detailed numerical algorithms: if good C/Fortran procedures are publicly available, should we incorporate these or try to convert them to Java and the structure above?
2. So far, we are working with one or more data structures following generally on the ``data frame'' structure in S. What is needed to supplement or as an alternative?