I downloaded administrative shapefiles from www.gadm.com/country (China). The Python (Basemap) code below tells me that there are 2365 polygons in it. However the .csv and .dbf files included with the shapefiles tell me that there are only 344 administrative divisions at this level (which corresponds to real geography).

Forgive my lack of experience with shapefiles, but what could explain this discrepancy? Is it possible that this would be caused by small islands being considered a different polygon? If so, are there specific techniques to get rid of them?

My worry is that this increases 8-fold the computation time of point-to-polygon spatial joints, since there are 8 times more polygons.

from mpl_toolkits.basemap import Basemap

m = Basemap(projection='merc',llcrnrlat=15,urcrnrlat=55,llcrnrlon=70,urcrnrlon=140,lat_ts=20,resolution='c')
m.readshapefile(china_shapefile, 'provinces', drawbounds=True)

print len(m.provinces)
  • 1
    Pure guesses here: could the shapefile also include administrative areas on a more detailed level as well? Do you get any hints by looking at the attribute table, maybe a column that displays the type of area? What does the polygons look like - like you would expect or covering a larger area/more detailed etc? When it comes to computation time, 2000 polygons is still a pretty decent amount and isn't likely to be much more of a problem in overlay analyses compared to 300. – Martin Nov 30 '15 at 9:49
  • 1
    Summing up my last comment: computation time isn't what you should worry about, it's whether you have relevant data for your analysis. – Martin Nov 30 '15 at 9:51
  • It seems that the detail level is the one expected, but there are indeed a lot of islands, which I know assume are counted as separate polygons (I wasn't sure if a polygon had to be fully connected or not). I don't know much about how the inner workings of these spatial join algorithms, but do you mean that getting rid of the small islands (90% of polygons) won't change much in speed, because although they are 90% of polygons, they are a much smaller fraction of the number of points that make up the polygons? Is there any other way to speed this up? – Alexis Eggermont Nov 30 '15 at 10:15
  • 1
    If you look at the attribute table, how many features have you got there? That should be what counts. You sure could have islands and even overlapping polygons (that can be really hard to spot - run Intersect with only one input). A polygon is basically defined as a surface, and can have holes in it, so if the outline isn't continuous all way round it just can't be a polygon. I doubt that it will make much difference in computation time, but it is of course hard to say without testing or at least seeing the data and knowing exactly what analysis you want to perform. – Martin Nov 30 '15 at 11:27
  • All solved. The islands were indeed the extra polygons, and I was able to multiply processing speed by a factor of 1000 by simplifying the shapefiles using mapshaper.org. The devil is in the detail level. Thanks for the help. – Alexis Eggermont Nov 30 '15 at 12:29

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.