I am a full-time consultant and provide services related to the design, implementation and deployment of mathematical programming, optimization and data-science applications. I also teach courses and workshops. Usually I cannot blog about projects I am doing, but there are many technical notes I'd like to share. Not in the least so I have an easy way to search and find them again myself. You can reach me at erwin@amsterdamoptimization.com.

Thursday, July 30, 2009

I have a .lp file which takes very long time (more than 40 hours without finish) when it is run on CPLEX's Interactive Optimizer version 10.2.

As I only have CPLEX v.10.2, I was wondering if someone could run the .lp file for me on the latest version of CPLEX (or Gurobi or other solvers), and then tell me the results.

After you finish running the case, please tell me: (1) The software you use and the software version. (2) How long does it take to finish? (3) The spec. of the computer that are used to run the .lp file, and whether the multi-thread function is turned on. (4) The final values of all the variables. (5) Why the .lp file takes so long?

This is a MIP model. An LP file (or MPS file) is a bad way to look at a MIP problem and give advice about better formulations. Difficult MIP models should always be implemented in a modeling language so you can look at the model, actually understand it and think about it.

Another good example how not to disseminate models is http://www.gamsworld.org/performance/index.htm. The models are totally obscured so they loose any meaning, and don’t allow any reasoning why a certain model gives performance problems for certain solvers.

When we want to learn something about the performance of a model I would argue that it is not helpful to use a format (LP/MPS) that hides the structure and that prohibits reformulations. This is especially the case for discrete models, but also for non-linear models.

1. I have this maximum likelihood problem. When I maximize the likelihood function, I get some results; when I take the log of the likelihood function and maximize it, I get different results. I suppose this is due to the other problem I seem to have, which is too many local maxima. Is this the typical case? Or is there something fundamentally wrong in the way I am posing the problem?

2. Econometricians and statisticians like to use log-likelihood functions instead of likelihood functions. The official argument is that "it is easier to manipulate/optimize". Yes, it is easier to derive, for humans. But it transforms a function that is bounded into something unbounded. If the optimization is going to be made by a computer does that really help or does that make it more difficult for the computer/algorithm? I want to do whatever is the more reliable technique. I hope this is not one of those "case-by-case" things.

Some points:

Taking the log can make the likelihood function better behaved. Some of the very large values that can be formed in sub-expressions can often be prevented this way.

Many good NLP solvers allow for bounds on the variables that allow you to steer the NLP solver away from the regions where the function can not be evaluated.

Use an alternative method to find a good starting point, e.g. a method of moments estimate.

Some of these problems can be difficult. Mixture models are a good example.

Sunday, July 26, 2009

he Excel plug-in allows only for one model. Data Envelopment Analysis models however consist of a number of small LP models: one for each DMU (Decision-Making Unit). However the individual LP models are very small. So it is reasonable to combine all small LP’s into one bigger LP. Basically we form:

Thursday, July 23, 2009

GDXXRW is a tool to exchange (read/write) data between GDX files and Excel spreadsheet files. The option CheckDate is helpful to speed up thing when importing data from Excel. If the GDX file is newer than the XLS file it will skip the actual data transfer. As this is COM based this step is somewhat expensive even fro small amounts of data (there is some overhead involved in starting Excel).

Thursday, July 16, 2009

When running a large NLP in the GAMS/IDE with the new 23.1 version, I noticed that the log from KNITRO is buffered. This is not ideal: a log line should appear right away. Also the log from XA is completely disappeared since a few versions. In both cases the result is staring at a blank window when the job is running. The log is really the user-interface for a solver, so it would be good if a little bit more care is devoted to this.

