Calculating distance to nearest border between many points and many candidate polygons

[Shir Bloch]

Hi everyone, seeking help for a spatial RD. [I couldn't find the answer I was looking for in existing posts, but if I missed something please let me know.]

I have data with coordinates for a bunch of structures in the US. I also have shape files for all the census tracts in the US, as well as an eligibility indicator (some tracts are eligible, others are not). For each of my structures, I want to calculate the distance to the nearest eligibility border (i.e., where an eligible tract meets an ineligible tract), to see if these structures are placed selectively on the eligible side of the border. Is there a way to do this in Stata? There may be an efficiency concern here as I have about 130k structures, and 84k tracts. Is there a way to do this without calculating the distance between each structure and each tract's nearest border (and then minimizing this distance across all estimated structure-tract pairs for each structure)?

The difference between this task and eg the thread about the smallest distance from a point to a border (and the helpful reply about the Texas/Mexico border) is that my eligible tracts are many and not necessarily contiguous, so I am not just finding the shortest distance to one single polygon but rather to a collection of them.

Thank you in advance! Looking forward to any response.

[Raymond Guiteras]

I have done something similar using user-written geoframe and geonear. However, it is a very brute-force approach -- to find the nearest neighbor to a point in file A among polygons in geoframe B, treat the shape frame for B containing the polygons as a series of points and use geonear to find the point that is the nearest neighbor. Then you use frget and frlink to get the IDs of the nearest neighbor and the distance.

I imagine there is a more efficient way to do this in Stata, and almost certainly using dedicated GIS tools. However, I've put an excerpt of my code below in case it is helpful. I have a set of points that are reported to be in a village with id PC11_ID. However, geoframe spjoin has told me that these points are in fact not in that village. I use geonear to find the nearest point of the polygon defining the village with id PC11_ID. The main difference with your problem is that I am only looking for the nearest point of one polygon, while you will be looking for the nearest point of all tracts. I don't think that will be a problem conceptually but it might be very slow. It may actually make your code simpler if you don't need to do all the subsetting that I did to try to cut down on the sample size.

//# distance for unmatched observations 
frame create distances_unmatched 
levelsof pc11_id 
foreach PC11_ID in `r(levels)' {
  di `"** pc11_id "`PC11_ID'"  "'

  count if pc11_id=="`PC11_ID'" 
  geoframe select if pc11_id=="`PC11_ID'", into(temp_spot) 

  tempfile temp_shp

  frame village {
    geoframe select if pc11_id=="`PC11_ID'", ///
      into(temp_village)
 
  }
  frame temp_village_shp : save `temp_shp'

  frame temp_spot {
    keep response_id obs_num plot_latitude plot_longitude match_village_spjoin
    geonear obs_num plot_latitude plot_longitude ///
      using `temp_shp', neighbors(shape_order _Y _X) nearcount(1)
    replace km_to_nid = 0 if match_village_spjoin==1

    keep obs_num response_id km_to_nid
  }
  frame distances_unmatched : xframeappend temp_spot

  frame drop temp_spot temp_village temp_village_shp 
  rm `temp_shp'

}

frame distances_unmatched : desc, f

frlink 1:1 response_id , frame(distances_unmatched)
frget km_to_nid, from(distances_unmatched)