Main Content

Create Choropleth Map of Population Density

This example shows how to create a choropleth map of population density for the six New England states in the year 2000.

Import low-resolution U.S. state boundary polygons, setting the map limits for the New England region.

MapLatLimit = [41 48];
MapLonLimit = [-74 -66];

NEstates = shaperead('usastatelo', 'UseGeoCoords', true, ...
   'BoundingBox', [MapLonLimit' MapLatLimit']);

Set up map axes with a projection suitable to display the New England states.

axesm('MapProjection', 'eqaconic', 'MapParallels', [],...
  'MapLatLimit', MapLatLimit, 'MapLonLimit', MapLonLimit,...
  'GLineStyle', '-')
geoshow(NEstates, 'DisplayType', 'polygon', 'FaceColor','green')

Identify the maximum population density for New England states.

maxdensity = max([NEstates.PopDens2000])
maxdensity = 1.1345e+03

Create an autumn colormap for the six New England states, and then use the flipud command to invert the matrix.

fall = flipud(autumn(numel(NEstates)));

Make a symbol specification structure, a symbolspec, that assigns an autumn color to each polygon according to the population density.

densityColors = makesymbolspec('Polygon', {'PopDens2000', ...
   [0 maxdensity], 'FaceColor', fall});

Display the map.

geoshow(NEstates, 'DisplayType', 'polygon', ...
   'SymbolSpec', densityColors)
title ({'Population Density in New England in 2000', ...
   'in Persons per Square Mile'})

Add a colorbar. You can also experiment with other colormaps.

caxis([0 maxdensity])
colormap(fall)
colorbar