Fun with maps in Stata

Today I was toying with Kevin Crow’s -shp2dta- and Maurizio Pisati’s -spmap-. -shp2dta- transforms spatial datasets in shapefile format into Stata data format and -spmap- graphs data onto maps. For starters, I browsed through Stata’s FAQ on -spmap- and Friedrich Huebler’s blog post “Guide to creating maps with Stata”, where instructions are clearly detailed. And, of course, I also checked out -help spmap- and -help shp2dta-.

It is surprising to find how easy it is to represent data onto maps using Stata. I know people who hire GIS specialists or buy expensive GIS software to draw the simplest of maps while they are armed with Stata. I have nothing against GIS experts or GIS software packages. Sure they are needed for some sophisticated GIS applications. But for very basic maps, I think -spmap- is more than enough.

The availability of basic geo-referenced database is also becoming less and less of a constraint. Global Administrative Areas (GADM), for example, provides free access to GIS datasets for almost all countries. PhilGIS* also provides free access to GIS data on the Philippines.

Using free shapefile datasets from GADM and provincial poverty data** from the Philippine National Statistics Coordination Board (NSCB), here is my first map:

The map shows poverty incidence in the Philippines by province for 2003 and 2009. Each shade of red represents a quartile (the default), where the darkest shade corresponds to the provinces with the highest incidence rates.

Here’s how I created the map:
cd c:/map    /* Change the directory to where the shapefiles are */

/* Convert shapefiles into Stata data */
shp2dta using PHL_adm1, database(phdb) coordinates(phxy) ///
    genid(id) replace

/* This created 2 Statadata files phdb.dta (master dataset) 
and phxy.dta (basemap dataset), where the variable id in the master
dataset corresponds to the variable _ID in the basemap dataset. */

use phdb.dta, clear /* Load the master dataset */
merge 1:1 ID_1 using prov.dta 
             /* Merge poverty data into the masterfile */
assert _merge==3 /* Make sure all provinces are merged */
drop _merge

/* Draw the maps */
local year 2003 2009
foreach y of local year{
    cap graph drop pov`y'
    spmap pov`y' using phxy.dta, ///
    id(id) fcolor(Reds)  ///
    title("Poverty Incidence `y'") ///
    name(pov`y') nodraw
graph combine pov2003 pov2009

  • PhilGIS has incorporated Philippine Standard Geographic Classification (PSGC) codes into the barangay (village) level dataset. This is very helpful. I have not used PhilGIS data yet because Stata returns the following error when I apply -shp2dta- using the data from PhilGIS:

st_addvar():  3300  argument out of range
read_dbf():     –  function returned error
<istmt>:     –  function returned error

I am waiting for a response from Stata Technical support. A previous Statalist post suggests that the error has something to do with strings exceeding the limit.

** I had to create a Stata data file out of the table in NSCB’s website and add the unique IDs for each province so that I can merge it later with the ‘master’ dataset. Data for Surigao del Norte and Maguindanao were used for Dinagat Islands and Shariff Kabunsuan, respectively, as they were not reported separately for these newly created provinces.

16 Responses

  1. Hi, which region corresponded to the region codes? (for NCR, ARMM etc at least)

  2. A doubt:
    I’m doing all that (with selecting the correct directory and other things necessaries.), but always occurs a error,
    >>>25: point, polyine, or polygon shapefile required”
    I thought that I’ve ever had the ‘shapefile’…

  3. […] Stata Daily. 2011. Fun with Maps in Stata. […]

  4. 1. Is there a way to assign points from one shapefile to locations within polygons from another shapefile in Stata?
    2. Is there a way to project the shapefiles in Stata?

    I am trying to fit points geocoded with geocode3 in Stata to spatial administrative units, but it seems like the projection from Google’s API is not really WGS 1984 Web Mercator. Without solving the projection problem, geocode3 is useless. Any advice?

  5. […] of Statalist for about a week now. I needed some help with the -spmap- program, which is really cool if you’re into making maps of geographic data, which is, you know, not a crazy idea for […]

  6. I had the same problem with

    st_addvar(): 3300 argument out of range

    It went away when I called

    set mem 1g

    prior to using spmap.

  7. I was wondering if you received any feedback on the error you were getting with the PhilGIS data?

    st_addvar(): 3300 argument out of range
    read_dbf(): – function returned error
    : – function returned error

    I am getting the same error when I try to load a relatively large (lots of polygons) shapefile into spmap.

    Any information would be greatly apprecaited.


  8. […] Fun with maps in Stata […]

  9. […] onto maps using Stata — and see another awesome graph — see Mitch Abdon’s “Fun with maps in Stata” over at the Stata […]

  10. I think one reason why we are so behind with what’s new is that for a long time we (at least the people I know) are not used to using Stata with the internet.* I learned the basics of Stata in 2002. That time (and many years after that) -help- (not even manuals) and friends are the only way I know how to explore Stata’s capabilities. In fact, I only learned about Statalist last year!

    *Either internet access was expensive or we were cautious not to expose our (then) pirated Stata. Sorry StataCorp, but there were copies of Stata floating around campus then (and maybe until now).

  11. Well, spmap was around long before 2008. It’s been mentioned on the Statalist, in blogposts like Friedrich’s, and it’s even in the FAQ on the Stata website.

    But still, you often need for a GIS software, e.g. to get diffrent projections, to do spatial analysis and so on.

  12. I guess Im one of those who “hired GIS specialists or buy expensive GIS software to draw the simplest of maps”… but guess what?? who knew about this stata module in 2008??
    good to know that we have this now in stata… only prob is it’s still dependent on available GIS/shape maps mostly from govt sources.

Leave a Reply