I'm trying to use Mathematica to make a map similar to the one here


i.e. a US map where each state is represented by a square of a given size, located as close as possible to its true geographical position while not overlapping with other squares. I get the data for the squares as follows:

usa = Import[
transform[s_] := StringTrim[s, Whitespace ~~ "(" ~~ ___ ~~ ")"]
polygons =
     "PlacemarkNames" /. usa[[1]]] -> ("Geometry" /. usa[[1]])];
(* Remove Alaska and Hawai *)
polygons =
  Join[{polygons[[1]]}, polygons[[3 ;; 10]], polygons[[12 ;;]]];
centers = First /@ polygons[[All, 2, 1]];
area[poly_] :=
   Flatten[First@poly Map[({1, -1} Reverse[#] &),
areas = Table[Total[area /@ polygons[[i, 2, 2 ;;]]], {i, 1, 48}];
sizes = Sqrt /@ areas;

which uses the true area of the state as square area, which is enough for now. I can plot them:

Graphics[Table[{RGBColor[RandomReal[], RandomReal[], RandomReal[]], 
   Rectangle[{centers[[i, 1]] - sizes[[i]]/2, 
     centers[[i, 2]] - sizes[[i]]/2},
    {centers[[i, 1]] + sizes[[i]]/2, 
     centers[[i, 2]] + sizes[[i]]/2}]}, {i, 1, 48}]]


Starting now, I demonstrate how far I've managed to go on my own. I want to move the squares to non-overlapping positions. I try to minimize the following score function:

dist[c1_, c2_, s1_, s2_] := Module[{dx, dy},
   dx = Abs[c1[1] - c2[1]]; dy = Abs[c1[2] - c2[2]];
   Max[0, dx - s1 - s2] + Max[0, dy - s1 - s2]
res = FindMinimum[
  Sum[Sum[dist[c[i], c[j], sizes[[i]], sizes[[j]]], {j, i + 1, 48}], {i, 1, 48}], 
  Flatten[Table[{c[i][j], centers[[i, j]]}, {i, 1, 48}, {j, 1, 2}], 1]]

but plotting this gives terrible results:

Graphics[Table[{RGBColor[RandomReal[], RandomReal[], RandomReal[]], 
    Rectangle[{c[i][1] - sizes[[i]]/2, c[i][2] - sizes[[i]]/2},
     {c[i][1] + sizes[[i]]/2, c[i][2] + sizes[[i]]/2}]}, {i, 1, 48}] /. res[[2]]]

I have not been able to see what is the issue with my code, because the score function clearly shouldn't be zero for this particular solution, but I am at loss to understand why. Any hint about that, or advice on tackling the problem in a different way, is welcome!

Doesn't work

  • 4
    $\begingroup$ Looks like a job for belisarius' answer to spread out rectangles to make them non-overlapping. Perhaps also subject it to a shape constraint (like in Heike's answer to the word cloud question) $\endgroup$ – rm -rf Sep 10 '12 at 21:00
  • $\begingroup$ @R.M yeah, I saw the word-cloud question (and answers). I'd prefer something that clearly uses Mathematica's optimizers rather than a specially crafted algorithm, as it is also intended as a demonstration of Mathematica's power to solve real-life tasks… $\endgroup$ – F'x Sep 10 '12 at 21:03
  • $\begingroup$ Did you see the link to belisarius' answer? There's a link to the code in the comments. Looks like it's a pretty straightforward application of that... I'll leave this for him to answer it for after all, he did all the hard work :) $\endgroup$ – rm -rf Sep 10 '12 at 21:06
  • $\begingroup$ I had the same buttonbar problem tonight and I saw this complaint elsewhere too. $\endgroup$ – Sjoerd C. de Vries Sep 10 '12 at 21:09
  • 2
    $\begingroup$ @R.M The code in that answer is crap, because it was done (almost) procedurally on purpose. I'll see if I can remember what I did and recast it functionally $\endgroup$ – Dr. belisarius Sep 11 '12 at 5:35

First set up the lower 48 lists, along with a set of weighted offset vectors. I'll use these for initial points in a constrained optimization. They help insofar as they can be adjusted to satisfy, or almost satisfy, the inequality constraints (which, inexplicably, I call "penalty" below).

sl48 = Take[sizes, 48];
cenl48 = Take[centers, 48];
middle = Total[cenl48]/48;
wmean = sl48.cenl48/Total[sl48]
offsets = Map[# - wmean &, cenl48];

Our function will measure how far the new centers are from the initial ones. This is needed because otherwise we tend to get a result that is a tight fit of weighted squares, but bears little resemblance to the USA map.

f[centers : {{_?NumberQ, _?NumberQ} ..}, 
  orig : {{_?NumberQ, _?NumberQ} ..}] := With[
  {diffs = centers - orig},
  Total@Map[#.# &, diffs]]

overlap[c1_, c2_, s1_, s2_] := Min[(s1 + s2)/2 - Abs[c1 - c2]] <= 0

penalty[centers_, wts_] :=

   overlap[centers[[i]], centers[[j]], wts[[i]], wts[[j]]], {i, 
    47}, {j, i + 1, 48}]]

cvars = Array[c, {48, 2}];

 res = FindMinimum[{f[cvars, cenl48], penalty[cvars, sizes]}, 
    Transpose[{Flatten[cvars], Flatten[cenl48 + offsets]}], 
    MaxIterations -> 50, Method -> "InteriorPoint"];]

(* Out[344]= {136.37, Null} *)

We can plot this as follows.

colors = Table[
   RGBColor[RandomReal[], RandomReal[], RandomReal[]], {48}];
rex2 = Transpose[{(cvars /. res[[2]]) - sl48/2, (cvars /. res[[2]]) + 
Graphics[Transpose[{colors, Rectangle @@@ rex2}]]

enter image description here

  • 2
    $\begingroup$ I guess I should mention that this can probably be done as well or better, and much faster, if recast as a linear programming problem using Manhattan distances rather than Euclidean in the objective function (and slightly changing the constraints an explicitly linear formulation). $\endgroup$ – Daniel Lichtblau Sep 12 '12 at 15:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.