The is already implemented in the GEOS library (`GEOSNearestPoints_r`) and is now wrapped in the `rgeos` package (>= 0.3-21), see `gNearestPoints`: R> gDistance(rect, poly.lake) [1] 33806 R> gNearestPoints(rect, poly.lake) class : SpatialPoints features : 2 extent : 1036613, 1070419, 825354, 825354 (xmin, xmax, ymin, ymax) coord. ref. : NA R> coordinates(gNearestPoints(rect, poly.lake)) x y 1 1036613 825354 2 1070419 825354