Draw maps in R
Use R to draw maps for your paper
Very often I have found researchers using google map or low resolution map to show population distributions (or survey regions). They usually don’t appear that great in a research journal atleast when printed out.
This is a step by step approach to draw maps in R. I demonstrate how to use functions in PBSmapping package of R. You can draw any map with two small piece of R code. The following are required (Assuming you have Windows):
1. R installed in your computer
2. PBSmapping – R package (and all its dependencies)
3. GSHHS Database file
Import GSHHS Database
GSHHS refers to “Global Self-consistent, Hierarchical, High-resolution Shoreline”. They present a high-resolution shoreline data set maintained by NOAA and are available to download for free. More details are available at : http://www.soest.hawaii.edu/wessel/gshhs/
The shorelines come in five resolutions:
- full resolution: Original data resolution.
- high resolution: About 80 % reduction.
- intermediate resolution: Another ~80 % reduction.
- low resolution: Another ~80 % reduction.
- crude resolution: Another ~80 % reduction.
They are each organized into 4 hierarchical levels:
- L1: boundary between land and ocean.
- L2: boundary between lake and land.
- L3: boundary between island-in-lake and lake.
- L4: boundary between pond-in-island and island.
STEP 1 (Download data)
Find the root directory of PBSmapping R package in your Windows. For example in my computer, the location is ‘C:\Users\Pazhayamadom\Documents\R\win-library\2.15\PBSmapping’
Download ‘GSHHS database in binary format’ (zip file) from the website.
Extract the contents into the root directory of PBSmapping.
The contents include database files (Data for world map) – gshhs_f.b, gshhs_h.b, gshhs_i.b ….etc.
STEP 2 (Load PBSmapping package)
Use the following code to load PBSmapping package in R
STEP 3 (Import data into R)
To import the required data (for map), you should have
1. xlim= range of X-coordinates or Longitudes (for clipping). The range should be between 0 and 360.
2. ylim= range of Y-coordinates or Latitudes (for clipping).
3. maxlevel= maximum level of polygons to import: 1 (land), 2 (lakes on land), 3 (islands in lakes), or 4 (ponds on islands).
The importGSHHS( ) function can be used for importing data (demonstrated with Map of India)
In the above function, there are 4 inputs or arguments seperated by commas. The data file for high resolution map is specified as gshhs_h.b with the source path or location (use gshhs_f.b file for full resolution). xlim are the minimum and maximum longitudes required. Similarly, ylim are the minimum and maximum latitudes required. maxlevel=4 means to import data for lakes, island and ponds along with land.
This process may take some time to read required data from the file and you should be getting the following message:
importGSHHS status: --> Pass 1: complete: 1493 bounding boxes within limits. --> Pass 2: complete. --> Clipping...
STEP 4 (Draw the map)
The final step is to plot the data which is now imported from the database file.
The plotMap ( ) function is used for this purpose.
R have the option of saving pictures in the format of JPEG, PNG, TIFF, Metafile, Postscript, Bmp and PDF. I prefer Metafile for MS Word, PNG for Latex and Postscript for graphical works (the best choices on my experience).
There are many optional arguments (inputs) for the functions mentioned above. I leave this for you to explore. But a few examples are given below:
To get grey shades for the land, use the plotMap( ) function with an extra argument col=’beige’.
plotMap (polys, col="beige")
If you would like the ocean (background) with a blue shade, use the following code:
plotMap (polys, col="beige", bg="lightblue")
Hope this was helpful. Thanks for reading the post.
P.S. – Have fun with changing the colour of your choice 🙂
Add Rivers and Borders
Here we import the database for rivers and borders (wdb_rivers_f.b and wdb_borders_f.b ).
The above function add rivers and borders to the map. This works only if the map is already plotted using plotMap( ) function.
You can add labels into the map using text ( ) function. Again this will work only if the map is already plotted.
text (70, 10, "Arabian Sea")
Here 70 and 10 are the x and y coordinates in the map / graph.