23
$\begingroup$

From an external simulation program, I have lots of particle tracks which follow random walks, all of which start at or near the origin $(0,0)$. This means that when I plot the walks in Mathematica, it looks a bit of a mess because all the tracks overlap with each other.

For example:

randomwalks = Table[
    Accumulate@RandomVariate[StableDistribution[1, 1.7, 0, 0, 2], {512, 2}], 
    {i, 3}];

ListLinePlot[randomwalks, PlotRange -> Full, PlotTheme -> "Minimal", 
    Axes -> None, Frame -> True, AspectRatio -> 1];

This gives the image on the left - the fact that all the walks start near the origin means that they overlap and it's hard to see what is going on.

enter image description here

What I'd like to plot instead is something like the image on the right, where the random walks have been manually translated to stop them overlapping.

Is there a way to achieve this translation automatically with Mathematica?

The direction of each jump in the random walk is important, so I want to keep this in the visualization (hence translations only and no rotations).

Note also the choice of a Lévy-stable distribution for the random walk rather than the normal distribution. This typically gives walks that are more spread-out with longer jumps (due to the fat distribution tails). In turn this leads to a more diverse spread of shapes to try and fit together, and means arranging the shapes in a simple grid (e.g. via Column[Graphics /@ Line /@ randomwalks] is not always as successful as I'd like.


My first thought was to calculate the bounding rectangle for each track, and then use the code in An algorithm to space out overlapping rectangles? to space out those bounding rectangles. However, this might leave large regions of whitespace due to the nature of the random walks, so perhaps a bounding polygon might be better to give a tighter arrangement of the tracks.

Given I have the track data itself, is there perhaps a way to do it without using the bounding rectangles? As suggested in a comment, perhaps using the convex hull of each path and then packing the resulting polygons?

regions = ConvexHullMesh[#] &/@ randomwalks
centroids = RegionCentroid[#] &/@ regions

Once we have the centroids, one approach would be to minimize the Euclidean distance between each centroid whilst ensuring the regions don't overlap - but that would need a function such as RegionOverlapQ[] or similar.

$\endgroup$
  • 5
    $\begingroup$ One idea: Take the convex hull of each path, and then your problem is reduced to packing convex polygons. $\endgroup$ – Joseph O'Rourke Sep 10 '15 at 11:23
  • 5
    $\begingroup$ Caveat: Optimal packing of convex polygons into a rectangular box is NP-hard, so you will be forced to use heuristics. (Mma has ConvexHull[ ].) $\endgroup$ – Joseph O'Rourke Sep 10 '15 at 11:34
  • 2
    $\begingroup$ Why not just Column[Graphics /@ Line /@ randomwalks] or something? $\endgroup$ – Kuba Sep 10 '15 at 11:41
  • 2
    $\begingroup$ @Kuba because I fancied a challenge? :-) $\endgroup$ – dr.blochwave Sep 10 '15 at 11:43
  • 2
    $\begingroup$ OMG My code there is mess $\endgroup$ – Dr. belisarius Sep 10 '15 at 14:21
13
$\begingroup$

One approach is to rasterize the lines separately and use correlation to measure overlap. The lines are positioned one by one, as close to the centre as possible without overlapping any previously positioned lines. Heike used this in her answer to the word cloud question and I borrowed the idea here.

Here's a quick & dirty implementation:

randomwalks = 
  Table[Accumulate@
    RandomVariate[StableDistribution[1, 1.7, 0, 0, 2], {512, 2}], {i, 15}];

r = 2.05 Max@Abs@randomwalks;

bin[walk_] := Binarize@ColorNegate@Erosion[
    Rasterize[Graphics[Line@walk, PlotRange -> r], ImageSize -> 500], 
    DiskMatrix[3]]

radial = RadialGradientImage[Black, {500, 500}];

iter[im1_, im2_] := Module[{corr, minx, miny},
  corr = ImageAdd[radial, ImageCorrelate[im1, im2]];
  {minx, miny} = 250 - PixelValuePositions[corr, Min@ImageData@corr][[1]];
  Sow[-{minx, miny} r/250];
  ImageAdd[im1, Image@RotateRight[ImageData@im2, {miny, -minx}]]]

offsets = Reap[Fold[iter, bin /@ randomwalks]][[-1, 1]] ~Prepend~ {0, 0};

ListLinePlot[
 Transpose[(Transpose[randomwalks, {1, 3, 2}] + offsets), {1, 3, 2}], 
 PlotRange -> Full, PlotTheme -> "Minimal", Axes -> None, 
 Frame -> True, AspectRatio -> 1]

enter image description here

$\endgroup$
  • 1
    $\begingroup$ Very much different from mine, but I am confident there are many different approaches on this problem (as already outlined on the question). $\endgroup$ – kirma Sep 12 '15 at 15:41
  • $\begingroup$ This one has ended up being more appropriate for my needs, although @kirma's answer is a very interesting idea nonetheless! $\endgroup$ – dr.blochwave Sep 21 '15 at 12:29
11
$\begingroup$

This solution is primarily a proof of concept; its run-time complexity makes it impractical for large amounts of items. Large amount means more than six, in this case.

Nonetheless, the idea is the following: individual plots are converted to polynominoes to be laid out on a rectangular grid. Every possible alignment on a finite grid is encoded to a Boolean equation representing legal layouts (every piece is present on exactly one position, and the pieces don't overlap), and a satisfying layout instance is searched for this equation. Result is baked back to datasets as offsets, which can be plotted now without overlaps.

ClearAll[layOutDatasets];

layOutDatasets[datasets_List, numinstances_Integer, 
   resolution_Integer, maxdim_List] :=
  Module[
   {normalizeDatasets, generateTile, generateAlignments, 
    generateOrientations, atMostOne, exactlyOne},

   normalizeDatasets[
     dss_List] := ((With[{offset = 
            First@CoordinateBoundingBox@#}, # - offset & /@ #]) & /@ 
       dss)/Max[Max@Abs[Subtract @@ CoordinateBoundingBox@#] & /@ dss];

   generateTile[ds_List] :=
    ImageData[
        ImageCrop@
         Binarize[ImageResize[#, resolution + 2], {0.99, 1}]] /. {0 ->
          True, 1 -> False} &@
     ImageReflect@
      Rasterize[
       Style[Graphics[{Thickness[1/(2 resolution)], Line@ds}, 
         PlotRange -> {{-1/resolution, 
            1 + 1/resolution}, {-1/resolution, 1 + 1/resolution}}], 
        Antialiasing -> False], ImageSize -> 8 (resolution + 2)];

   generateAlignments[patt_List] :=
    With[{maxoff = maxdim - Dimensions@patt},
     Flatten[Table[{i, j}, {i, 0, maxoff[[1]]}, {j, 0, maxoff[[2]]}], 
      1]];

   generateOrientations[patt_List] :=
    PadRight[PadLeft[patt, Dimensions[patt] + #, False], maxdim, 
       False] & /@ generateAlignments[patt];

   atMostOne[v_List] := BooleanCountingFunction[{0, 1}, Length@v] @@ v;
   exactlyOne[v_List] := BooleanCountingFunction[{1}, Length@v] @@ v;

   Module[{normdatasets, pieces, placemasks, sols, offsols},
    normdatasets = normalizeDatasets[datasets];
    pieces = generateTile /@ normdatasets;

    placemasks = generateOrientations[#] & /@ pieces;

    With[{vars = 
       Flatten@MapIndexed[c@First@#2 /@ Range@Length@#1 &, placemasks],
      varmaps = 
       Flatten@MapIndexed[
         c[#2[[1]]][#2[[2]]] -> Reverse@#1/resolution &, 
         generateAlignments[#] & /@ pieces, {2}]},

     sols = 
      SatisfiabilityInstances[
       And @@ (atMostOne /@ 
           Flatten[
            Transpose[
             Flatten[
              MapIndexed[#1 && c[#2[[1]]][#2[[2]]] &, 
               placemasks, {4}], 1], {3, 1, 2}], 1]) &&

        And @@ MapIndexed[
          With[{constants = c@First@#2 /@ Range@Length@#1},
            exactlyOne@constants] &, placemasks],
       vars, numinstances];

     offsols = (Pick[vars, #, True] /. varmaps) & /@ sols;

     MapThread[
        Function[{set, off}, # + off & /@ set], {normdatasets, #}] & /@
       offsols]]];

An example of use:

ListLinePlot[#, Axes -> False, AspectRatio -> Automatic] &@
 First@layOutDatasets[
   Table[Accumulate@
     RandomVariate[StableDistribution[1, 1.7, 0, 0, 2], {512, 2}], {i,
      5}], 1, 5, {15, 15}]

enter image description here

Another run:

enter image description here

EDIT: Note that there is no reason this method would be limited to a rectangular grid, regular alignments or even specific rotations. Generating more varied alternatives is always an option (although it increases size of the search space, which increases time complexity exponentially), and generating random alignments instead of regular lattice alignment is also an option. Hexagonal lattice would be quite interesting, but "rasterizing" graphics to such a grid is somewhat cumbersome.

$\endgroup$
  • $\begingroup$ This looks really cool - I'll give it a go later on and report back properly! $\endgroup$ – dr.blochwave Sep 12 '15 at 11:29
  • $\begingroup$ @blochwave It's more than a bit unpolished and has hacks all around - but I decided to dump it here anyway, because polishing might have never happened... The core idea is pretty simple, though. $\endgroup$ – kirma Sep 12 '15 at 11:33
  • 1
    $\begingroup$ you never know, community action might help polish it! $\endgroup$ – dr.blochwave Sep 12 '15 at 11:36

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.