Skip to content

basic polyfon shapefile subsetting

2 messages · Chris Fowler, Roger Bivand

Original
Chris Fowler
# 
I swore I would figure this out without going to the list, but I just 
don't understand the data structure well enough to make this work. Easy 
pickings for many of you on this list.

I have an ArcGIS shapefile "base.shp" that has 23034 polygons (U.S. 
census tracts). I want to end up with a subset of this file that has no 
'islands' and a neighbor list to accompany it. What am I doing wrong?

baseShp<-readShapePoly("Export_Output.shp")	#read in a shapefile
Q1B<-poly2nb(baseShp, row.names = baseShp$FIPS)	#1st order Queen
summary(Q1B)	#731 with 0 neighbors we need to remove these
x<-card(Q1B)	#create a vector with the number of neighbors
baseShp$numNeighbors<-x	#assign that vector as an additional data column
length(baseShp$numNeighbors)	#check to make sure all 23034 have value
summary(baseShp$numNeighbors)	#...and that those values are what I 		
				#expect
submap <- baseShp[baseShp$attr.data$numNeighbors >=1] #pull out those
				# with at least 1 neighbor
plot(submap)		#something plots....
length(submap$numNeighbors)  #0, no data transferred
subm<-as.data.frame(submap)  # convert to data frame for easier looking
str(subm)  #we still have 23034 observations and now we have no data


So, why doesn't my subseting code bring data with it and why doesn't it 
actually do the subsetting?

Thanks in advance,

Chris
Roger Bivand
#
On Fri, 23 Oct 2009, Chris Fowler wrote:

            

      
      
      
    

        
Use

baseShp[baseShp$numNeighbors >=1, ]

There are two mistakes, one the missing comma, the other the use of a 
non-existing column, which will result in baseShp[NULL >=1] being 
evaluated.

Roger