30
$\begingroup$

I'd like to plot county data as in this D3 example:

enter image description here

Is there a way to do this with GeoGraphics, or at least to get path polygons per county and then plot them with different colors?

The first thing I tried was this post which didn't work. Here's a snippet of my data:

(*county    value*)
49005   0.5742412852924035
49007   0.8317715755644481
49009   0.5453979492171142
49011   0.28254890319856485
49013   0.7479376357589627
49015   0.18711805968098272
49017   0.7457292949804695
49019   0.567573145972929
49021   0.16246935549286112
49023   0.9195138046651179
49025   0.4367034052124463
49027   0.6982473831931622
49029   0.8338046457979154
49031   0.9912431029895752
49033   0.7835820925255039
49035   0.1239835744258424
49037   0.6332753189833022
49039   0.1020405619672482
49041   0.6598660862197654
49043   0.292093295985822

Sadly these are missing:

EntityValue[
 Entity["AdministrativeDivision", {"AbbevilleCounty", "SouthCarolina",
    "UnitedStates"}], 
 EntityProperty["AdministrativeDivision", "FIPSCode"]]
$\endgroup$

1 Answer 1

30
$\begingroup$

Updated text

I made a tutorial out of this answer and posted it on Wolfram Community. It's slightly different, and may be worth to check out.

Data from Wolfram Knowledgebase

The most straightforward way to plot data about counties is to use the built-in county entities and GeoRegionValuePlot.

counties = Cases[
   EntityList@Entity["AdministrativeDivision", {"Country" -> Entity["Country", "UnitedStates"]}],
   Entity[_, {_, _, _}] (* separates states from counties *)
   ];

GeoRegionValuePlot[counties -> "ForeignBornFraction"]

US county region value plot

Here ForeignBornFraction is an entity property to adminstrative divisions.

Creating a US county entity store

Wouldn't it be nice if we could visualize our own data like this?

GeoRegionValuePlot[counties -> "CustomUnemploymentData"]

In the just released Mathematica 11 we can achieve that goal with entity stores. Since the built-in county entities don't have county FIPS codes and we need FIPS codes to correlate the entities with our data, I'll take the FIPS codes that JasonB shared in his answer, and make my own entity system.

