Review: Multidimensional Scaling Using Majorization: SMACOF in R¶

Citation

de Leeuw, J., & Mair, P. (2009). Multidimensional scaling using majorization: SMACOF in R. Journal of Statistical Software, 31(3), 1–30.
DOI

Abstract¶

In this paper we present the methodology of multidimensional scaling problems (MDS) solved by means of the majorization algorithm. The objective function to be minimized is known as stress and functions which majorize stress are elaborated. This strategy to solve MDS problems is called SMACOF and it is implemented in an R package of the same name which is presented in this article. We extend the basic SMACOF theory in terms of configuration constraints, three-way data, unfolding models, and projection of the resulting configurations onto spheres and other quadratic surfaces. Various examples are presented to show the possibilities of the SMACOF approach offered by the corresponding package.

de Leeuw and Mair (2009) present the majorization derivation of SMACOF, an acronym for Scaling by MAjorizing a COmplicated Function. The original derivation is due to de Leeuw (1977), which predates digital object identifiers; this paper is the citable exposition. Sections 2 and 3.1 contain the material Hifuku uses. The remainder describes the R package and extensions to constrained, three-way and quadratic-surface scaling, which are not used.

Majorization¶

Section 2 states the principle. To minimize an awkward function \( f(x) \), construct a simpler surrogate \( g(x, y) \) that majorizes it, meaning \( g(x, y) \ge f(x) \) for all \( x \), and that touches it at the supporting point, \( f(y) = g(y, y) \). Minimizing the surrogate then yields the sandwich inequality

\[f(x^\star) \le g(x^\star, y) \le g(y, y) = f(y)\]

so each iteration cannot increase \( f \). The procedure is to choose a starting value, minimize the surrogate, and stop when the decrease falls below a tolerance. The paper notes that majorization is a prescription for constructing algorithms rather than an algorithm itself, and is related to EM and to MM methods more generally.

Stress, and the role of the weights¶

Section 3.1 defines stress for an \( n \times p \) configuration \( X \) against a dissimilarity matrix \( \Delta \):

\[\sigma(X) = \sum_{i<j} w_{ij}\,\bigl(\delta_{ij} - d_{ij}(X)\bigr)^2\]

The weight matrix \( W \) is symmetric, non-negative and hollow. Two of its stated properties matter here.

The first is its documented use: the paper notes that \( W \) can impose a missing-value structure, with \( w_{ij} = 1 \) where \( \delta_{ij} \) is known and \( w_{ij} = 0 \) where it is missing. This is exactly how Hifuku uses it in chart selection. A candidate set of anchor triangles retains some pairs of anchors and omits others, and the retained set is expressed as a weight matrix rather than as a reduced problem.

The second is the normalization, equation 5:

\[\sum_{i<j} w_{ij}\,\delta_{ij}^2 = \frac{n(n-1)}{2}\]

Together these explain a result that is easy to misread. Because stress sums only over pairs with positive weight, removing edges removes terms, so raw stress falls as a candidate edge set is thinned. The decrease is an artifact of the edge count, not evidence of a better embedding. Hifuku's reconstruction_error therefore scores a re-embedded configuration against every pairwise distance, including those the candidate set dropped, and normalizes so that sets of different size are comparable.

The Guttman transform¶

The derivation decomposes stress into a constant, a convex quadratic, and a concave term. With \( A_{ij} = (e_i - e_j)(e_i - e_j)' \), the paper defines

\[V = \sum_{i<j} w_{ij} A_{ij}, \qquad B(X) = \sum_{i<j} w_{ij}\,s_{ij}(X)\,A_{ij}\]

where \( s_{ij}(X) = \delta_{ij} / d_{ij}(X) \) when \( d_{ij}(X) > 0 \) and zero otherwise, giving

\[\sigma(X) = \frac{n(n-1)}{2} + \operatorname{tr} X'VX - 2\operatorname{tr} X'B(X)X\]

The Cauchy-Schwarz inequality supplies \( \rho(X) \ge \operatorname{tr} X'B(Y)Y \), which linearizes the concave term and produces a quadratic majorizer. Setting its derivative to zero and solving with the Moore-Penrose inverse of \( V \) gives the update

\[\bar{X} = V^{+} B(Y)\,Y\]

known as the Guttman transform. Where all weights are one, \( V \) reduces to \( 2n(I - n^{-1}\mathbf{11}') \) and the update simplifies to \( \bar{X} = n^{-1} B(Y) Y \).

Majorization guarantees a non-increasing sequence of stress values, at a linear convergence rate.

Relevance to Hifuku¶

Hifuku uses SMACOF in one place: scoring a candidate set of anchor triangles during chart selection. The weighted formulation is what makes the question expressible at all, since a candidate set is naturally a pattern of retained and omitted anchor pairs.

The paper also documents a limitation that shaped the implementation. It notes that local minima are more likely in low-dimensional solutions, while full-dimensional scaling has no local minimum problem. Two-dimensional embeddings are the worst case. Hifuku therefore uses classical multidimensional scaling, which is deterministic and has a closed form, wherever a fixed two-dimensional configuration is required, and reserves SMACOF for the weighted problem, where no closed form exists.

Section	Content
Section 1	Introduction and MDS taxonomies
Section 2	Majorization: the sandwich inequality and the iterative scheme
Section 3.1	Stress, the weight matrix, the stress decomposition, and the Guttman transform
Section 3.2	SMACOF with restrictions on the configuration
Section 4	Nonmetric SMACOF variants
Section 5	Extended SMACOF: projection onto quadratic surfaces
Section 6	The R package smacof, with worked examples
Section 7	Discussion