The assignment to a is somewhat expensive: we can do this faster (see http://yetanothermathprogrammingconsultant.blogspot.com/2009/07/fill-parameters-with-random-sparse-data.html). The memory use seems correct: 407 MB for 12501427 elements, which is about 34 bytes per nonzero element. We use the 64-bit version of GAMS. With the 32-bit we will see smaller numbers: 304 MB or about 25 bytes per element. (That explains why 64-bit versions are not faster as 32-bit version: they have to drag more Mb around as pointers use up more memory: 8 bytes vs. 4 bytes). This is all ok.

The assignment of b1 is also ok. It is fast and only adds 500 elements, so it does not even register in the total memory usage.

The assignment of b2 is strange. It is understandable it is slower: the summation is using indices in the “wrong” order. The sudden spike in memory usage is less obvious. The additional hefty 806 MB (approximately twice what is needed to store parameter a itself) is used internally to create faster access to parameter a. This automatic reordering option can be turned off with OPTION SYS11=1;. If we use this option the profile report looks like:

This report has the memory profile we expected. Note that the assignment b2 is actually faster without reordering. The assignment b3 however shows why reordering though expensive (in memory and CPU time) can really pay off. Note: the “Other” in the table is related to the handling of the OPTION SYS11.

So if you see unexpected large spikes in memory usage, this may well be related to the automatic reordering of symbols that can take place when assignments are not in the “correct” index order.

The big advantage of an LP model is that we can extend this to a MIP for more complicated modeling situations. Although some solvers support MIQP models, the MIQP algorithms are often not as fast as MIP algorithms. (Microsoft Solver Foundation does have a QP algorithm, but as it is based on an interior point algorithm it is not really suited to be used inside a branch-and-bound framework, so there is no standard MIQP capability).

Everything else can stay the same. The results are shown in the screen dump below. As you can see the LP offers a reasonable approximation: the portfolios for different trade-offs between risk and return are quite similar. (Click on the picture to enlarge).

This exercise shows that beyond using the MSF Excel plug-in, the .NET API’s in combination with the modeling language OML are a powerful tool to write applications and that these solutions can be easily integrated in an Office Environment (Access, Excel).

In preparation for a presentation on different modeling languages for a project at a financial firm, I implemented a demo QP/LP portfolio model with an Excel front-end and a .NET/C# DLL.

MS Solver Foundation has a great tool to develop Excel based optimization models: the Excel plug-in. One of the restrictions however is that only one single model can be handled this way. Sometimes we need to solve a series of models or different variants of a model. In this post I will show this can be achieved by creating a .NET DLL that uses the MS Solver Foundation API’s on one hand, and interfaces with Excel through a COM interface on the other hand.

Note that there are other possible approaches such as VSTO to handle this. VSTO allows you to create .Net applications inside an Office environment. The Microsoft Solver Foundation Excel plug-in is an example of a VSTO application.

We used a portfolio problem as an example. Basically we have a multi-criteria problem:

where R are the mean adjusted returns. This separable formulation suggests an approximate LP formulation:

risk ≈ Σ |wt|/T

The OML model for the QP formulation looks like:

Afterwards we can return x[i], Return and calculate the Risk term as (Obj+Return)*NumT/Lambda. We wanted to keep the quadratic term in the objective for obvious reasons. Note that we need to substitute out {0} and {1} later so it holds values for the scalars Lambda and NumT. The substitution is done by String.Format, and we need to protect the other constructs like {t,T} by adding another pair of curly brackets: {{t,T}}.

There is a little bit of extra code to retrieve the solution. At the end we collect solution data as:

In Excel it is easy to display this graphically:

The architecture of this demo app is as follows (click to enlarge): The database not only stores the data, but also does some of the data manipulation: e.g. the mean adjusted returns are calculated in Access through queries.

Monday, July 6, 2009

In http://www.iging.com/optimization/ yet another input format for LP’s is proposed. I don’t think this particular format will work in practice: it stores all elements, even the zero ones. Large LP’s have typically a large number of nonzero elements. Lets take the medium sized LP indus89 (m=2,726, n= 6,570) from the GAMS model library:

Anything dense will never work for large models (“dense does not make sense”). I doubt this will become a standard format. Besides the inability to store large problems, it seems that ASCII formats are more popular, not in the least because they are platform independent and because you can actually look at the content.

Thursday, July 2, 2009

For testing I need sometimes very large, but very sparse data to play with. To generate random data, the easiest is something like:

a(i,j,k)$(uniform(0,1)<fracsparse) = uniform(-100,100);

However this will loop over every (i,j,k). If very sparse this can be become relatively expensive. It is faster in that case to loop over the expected number of nonzero elements, and place them randomly in the parameter. For this we use the lag operator for indices.

The two methods are not exactly the same: the second method can never generate more than nn nonzeroes, while the first method can. However for testing algorithms in GAMS the generated data is often “random enough” to give an indication about performance etc.

$ontext

fill a(i,j,k) and b(i,j,k) with random sparse data.Both sparsity pattern as the data should be random. We demonstrate two methods. The second method is efficient ifthe parameters are *very* sparse. $offtext

$set n 500 set i /i1*i%n%/; scalar n /%n%/; alias (i,j,k);

scalar fracsparse 'fraction of nonzero elements' /0.0001 /;

** method 1* this is slow because we loop over (i,j,k)* parameter a(i,j,k); a(i,j,k)$(uniform(0,1)<fracsparse) = uniform(-100,100);