fipsRules = Import["http://pastebin.com/raw/YSsKAbjn"]; (* JasonB's FIPS codes file, see his answer for details. *)
testdata = Rest@Import["https://gist.githubusercontent.com/mbostock/4060606/raw/25385f68a3be5c9dbe36af27fc2498fb2aab6bc0/unemployment.tsv", "TSV"];

getEntity[fipsRules_][{fips_, unemployment_}] := 
 With[{entity = (IntegerString[fips] /. fipsRules)}, 
  If[MatchQ[entity, _String],
   Nothing, (* If the FIPS code does not exist in fipsRules *)
   getEntityIdentifier[entity] -> <|
     "Label" -> Style[entity["Name"], Bold],
     "ParentEntity" -> entity,
     "FIPSCode" -> fips,
     "Unemployment" -> unemployment
     |>
   ]]

getEntityIdentifier[Entity[_, {county_, state_, country___}]] := county <> "-" <> state
getEntityIdentifier[Entity[_, {"DistrictOfColumbia", _}]] := "DistrictOfColumbia"

countyEntities = Apply[Association, getEntity[fipsRules] /@ testdata];

So far we've imported both the unemployment data and the FIPS codes and used that to create an Association in the form dictated by EntityStore. Each entity in our custom entity store has the properties ParentEntity, FIPSCode and Unemployment. Label is also present, and is used for presentation.

ParentEntity is the built-in county entity and I will use that and computed properties to effect a sort of inheritance.

store = EntityStore["USCounty" -> <|
     "Entities" -> countyEntities,
     "Properties" -> <|
       "ParentEntity" -> <|
         "Label" -> "Parent entity"
         |>,
       "FIPSCode" -> <|
         "Label" -> "County FIPS code"
         |>,
       "Polygon" -> <|
         "DefaultFunction" -> (#["ParentEntity"]["Polygon"] &),
         "Label" -> "Polygon"
         |>,
       "Position" -> <|
         "DefaultFunction" -> (#["ParentEntity"]["Position"] &),
         "Label" -> "Position"
         |>,
       "Latitude" -> <|
         "DefaultFunction" -> (#["ParentEntity"]["Latitude"] &),
         "Label" -> "Latitude"
         |>,
       "Longitude" -> <|
         "DefaultFunction" -> (#["ParentEntity"]["Longitude"] &),
         "Label" -> "Longitude"
         |>,
       "HasPolygon" -> <|
         "DefaultFunction" -> (#["ParentEntity"]["HasPolygon"] &),
         "Label" -> "Has polygon?"
         |>,
       "Name" -> <|
         "DefaultFunction" -> (#["ParentEntity"]["Name"] &),
         "Label" -> "Name"
         |>
       |>,
     "EntityClasses" -> <|
       "Alaska" -> <|
         "Entities" -> ("ParentEntity" -> (MatchQ[#, Entity[_, {_, "Alaska", _}]] &)),
         "Label" -> Style["Alaska", Bold]
         |>,
       "Hawaii" -> <|
         "Entities" -> ("ParentEntity" -> (MatchQ[#, Entity[_, {_, "Hawaii", _}]] &)),
         "Label" -> Style["Hawaii", Bold]
         |>,
       "Mainland" -> <|
         "Entities" -> ("ParentEntity" -> (MatchQ[#, Entity[_, {_, Except["Alaska" | "Hawaii"], _}]] &)),
         "Label" -> Style["Mainland", Bold]
         |>
       |>
     |>];

The inheritance that I'm talking about is the properties that use DefaultFunction to compute their values. Besides the properties I also define three classes of counties, one for Alaska, one for Hawaii and one for the US mainland. The reason for this is that later I'll want to manipulate the three areas separately; I'll want to rescale and translate Alaska and Hawaii to the positions which they have in the map in the original post.

In order to be able to use the entity store you'll want to add it to $EntityStores:

AppendTo[$EntityStores, store]

I would also recommend adding as a resource object:

ro = ResourceObject[<|
    "Name" -> "USCounty",
    "ResourceType" -> "DataResource",
    "Content" -> store
    |>];
ResourceData[ro]; (* Needed to call this, otherwise the next line did not work. *)
ro["ResourceLocations"]

The advantage of this is that resource objects should be, I think, that they are persistent, i.e. next time you start up Mathematica you can get ahold of your data with

ResourceObject["USCounty"]

and add it to $EntityStores using

AppendTo[$EntityStores, ResourceObject["USCounty"]]

But right now this doesn't seem to work. It works within one session, but in a new session the resource object has disappeared. This feature as well is new in Mathematica 11. We can now do things like

EntityList["USCounty"] // Shallow

Mathematica graphics

EntityList@EntityClass["USCounty", "Hawaii"]

Mathematica graphics

Back to GeoRegionValuePlot

With the entity store in place we can now use GeoRegionValuePlot to plot our custom data like this:

GeoRegionValuePlot[EntityList[EntityClass["USCounty", "Mainland"]] -> "Unemployment"]

Mathematica graphics

With a few options it looks almost the same as the image in the question:

getColor[val_] := Which[
  0 < val < 0.01667, RGBColor[{247, 251, 255}/255],
  0.0166 < val < 0.0333, RGBColor[{222, 235, 247}/255],
  0.0333 < val < 0.05, RGBColor[{198, 219, 239}/255],
  0.05 < val < 0.0666, RGBColor[{158, 202, 225}/255],
  0.0666 < val < 0.0833, RGBColor[{107, 174, 214}/255],
  0.0833 < val < 0.1, RGBColor[{66, 146, 198}/255],
  0.1 < val < 0.1166, RGBColor[{33, 113, 181}/255],
  0.1166 < val < 0.1333, RGBColor[{8, 81, 156}/255],
  0.1333 < val < 0.15, RGBColor[{8, 48, 107}/255],
  True, RGBColor[{8, 48, 107}/255]
  ]

GeoRegionValuePlot[
 EntityList[EntityClass["USCounty", "Mainland"]] -> "Unemployment",
 ColorFunction -> getColor,
 ColorFunctionScaling -> False,
 PlotLegends -> None,
 PlotStyle -> Directive[EdgeForm[None]],
 GeoBackground -> None,
 GeoProjection -> {
   "LambertAzimuthal",
   "Centering" -> GeoPosition[{30, -195/2}]
   },
 PlotRange -> {{-0.37, 0.38}, {-0.13, 0.38}},
 ImageSize -> 1000
 ]

Mathematica graphics

The last details

The previous map still lacks the white borderlines, Hawaii, and Alaska. For Hawaii and Alaska we need to use Inset, so we have to use GeoGraphics. Hence, below is a GeoGraphics version of the previous plot.

getPolygon[entity_] := {
  GeoStyling[None],
  getColor[entity["Unemployment"]],
  Polygon[entity]
  }

alaska = GeoGraphics[
   getPolygon /@ EntityList[EntityClass["USCounty", "Alaska"]], 
   GeoBackground -> None];
hawaii = GeoGraphics[
   getPolygon /@ EntityList[EntityClass["USCounty", "Hawaii"]], 
   GeoBackground -> None];

map = GeoGraphics[{
    getPolygon /@ EntityList[EntityClass["USCounty", "Mainland"]],
    Inset[alaska, {-0.3, -0.06}, {0, 0}, 0.3],
    Inset[hawaii, {-0.15, -0.06}, {0.1569, -0.0436}, 0.3]
    },
   GeoBackground -> None,
   GeoProjection -> {
     "LambertAzimuthal",
     "Centering" -> GeoPosition[{30, -195/2}]
     },
   PlotRange -> {{-0.37, 0.38}, {-0.13, 0.38}},
   ImageSize -> 1000
   ]

Mathematica graphics

To add the borderlines I first make a black and white mask, where states are black and the background is white. The state borderlines that I want to draw are white as well; I then use ImageAdd to put the images together. The counties will be unaffected because we'll only be adding black to them, except for on the borderlines.

states = Cases[
   EntityList@Entity["AdministrativeDivision", {"Country" -> Entity["Country", "UnitedStates"]}],
   Entity[_, {_, _}]
   ];

alaska = GeoGraphics[
   Polygon /@ EntityList[EntityClass["USCounty", "Alaska"]], 
   GeoBackground -> None];
hawaii = GeoGraphics[
   Polygon /@ EntityList[EntityClass["USCounty", "Hawaii"]], 
   GeoBackground -> None];

borders = GeoGraphics[{
    FaceForm[Black],
    EdgeForm[White],
    GeoStyling[None],
    Polygon /@ states,
    Inset[alaska, {-0.3, -0.06}, {0, 0}, 0.3],
    Inset[hawaii, {-0.15, -0.06}, {0.1569, -0.0436}, 0.3]
    },
   GeoBackground -> None,
   GeoProjection -> {
     "LambertAzimuthal",
     "Centering" -> GeoPosition[{30, -195/2}]
     },
   PlotRange -> {{-0.37, 0.38}, {-0.13, 0.38}},
   ImageSize -> 1000
   ];

ImageAdd[map, borders]

Mathematica graphics

$\endgroup$
9
  • $\begingroup$ Is there an easy way to go from the Fips code to names and back? $\endgroup$
    – M.R.
    Aug 7, 2016 at 1:45
  • $\begingroup$ @M.R. I don't understand what the situation you are thinking about is, the way I see it you have them both at the same time so you don't don't have to go from one to the other. $\endgroup$
    – C. E.
    Aug 7, 2016 at 8:31
  • $\begingroup$ My situation is that I have a 2 column csv that contains (FIP, Value) pairs and I'd like to plot them but I don't know how to go from FIP to county in Mathematica's geo plots. $\endgroup$
    – M.R.
    Aug 7, 2016 at 20:41
  • $\begingroup$ Trying to plot data in census blocks as well, I'll add some example data to the question. $\endgroup$
    – M.R.
    Aug 7, 2016 at 21:25
  • $\begingroup$ see my updated question... $\endgroup$
    – M.R.
    Aug 7, 2016 at 21:29

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.