```{r equi_shepardplot, fig.align = "center", echo = FALSE} par(mfrow=c(1,2), pty = "s") plot(h2a$delta, h2a$dist, xlab = "dissimilarities", ylab = "distances", col = "RED", xlim = c(0, 1.5), ylim = c(0, 1.5)) abline (0, 1, col = "BLUE", lwd = 2) plot(h2b$delta, h2b$dist, xlab = "dissimilarities", ylab = "distances", col = "RED", xlim = c(0, 1.5), ylim = c(0, 1.5)) abline (0, 1, col = "BLUE", lwd = 2) ```

To illustrate complementary slackness we give the values of the constraint functions $d_{ij}^2(X)-\delta_{ij}^2$, which are all non-positive, and the values of the Lagrange multipliers, which are all non-negative. We see that a constraint function that is equal to zero means a positive Lagrange multiplier, and a constraint function that is negative means a zero Lagrange multiplier. ```{r equi_complement} mprint (h2b$constraints, d = 6, f = "+") mprint (h2b$multipliers, d = 6, f = "+") ``` ##Dutch Political Parties 1967 As the next illustration we use data from @degruijter_67, with average dissimilarity judgments between Dutch political parties in 1967. The data are ```{r poldist-data, echo = FALSE} poldist <- structure( c( 5.63, 5.27, 4.6, 4.8, 7.54, 6.73, 7.18, 6.17, 6.72, 5.64, 6.22, 5.12, 4.59, 7.22, 5.47, 5.46, 4.97, 8.13, 7.55, 6.9, 4.67, 3.2, 7.84, 6.73, 7.28, 6.13, 7.8, 7.08, 6.96, 6.04, 4.08, 6.34, 7.42, 6.88, 6.36, 7.36 ), Labels = c("KVP", "PvdA", "VVD", "ARP", "CHU", "CPN", "PSP", "BP", "D66"), Size = 9L, call = quote(as.dist.default(m = polpar)), class = "dist", Diag = FALSE, Upper = FALSE ) print (poldist) poldata <- as.vector (poldist) polnames <- attr (poldist, "Labels") ``` ```{r below_comp, echo = FALSE, cache = FALSE} h1a <- smacof (poldata) h1b <- smacofBelow (poldata) ``` The optimal `smacof` solution in two dimensions needs `r h1a$itel` iterations to arrive at stress `r h1a$stress`. Approximation from below uses `r h1b$itel` iterations and finds stress `r h1b$stress`. At the solution there are `r length (which (abs (h1b$constraints) < 1e-6))` active constraints, but as the Shepard plots in figure `r figure_nums (name = "gruijter_regular_shepard", display = "n")` show, the active constraints now do not correspond with the largest dissimilarities, and the two configurations in figure `r figure_nums (name = "gruijter_regular_conf", display = "n")` are quite different. Although that could, of course, also be because they correspond with two different local minima, which mostly differ in the location of the D66 party (a neoliberal party, full of civilized protest).

```{r smacof_confplot, fig.align = "center", echo = FALSE} par(mfrow=c(1,2), pty = "s") plot(h1a$x, type = "n", xlab = "dim 1", ylab = "dim 2")#, xlim = c(-.15, .15), ylim = c(-.15, .15)) text(h1a$x, polnames, col = "RED") plot(h1b$x, type = "n", xlab = "dim 1", ylab = "dim 2")#, xlim = c(-.15, .15), ylim = c(-.15, .15)) text(h1b$x, polnames, col = "RED") ```

```{r smacof_shepardplot, fig.align = "center", echo = FALSE} par(mfrow=c(1,2), pty = "s") plot(h1a$delta, h1a$dist, xlab = "dissimilarities", ylab = "distances", col = "RED")# xlim = c(0, .25), ylim = c(0, .25)) abline (0, 1, col = "BLUE", lwd = 2) plot(h1b$delta, h1b$dist, xlab = "dissimilarities", ylab = "distances", col = "RED")# xlim = c(0, .25), ylim = c(0, .25)) abline (0, 1, col = "BLUE", lwd = 2) ```

In the next analysis we require that all distances are less than or equal to 8.13, the maximum dissimilarity. In order words, the maximum distances must be less than or equal to the maximum dissimilarity. ```{r gruijter_maxbound_comp, echo = FALSE} h1c <- smacofBelow (poldata, bnd = rep(max(poldata), length (poldata))) ``` The optimal solution under these constraints in two dimensions needs `r h1c$itel` iterations to arrive at stress `r h1c$stress`. At the solution there are `r length (which (abs (h1c$constraints) < 1e-6))` active constraints, i.e. six distances equal to 8.13.

```{r gruijter_maxbound_plot, fig.align = "center", echo = FALSE} par(mfrow=c(1,2), pty = "s") plot (h1c$x, type = "n", xlab = "dim 1", ylab = "dim 2") text(h1c$x, polnames, col = "RED") plot (h1c$delta, h1c$dist, xlab = "dissimilarities", ylab = "distances", col = "RED") abline (0, 1, col = "BLUE", lwd = 2) ```

And finally an analysis in which we require that distances between ARP, CHU, KVP (the Christian Democrat parties) are less than 0.1, and the distances between PvdA, PSP, CPN (the leftist parties) are also less than 0.1. This will effectively collapse them in the configuration plot. ```{r gruijter_infbound_comp, echo = FALSE} bnd <- rep (Inf, length (poldata)) bnd[c(3, 4, 22, 12, 13, 31)] <- 4 h1d <- smacofBelow (poldata, bnd = bnd, verbose = FALSE) ``` The optimal solution under these constraints in two dimensions needs `r h1d$itel` iterations to arrive at stress `r h1d$stress`. At the solution all `r length (which (abs (h1d$constraints) < 1e-6))` constraints are active.

