New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Scaling and universality in continuous length combinatorial optimization

Communicated by Peter J. Bickel, University of California, Berkeley, CA, August 13, 2003 (received for review June 1, 2003)
Abstract
We consider combinatorial optimization problems defined over random ensembles and study how solution cost increases when the optimal solution undergoes a small perturbation δ. For the minimum spanning tree, the increase in cost scales as δ^{2}. For the minimum matching and traveling salesman problems in dimension d ≥ 2, the increase scales as δ^{3}; this is observed in Monte Carlo simulations in d = 2, 3, 4 and in theoretical analysis of a meanfield model. We speculate that the scaling exponent could serve to classify combinatorial optimization problems of this general kind into a small number of distinct categories, similar to universality classes in statistical physics.
The interface of statistical physics, algorithmic theory, and mathematical probability is an active research field, containing diverse topics such as mixing times of Glaubertype dynamics (ref. 1 and many others), reconstruction of broadcast information (2), and probabilistic analysis of paradigm computational problems such as kSAT (3–5). In this article we introduce an article topic whose motivation is simpler than those.
Freshman calculus tells us that, for a smooth function attaining its minimum at x*, for x near x* the relation between δ = x – x* and ε = F(x) – F(x*) is ε ≈ (1/2)F″(x*)δ^{2}. If instead we consider a function on ddimensional space, sophomore calculus tells us that similarly for appropriate c. So in a sense the scaling exponent 2 is naturally associated with “smooth” or “regular” optimization problems.
Now consider a graphbased combinatorial optimization problem, such as the traveling salesman problem (TSP): each feasible solution has n constituents (edges) and associated continuous costs (lengths), the sum of which gives the overall solution cost. Compare an arbitrary feasible solution x with the optimal (minimal) solution x*, unique, which for generic lengths is unique, by the two quantities where s(n) expresses the rate at which the optimal cost scales in n. Define ε_{n}(δ) to be the minimum value of ε_{n}(x) over all feasible solutions x for which δ_{n}(x) ≥ δ. Although the function ε_{n}(δ) will depend on n and the problem instance, we anticipate that for typical instances drawn from a suitable probability model it will converge in the n → ∞ limit to some deterministic function ε(δ).
The universality paradigm from statistical physics suggests there may be a scaling exponent a such that and that the exponent should be robust under model details. In statistical physics, universality classes are typically defined by critical exponents that characterize the behavior of measurable quantities both near and at a phase transition. Although a is not a critical exponent here, and there is no phase transition, we suggest that it could play a similar role, categorizing combinatorial optimization problems into a small set of classes. If our analogy with freshman calculus is apposite, we expect that the simplest problems should have scaling exponent 2.
This approach may seem obvious in retrospect and fits within a longstanding tradition in the physical sciences (see Discussion). However, it has never been proposed or studied explicitly. In this article we report on three aspects of our program. For the minimum spanning tree (MST), a classic “algorithmically easy” problem solvable to optimality by greedy methods, we confirm that the scaling exponent is indeed 2. We then turn to two harder problems: minimum matching (MM) and the TSP. Under a meanfield model, our mathematical analysis methods combined with numerics show that the scaling exponent is 3 for both MM and TSP, independent of the pseudodimension defined below. For the Euclidean model the exponent is 2 in the (essentially trivial) 1D case, while Monte Carlo simulations suggest it is 3 in higher dimensions.
Models
In the Euclidean model we take n random points in a ddimensional cube whose volume scales as n. Interpoint lengths are Euclidean distances. To reduce finitesize effects, we take the space to have periodic (toroidal) boundary conditions when calculating the distances.
In the meanfield or random link model we imagine n random points in some abstract space such that the (^{n}_{2}) vertex pair lengths are i.i.d. random variables distributed as n^{l/}^{d}l, with probability density p(l) ∼ l^{d}^{–1} for small l. Here 0 < d < ∞ is the pseudodimension parameter and the distribution of small single interpoint lengths mimics that in the Euclidean model of corresponding dimension d, up to a proportionality constant. Both models are set up so that nearestneighbor distances are of order 1 and the scaling of overall cost in the optimization problems is s(n) = n.
A Simple Case: The MST
For the MST, given any reasonable model of interpoint lengths including the two models above, we expect a scaling exponent of 2. We will provide a rigorous account elsewhere, but the underlying idea is simple. The classical greedy algorithm gives the following explicit inclusion criterion for whether an edge e = (v_{1}, v_{2}) of a graph belongs in the MST. Consider the subgraph containing edges between any two vertices within length t of each other. Let perc(e) ≤ length(e) be the smallest t that keeps v_{1} and v_{2} within the same connected component. It is not difficult to see that e ∈ MST if and only if length(e) = perc(e).
Given a probability model for n random points and their interpoint lengths, define a measure μ_{n}(x) on x ∈ (0, ∞) in terms of the expectation For any reasonable model we expect an n → ∞ limit measure μ(x), with a density ν(x) = dμ/dx having a nonzero limit ν(0^{+}).
Now modify the MST by adding an edge e with length(e) – perc(e) = b, for some small b, to create a cycle: then delete the longest edge e′ ≠ e of that cycle, which necessarily has length(e′) = perc(e). This gives a spanning tree containing exactly one edge not in the MST and having length greater by b. Repeat this procedure with every edge e for which 0 < length(e) – perc(e) < β, for some small β. The number of such edges is nμ(β) ≈ nν(0^{+})β to first order in β, and as there is negligible overlap between cycles, each of the new edges will increase the tree length by ∼β/2 on average. So This construction must yield essentially the minimum value of ε for given δ, so the scaling exponent is 2.
Poisson Weighted Infinite Tree (PWIT)
We now consider the MM and TSP. In MM, we ask for the minimum total length L_{n} of n/2 edges matching n random points and study the normalized limit expectation lim_{n}_{→∞} (2/n)E[L_{n}]. Taking the meanfield model with d = 1 for simplicity, the limit value π^{2}/6 was obtained in ref. 6 by using the replica method from statistical physics. We work in the framework of ref. 7, which rederives this limit rigorously by doing calculations within an n = ∞ limit structure, the PWIT.
Briefly, the PWIT is an infinite degree rooted tree in which the edge weights (lengths) at each vertex are distributed as the successive points 0 < ξ_{1} < ξ_{2} <... of a Poisson process with a mean number x^{d} of points in [0, x], i.e., a process with rate increasing as d x^{d}^{–1}. In this way, the PWIT corresponds to the meanfield model at a given d (see ref. 8 for review).
Consider a matching on an instance of a rooted PWIT, as well as a matching on the same instance but with the root removed, as shown in Fig. 1. Introduce the variable Both lengths are infinite, so this is interpreted as a limit of finite differences. If X_{i} is the analogous quantity for the ith constituent subtree of the rootless PWIT instance and ξ_{i} the length of the root's ith edge, these variables satisfy the recursion [1]
Now take the {ξ_{i}} to be the Poissondistributed edge lengths and the {X_{i}} to be independent random variables from the same random process that produces X. Eq. 1 is then a distributional equation for X and can be shown (7) for d = 1 to have as its unique solution the logistic distribution [2]
The PWIT structure further leads to the following inclusion criterion. Consider an edge of length x in the tree, and the two subtrees formed by deleting that edge. The memoryless nature of the Poisson process allows us to consider each of these subtrees as independent copies of a PWIT, with their roots at the vertices of the deleted edge. It may be seen that including the edge in the optimal matching incurs a cost of x – X_{1} – X_{2}, where X_{1} and X_{2} are the X variables as defined above, but for the two subtrees. Thus, an edge of length x is present in the minimal matching if and only if [3]
The probability density function for edge lengths in the MM is then Here X_{1} and X_{2} are independent random variables distributed according to Eq. 2, from which the mean edge length can be calculated:
MeanField MM and TSP
The previous section summarized analysis from ref. 7; now we continue with new analysis. To study scaling exponents, we introduce a parameter λ > 0 that plays the role of a Lagrange multiplier. Penalize edges used in the optimal matching by adding λ to their length. Let us study optimal solutions to the MM problem on this new penalized instance. Precisely, on a realization of the PWIT, define Y and Z as where Y and Z differ in the definition of the edge lengths of the new tree: for Y, the edges penalized are those used by the original rooted optimal matching: for Z, they are those used by the original rootless optimal matching.
For the penalized problem the recursion Eq. 1 for X is supplemented by the following recursions for (X, Y, Z) jointly. Let i* be the value of i that minimizes ξ_{i} – X_{i}. Then where, as before, the {Y_{i}} and {Z_{i}} are independent random variables from the same random process producing Y and Z.
Moreover, we get an inclusion criterion, analogous to Eq. 3:an edge of length x is included if and only if In terms of the expected unique joint distribution for (X, Y, Z), the quantities δ and ε that compare the penalized solution (as a nonoptimal solution of the original problem) with the original optimal solution are and By the theory of Lagrange multipliers these functions ε(λ), δ(λ) determine ε(δ). We do not have explicit analytic expressions analogous to Eq. 2 for the joint distribution of (X, Y, Z) in terms of λ. However, we can use routine bootstrap Monte Carlo simulations to simulate the distribution and thence estimate the functions δ(λ) and ε(λ) numerically. And as indicated in refs. 7, 9, and 10 the meanfield MM and the meanfield TSP can be studied by using similar techniques; the TSP analysis is just a minor variation of the MM analysis. For instance, recursion Eq. 1 becomes where min^{[2]} denotes second minimum.
Table 1 reports numerical results showing good agreement with ε ∞ δ^{3} in both problems for d = 1. These numerics are compatible with independent MM results obtained recently (11), as well as with our direct simulations on meanfield TSP instances at n = 512. The same exponent 3 arises for other d.
Euclidean MM and TSP
We consider the d = 1 case where the scaling exponent can be found exactly and give numerical results for other cases. We restrict the discussion to the Euclidean TSP, although as for the meanfield model, MM is phenomenologically similar.
Take the Euclidean TSP in d = 1, with periodic boundary conditions. The optimal tour here is trivial (with high probability a straight line of length n) but nevertheless instructive to analyze. As before, add a penalty term λ to each edge used in the tour and consider how the optimal tour changes in this new penalized instance. When λ is small, changes to the tour will consist of “2changes” shown in Fig. 2 and will occur when an original edge length is < λ. A simple nearestneighbor distance argument gives the distribution of edge lengths in the original tour as p(l) ∼ e^{–}^{l}. Since two edges are modified in each 2change, The scaling exponent of 2 is not surprising, as the 1D TSP behaves very similarly to the 1D MST on penalized instances. Furthermore, it is consistent with the intuition that the “easiest” problems scale in this way. A similar argument applies to MM, and in both cases a more rigorous analysis yields bounds on δ and ε that confirm the exponent.
For d > 1, numerical results are shown in Fig. 3. These have been obtained by finding exact solutions to randomly generated n = 512 Euclidean instances in d = 2, 3, 4, using the concorde TSP solver available at www.math.princeton.edu/tsp/concorde.html. For each instance, the optimum was obtained on the original instance as well as on the instance penalized with a range of λ values. For each λ value, δ(λ) and ε(λ) were averaged over the sample of instances. The resulting numerics are closely consistent with a scaling exponent of 3 (in spite of suffering from some finitesize effects at smaller δ), suggesting that the meanfield picture gives the correct exponent in all but the trivial 1D case. In the language of critical exponents, this would correspond to an “upper critical dimension” of 2.
Discussion
The goal of our scaling study has been to address a new kind of problem in the theory of algorithms, using concepts from statistical physics. Traditionally, work on the TSP within the theory of algorithms (12) has emphasized algorithmic performance, rather than the kinds of questions we ask here. Rigorous study of the Euclidean TSP model within mathematical probability (13) has yielded a surprising amount of qualitative information: existence of an n → ∞ limit constant giving the mean edgelength in the optimal TSP tour (14), and large deviation bounds for the probability that the total tour length differs substantially from its mean (15). However, calculation of explicit constants in dimensions d ≥ 2 seems beyond the reach of analytic techniques. For the meanfield bipartite MM problem, impressive recent work (26, 27) has proven an exact formula giving the expectation of the finiten minimum total matching length, though such exact methods seem unlikely to be widely feasible.
On the other hand, there has been significant progress over the past 20 years in the use of statistical physics techniques on combinatorial optimization problems in general. Finding optimal solutions to these problems is a direct analog to determining ground states in statistical physics models of disordered systems (16). This observation has motivated the development of such approaches as simulated annealing (17), the replica method (18), and the cavity method (4). Condensed matter physics, particularly models arising in spin glass theory, has provided a powerful means to study algorithmic problems: at the same time, algorithmic results have implications for the associated physical models. It is instructive to consider our work in that context.
Researchers in the physical sciences have long been interested in the lowtemperature thermodynamics (18, 19) of disordered systems, investigating properties of nearoptimal states in spin glass models. Our procedure for studying nearoptimal solutions by way of a penalty parameter is similar to a method, known as εcoupling (20–22), used for calculating lowenergy excitations in spin glasses. Making use of this method, physicists have obtained quantities closely analogous to our scaling exponents for models of RNA folding (22). Furthermore, in the last year independent work (11) has explored εcoupling on MM, numerically identifying a different but related scaling exponent.
For the TSP, analytical and numerical studies were performed >15 years ago (23, 24) on the thermodynamics of the model, with overlap quantities calculated for nearoptimal solutions. The results have suggested that at low temperature T, the cost excess ε scales as T^{2} while the average fraction of differing edges between solutions (1 – q, with q being the “overlap fraction”) scales as T. This leads to ε ∼ (1 – q)^{2}, in apparent contradiction with our exponent of 3. However, at low temperatures, q represents overlaps between typical nearoptimal solutions, whereas our δ measures overlaps between a nearoptimal solution and the optimum. The different definitions of these two quantities could account for the discrepancy in scaling exponent: it is not surprising that 1 – q grows faster than δ as one considers solutions of increasing cost. At the same time, a possible implication of these results is that at low temperature, δ ∼ T^{2/3}. We are not aware of any direct theoretical arguments to explain this and consider it an intriguing open question.
It is also important to note that the underlying property δ → 0 as ε → 0 cannot always be taken for granted. This property is called asymptotic essential uniqueness (AEU) (7). AEU requires, among other things, that the optimum itself be unique. In principle, even if it is not, one could still analyze nearoptimal scaling by considering sufficiently local perturbations from a given optimum. It is natural to expect the resulting exponent to be independent of the specific optimum chosen. However, this may not be true in the event of what statistical physicists call replica symmetry breaking (RSB) (18, 19). AEU is a special case of replica symmetry, so while RSB implies the absence of AEU, the absence of AEU does not necessarily imply RSB. A current debate in condensed matter literature concerns whether or not lowtemperature spin glasses display RSB (20, 21, 25). It is generally believed that RSB is incompatible with unique nonzero values of various scaling exponents. Thus, the correct approach to analyzing nearoptimal scaling in such problems remains another open question.
One final example may serve to illustrate the diversity of possible applications for our type of scaling analysis, as well as an instance where the absence of AEU is surmountable. In oriented percolation on the 2D lattice, there are independent random traversal times on each oriented (up or right) edge. The percolation time T_{n} is the minimum, over all (^{2}^{n} n) paths from (0, 0) to (n, n), of the time to traverse the path. So (2n)^{–1}E[T_{n}] → t*, a time constant. It is elementary that there will be nearoptimal paths, with lengths T′_{n} such that n^{–1} (E[T′_{n}] – E[T_{n}]) → 0 and which are almost disjoint from the optimal path. So our ε(δ) analysis applied to paths will not be useful: even with a unique optimum, AEU will not hold. But we can rephrase the problem in terms of flows. A flow on the n × n oriented torus assigns to each edge a flow of size ∈ [0, 1], such that at each vertex, inflow equals outflow. Let t(δ) be the minimum, over all flows with mean flowperedge = δ, of the flowweighted average edge traversal time. In the n → ∞ limit, one can show that as δ → 0, t(δ) → t* where t* is the same limiting constant as before. We therefore expect a scaling t(δ) – t* ∼ δ^{a}. Mean field analysis gives the scaling exponent a = 2, and Monte Carlo study of the d = 2 case is in progress.
Conclusions
We have studied the scaling of the relative cost difference ε between optimal and nearoptimal solutions to combinatorial optimization problems, as a function of the solution's relative distance δ from optimality. This kind of scaling study, although well accepted in theoretical physics, is new to combinatorial optimization. For the MST, we have found ε ∼ δ^{2}. For the MM and TSP, in the 1D Euclidean case ε ∼ δ^{2} as well, while in both the meanfield model and higher Euclidean dimensions ε ∼ δ^{3}.
The scaling exponent may categorize combinatorial optimization problems into a small number of classes. The fact that MST is solvable by a simple greedy algorithm, and that the 1D case of the MM and TSP is essentially trivial, suggests that a scaling exponent of 2 characterizes problems of very low complexity. The exponent of 3 characterizes problems that are algorithmically more difficult. Of course, this is a different kind of classification from traditional notions of computational complexity: MM is solvable to optimality in O(n^{3}) time whereas the TSP is in the NPhard class. Rather, these exponent classes are reminiscent of universality classes in statistical physics, which unite diverse physical systems exhibiting identical behavior near phase transitions.
A key question in the study of critical phenomena is whether meanfield models correctly describe phase transition behavior in the geometric models they approximate. The TSP and MM do not involve critical behavior, but the fact that meanfield and geometric scaling exponents coincide for d ≥ 2 is significant. It provides evidence that in a combinatorial setting, the meanfield approach can give a valuable and accurate description of the structure of nearoptimal solutions.
Acknowledgments
We thank Mike Steele for helpful discussions. D.A.'s research is supported by National Science Foundation Grant DMS0203062. A.P. acknowledges support from Department of Energy Grant LDRD/ER 20030137 and the kind hospitality of the Institute for Pure and Applied Mathematics at the University of California (Los Angeles), where much of this work was conducted.
Footnotes

