Plotting scatterplot with non-dominated front in R

I have a dataset with two objectives (runtime and solution quality). Is there some way to draw a scatterplot in R that also draws a non-dominated front of the two objectives?

• Non-dominating or non-dominated? Oct 10 '11 at 21:00
• @Henry: "non-dominated front". Edited. Oct 10 '11 at 21:30
• @Chris can you please be more specific? It is not clear what you are looking for. Oct 11 '11 at 13:57
• @David: I'm sorry. Well, I want to draw a scatterplot in R where the non-dominated points are highlighted in some way. Oct 14 '11 at 16:08
• It's very easy to highlight a subset of points in R, by using points() on a subset of the data with different settings for point character or colour from the original plot you are adding the points to. So if you have a way of working out which are the "non-dominated points" (and sorry, I have no idea what this is) it is easy to highlight them. Jan 18 '12 at 0:14

Presumably "runtime" is better when it is low and "solution quality" is better when it is high. To make these variables more consistent in how they represent the attributes, first re-express them so that high numbers always correspond to better values. For this question, the reciprocal of runtime (which could be interpreted as number of runs per second) would be a good choice, but the negative of runtime would work fine too.

With this convention established, we may plot each record in the dataset by associating the re-expressed "runtime" value with the x coordinate and the "solution quality" value with the y coordinate, creating a scatterplot. One record $(x,y)$ dominates another $(x',y')$ when simultaneously $x\ge x'$, $y\ge y'$, and $(x,y) \ne (x',y')$. No rational actor (using solely information about "runtime" and "solution quality") would choose the dominated option $(x',y')$ when $(x,y)$--which is at least as good on both attributes and better on at least one--is available. Therefore such an actor would focus their interest on the non-dominated options and would be free to neglect all others. All points dominated by the solid red point are shown in gray. All possible points that it can dominate constitute the (closed) rectangle shown in pink. Because evidently no points in the scatterplot dominate the solid point, it must lie on the non-dominated frontier.

The set of non-dominated options forms the vertices of part of the "quasi-convex hull" or "non-dominated frontier" of the scatterplot. This is an analog of the better-known convex hull. Indeed, the quasi-convex hull will be invariant under all increasing monotonic re-expressions of the variables, $f(x)$ and $g(y)$, and one can always find such re-expressions for which the quasi-convex hull actually becomes the convex hull. This connection, and the fact that optimal construction of the convex hull of $n$ points requires $O(n\log(n))$ computation, suggests that we seek an algorithm for the quasi-convex hull that performs at least as well. Such algorithms exist. Here is one, given in R as requested, which is written to port easily to other computing platforms:

hull.qc <- function(xy) {
i <- order(xy[, 1], xy[, 2], decreasing=TRUE)
frontier <- rep(0, length(i))
k <- 0; y <- -Inf
for (j in i) {
if (xy[j, 2] > y) {
frontier[k <- k+1] <- j
y <- xy[j, 2]
}
}
return(frontier[1:k])
}

It performs a decreasing lexicographic sort of the coordinates (time: $O(n\log(n))$) and then scans through the first coordinate (time: $O(n)$), putting into effect a line scan algorithm. The points where a new larger value of the second coordinate is found are recorded. Their indexes within the input array xy are returned.

We can prove this algorithm is correct by induction. The initial point, chosen to have a maximal first coordinate and (among such points) a maximal second coordinate, obviously is not dominated. Using it we can eliminate all other points that it does dominate, including itself. The first non-eliminated point next encountered in the algorithm necessarily has a larger value of its second coordinate. It obviously is not dominated by any point because if it were, that other point would already have been encountered. Therefore the algorithm indeed selects non-dominated points and (by induction on the number of points) it finds them all, QED.

To illustrate these concepts and this algorithm, here is plot of a dataset of $16$K records (all with integral values, shown jittered) along with its quasi-convex hull (a dark line at the upper right), its vertices marked as requested in the question, and the "dominance rectangles" of those vertices colored in. The code to do the calculations, make a scatterplot, and mark the vertices of its quasi-convex hull appears below. It includes a slightly faster R-centric version of hull.qc. It processes a million points in approximately two seconds on this machine.

#
# A faster solution (for R).
#
hull.qc <- function(xy) {
i <- order(xy[, 1], xy[, 2], decreasing=TRUE)
y <- xy[i, 2]
frontier <- which(cummax(y) <= y)
#
# Eliminate interior points along edges of the hull.
#
y.0 <- y[frontier]
frontier <- frontier[c(TRUE, y.0[-1] != y.0[-length(y.0)])]
return(i[frontier])
}
#
# Generate data.
#
library(MASS)
set.seed(17)
n <- 2^14
xy <- floor(mvrnorm(n, c(1,1), matrix(c(1, -1/2, -1/2, 1), 2))^2)
#
# Do the work.
#
system.time(frontier <- hull.qc(xy))
xy.f <- xy[frontier, , drop=FALSE]
#
# Visualization.
#
plot(xy, xlab="X", ylab="Y", main="Quasiconvex Hull")
points(xy.f, pch=19, col="Red")