Bounding data by two parallel lines with minimum distance between them I have a set of data samples that approximately follow a straight line in 2D. I need to find two parallel lines that are spaced as close as possible such that all of the samples lie between the lines.
I'll need to justify how these lines were chosen to meet this requirement, so I'd like to understand the solution in a somewhat formal way. (ButI'll take what I can get for now!)
Previously, I used linear regression to get a best fit a line, and then found two lines that are parallel to the fit line. I noticed that these two lines are not equally distant from the best fit line. This led me to wonder whether a tighter pair of parallel lines might be found (with rigor) to bound the data.
I'm pressed for time on this, so I am hoping someone can point me toward where I can find an answer to this. (Or to provide an answer would be even more helpful!)
Responding to the request for clarification:
The data points have both systematic (non-random) and random deviation from some ideal line. I’m not interested in estimating the parameters of that ideal line, though.  We can assume data is unitless (it is actually the same units on each axis). Distance is Euclidean and I am only interested in minimizing the distance between the bounding lines such that all samples lie between them.
 A: It is immediate (from the definitions) that the lines must pass through extremal points of the point set.  Because at least one of them must contain at least two of the points (provided there is more than one point!), at least one of them will be determined by an edge of its convex hull.  Thus, for $n$ points with $m$ extremals, there is an algorithm requiring $O(n(\log(n)+m))$ time consisting of finding the convex hull and testing its edges.  
(This is very nearly the same problem as finding a minimum-area bounding rectangle for the points and the solution below is a modification of the one I posted at https://gis.stackexchange.com/a/22934/664.)
It's more convenient in an array-oriented language like R to be a little less than optimal: for each edge of the convex hull, rotate the hull so that edge lies on the x-axis and all points lie in the upper half plane.  Record the largest y-coordinate: this is the width of the point cloud orthogonal to that edge.  An edge for which this width is minimal determines one line and the corresponding width determines the oriented distance between the lines, thereby specifying the other line.
Here are some randomly-generated examples.  The red line is determined by the optimal edge and the gray line is determined by the oriented distance from the red line: typically it passes through just one extremal point.

Here is the R code to generate a solution given an $n\times 2$ array of coordinates p.  It assumes there is more than one distinct point in p.
span <- function(p) {
  a <- chull(p)                     # Indexes of extremal points, negatively oriented
  e <- p[c(a[-1], a[1]), ,drop=FALSE] - p[a, ] # Edge vectors
  e <- e / sqrt(rowSums(e^2))       # Unit edge dirction vectors
  n <- cbind(e[, 2], -e[, 1])       # Unit normal vectors
  w <- apply(tcrossprod(n, p[a,]), 1, function(x) max(x) - min(x))  # Widths
  i <- which.min(w)                 # Index (into `a`) of best edge
  list(origin=p[a[i],], direction=e[i,], normal=n[i,], width=w[i])
}

It returns a point on one line (origin), a unit direction vector for the line (direction), a unit normal pointing towards the second line (normal), and the distance between those lines (width).
Here is an example of its use.
n <- 4
set.seed(17)
par(mfrow=c(1,n))
invisible(replicate(n, {
  # Create sample data.
  p <- matrix(rnorm(9*2), ncol=2)

  # Find the two lines.
  s <- span(p)

  # Convert the result into parameters for `abline`.
  slope <- s$direction[2] / s$direction[1]
  origin.1 <- s$origin
  origin.2 <- s$origin + s$width * s$normal

  # Display everything.
  plot(p, asp=1)
  if (s$direction[1] != 0) {
    abline(a=origin.1[2]-slope*origin.1[1], b=slope, col="Red")
    abline(a=origin.2[2]-slope*origin.2[1], b=slope, col="Gray")
  } else {
    abline(v = c(origin.1[1], origin.2[1]), col=c("Red", "Gray"))
  }
}))
par(mfrow=c(1,1))

This solution is reasonably fast: on this nine-year-old workstation, it takes one second when given ten million normally-distributed points.