↵§ To whom correspondence should be addressed. Email: percus{at}lanl.gov.

Abbreviations: TSP, traveling salesman problem; MST, minimum spanning tree; MM, minimum matching; PWIT, Poisson weighted infinite tree; AEU, asymptotic essential uniqueness; RSB, replica symmetry breaking.
 Received June 1, 2003.
 Copyright © 2003, The National Academy of Sciences
References
 ↵
 ↵
 ↵
Coppersmith, D., Gamarkin, D., Hajiaghayi, M. T. & Sorkin, G. B. (2003) in Proceedings of the 14th Annual ACMSIAM Symposium on Discrete Algorithms, ed. FarachColton, M. (Association for Computing Machinery, New York), pp. 364–373.
 ↵
Mézard, M., Parisi, G. & Zecchina, R. (2002) Science 297, 812–815.pmid:12089451
 ↵
 ↵
 ↵
 ↵
Aldous, D. J. & Steele, J. M. (2004) in Probability on Discrete Structures, ed. Kesten, H. (Springer, Berlin), in press.
 ↵
Krauth, W. & Mézard, M. (1987) Europhys. Lett. 8, 213–218.
 ↵
Cerf, N. J., Boutet de Monvel, J., Bohigas, O., Martin, O. C. & Percus, A. G. (1997) J. Phys. I 7, 117–136.
 ↵