```{r gruijter_infbound_plot, fig.align = "center", echo = FALSE} par(mfrow=c(1,2), pty = "s") plot (h1d$x, type = "n", xlab = "dim 1", ylab = "dim 2", xlim = c(-5, 5), ylim = c(-5, 5)) text(h1d$x, polnames, col = "RED") plot (h1d$delta, h1d$dist, xlab = "dissimilarities", ylab = "distances", col = "RED") abline (0, 1, col = "BLUE", lwd = 2) ```

##Ekman Color Data ```{r ekman, echo = FALSE} data (ekman, package = "smacof") nms <- as.character (attr (ekman, "Labels")) ekman <- as.vector (1 - ekman) ``` ```{r ekman_run, echo = FALSE, cache = FALSE} h3a <-smacof (ekman) h3b <- smacofBelow (ekman) ``` The final example are the color data from @ekman_54. The optimal `smacof` solution in two dimensions needs `r h3a$itel` iterations to arrive at stress `r h3a$stress`. Approximation from below uses `r h3b$itel` iterations and finds stress `r h3b$stress`. At the solution there are `r length (which (abs (h3b$constraints) < 1e-6))` active constraints. Again the bounded solution is a shrunken version of the regular `smacof` solution, and figure `r figure_nums (name = "ekman_regular_shepard", display = "n")` shows that the active constraints correspond with the largest dissimilarities and distances.

```{r ekman_confplot, fig.align = "center", echo = FALSE} par(mfrow=c(1,2), pty = "s") plot(h3a$x, type = "n", xlab = "dim 1", ylab = "dim 2")# xlim = c(-.08, .08), ylim = c(-.08, .08)) text(h3a$x, nms, col = "RED") plot(h3b$x, type = "n", xlab = "dim 1", ylab = "dim 2")# xlim = c(-.08, .08), ylim = c(-.08, .08)) text(h3b$x, nms, col = "RED") ```

```{r ekman_shepardplot, fig.align = "center", echo = FALSE} par(mfrow=c(1,2), pty = "s") plot(h3a$delta, h3a$dist, xlab = "dissimilarities", ylab = "distances", col = "RED")# xlim = c(0, .15), ylim = c(0, .15)) abline (0, 1, col = "BLUE", lwd = 2) plot(h3b$delta, h3b$dist, xlab = "dissimilarities", ylab = "distances", col = "RED")# xlim = c(0, .15), ylim = c(0, .15)) abline (0, 1, col = "BLUE", lwd = 2) ```

Figure `r figure_nums("ekman_distances", display = "n")` shows that the two configurations in figure `r figure_nums ("ekman_regular_conf", display = "n")` may look to be proportional, they really are not. There are many deviations from proportionality if we compare the fitted distances in the regular and bounded `smacof` solutions.

```{r ekman_distances, fig.align = "center", echo = FALSE} plot(h3a$dist, h3b$dist, col = "RED", xlab = "regular distances", ylab = "bounded distances") ```

#Discussion If $p=n-1$ our multidimensional scaling problem is full-dimensional scaling or FDS (@deleeuw_groenen_mair_E_16e). In that case stress can be parametrized by using the scalar products $C$ and we must minimize $$ \mathbf{stress}(C):=\sum_{i=1}^n\sum_{j=1}^n w_{ij}(\delta_{ij}-\sqrt{\mathbf{tr}\ A_{ij}C})^2 $$ over $C\gtrsim 0$ and $\mathbf{tr}\ A_{ij}C\leq\alpha_{ij}$. This means that FDS with upper bounds means minimizing a convex function over a convex set, in fact solving a semidefinite programming problem. It also means that the concept of Gower rank starts playing a role (@deleeuw_E_16k). MDS from below can also be used with different loss functions. Minimizing $$ \sum_{i=1}^n\sum_{j=1}^nw_{ij}|\delta_{ij}-d_{ij}(X)| $$ over $d_{ij}(X)\leq\delta_{ij}$ is the same as maximizing the weighted sum of the distances, a convex function. Thus we can use tangential majorization (@deleeuw_B_16b), which means the majorization function is linear and has no quadratic term. Alternatively, the problem of minimizing $$ \sum_{i=1}^n\sum_{j=1}^nw_{ij}|\delta_{ij}^2-d_{ij}^2(X)| $$ over $d_{ij}(X)\leq\delta_{ij}$ means maximizing a convex quadratic form. Although the objective function is quadratic, we still need tangential majorization to get a convex objective to minimize. We do not address the problem in this paper of minimizing stress under constraints of the form $\beta_{ij}\leq d_{ij}(X)\leq\alpha_{ij}$. This is a fundamentally more complicated problem, because the lower bounds make the solution set non-convex. In fact, using lower bounds may make the solution set empty, because there may not be $p$-dimensional Euclidean distances satisfying the constraints. It will be of some interest to apply our algorithm to unfolding examples, where we could have upper bound constraints on the row distances, the column distances, or the row-column distances. Of course in the unfolding case lower bound constraints are at least as interesting, because they can in principle be used to avoid the trivial solutions that are so common in unfolding, but we have argued that computationally lower bound constraints are far more difficult to handle. #Appendix: Code ##below.R ```{r file_auxilary, code = readLines("below.R")} ``` ##auxilary.R ```{r file_auxilary, code = readLines("auxilary.R")} ``` ##mdsUtils.R ```{r file_auxilary, code = readLines("mdsUtils.R")} ``` ##smacof.R ```{r file_auxilary, code = readLines("smacof.R")} ``` ##smacofBelow.R ```{r file_auxilary, code = readLines("smacofBelow.R")} ``` #References