Ratiéville, M. (2002) Ph.D. thesis (Université Pierre et Marie Curie, Paris and Università degli Studi di Roma “La Sapienza,” Rome).
 ↵
Gutin, G. & Punen, A. P., eds. (2002) The Traveling Salesman Problem and its Variations (Kluwer, Dordrecht, The Netherlands).
 ↵
Steele, J. M. (1997) Probability Theory and Combinatorial Optimization (Society for Industrial and Applied Mathematics, Philadelphia).
 ↵
 ↵
 ↵
 ↵
Kirkpatrick, S., Gellatt, C. D. & Vecchi, M. P. (1983) Science 220, 671–680.
 ↵
Mézard, M., Parisi, G. & Virasoro, M. A. (1987) Spin Glass Theory and Beyond (World Scientific, Singapore).
 ↵
Dotsenko, V. (2001) Introduction to the Replica Theory of Disordered Statistical Systems (Cambridge Univ. Press, Cambridge, U.K.).
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
Nair, C., Prabhakar, B. & Sharma, M. (2003) in Proceedings of the 44th Annual IEEE Symposium on Foundations of Computer Science (Institute of Electrical and Electronics Engineers, Los Alamitos, CA), in press.
 ↵
Linusson, S. & Wästlund, J. (2003) Technical report LiTHMATR200303 (Linköping University, Linköping, Sweden).
Citation Manager Formats
Sign up for Article Alerts
Jump to section
You May Also be Interested in
More Articles of This Classification
Physical Sciences
Related Content
 No related articles found.
Cited by...
 No citing articles found.