4

I have a polygon shapefile with about a dozen polygons (Area.shp). I have an equation that calculates a value for points (or small areas) as a function of distance to the edge of these polygons. I would ultimately like to apply this function to all points within each of these polygons and sum them (within each polygon) to get a single value for the polygon.

It seems like one of those things that shouldn't be too difficult, but I can't seem to get a working solution. I had an idea of converting the shapefile to raster then raster to points then taking the distance of each point to the nearest polygon edge, then apply the equation. It seemed like a reasonable, if clunk solution but unfortunately, I get all sorts of errors during various conversions using ArcMap 10.

Does anyone have any suggestions on how to do this? I am happy to use ArcGIS, QGIS, GRASS, Python, or export to R or any other solution. I am a novice in all of these systems (except R) but am willing to use whatever works.

[UPDATE] Additional Information: When I was thinking of points, it was really as a way to approximate numerical integration. I have an equation where the number of animals within a particular patch of habitat (polygon) is a function of distance to the edge of the habitat (nearest edge of the polygon). The equation is a 4pt logistic function 

N(i) = C + A/(1 + e^-(x-D)/B)

where all values are known (previously estimated) constants except x which is distance to the edge. N(i) is the number of animals in a small area, roughly 2 m^2. I'm interested in how the size and shape of habitat patches (polygons) affects the total N (sum of N(i) within a polygon) in a patch. I have a variety of patches identified from GIS layers as a single shape file with about a dozen polygons. For the purpose of this study, I can assume that no polygons overlap.

Improve this question
3
  • Summing values at "all points within" a polygon makes no sense because there are uncountably many such points. Perhaps you are seeking to integrate a function of distance over each polygon? Integrals are computed with focal means. As a practical matter, how you go about this computation depends on (a) whether any of the polygons overlap and (b) how the "equation" is calculated. Could you please supply those details?
    – whuber
    Jul 16, 2013 at 20:23
  • I added the information to the question. I'm sorry about being so vague. I haven't done much GIS work beyond creating maps, so I don't even really have the appropriate vocabulary or expertise to adequately specify the question. I'm happy to add more information as needed. Thanks for your patience.
    – djhocking
    Jul 16, 2013 at 20:54
  • 3
    You did great with the edit: the secret is to present the question in your terms rather than in technical or GIS terms. That at least guarantees everyone at your end of the conversation will understand what we're talking about :-) and the rest of us can ask for any clarification we need.
    – whuber
    Jul 16, 2013 at 21:19

1 Answer 1

Reset to default
6

Integration is done by summing the values and multiplying the common cell area (equal to the square of the cellsize).

Here is an example region to illustrate.

Region

  1. Create the Euclidean distance grid to the complement of the polygons. To do this, convert the polygons to raster format. This will place NoData values at all cells outside the polygons. Use IsNull and SetNull to place NoData only at cells inside the polygons. The compute the Euclidean distance grid of that. When hillshaded it should look something like this, with peaks inside each polygon and "ridges" extending along their "skeletons":

    Distance grid

  2. Use math operations, such as those offered in the Raster Calculator, to compute the values of N. In the formula, "x" is the Euclidean distance grid and A..D are constants. Here is a picture of the resulting grid, in pseudo 3D perspective:

    N grid

  3. Using the original polygon grid as the zone grid, compute the zonal sum of the result of (2). It will be a table with one number for each grid. Multiply those numbers by the square of the zone grid's cell size: those are the desired values.

Improve this answer
6
  • The logistic function is not likely to be a realistic model of population density. But you will find that out for yourself if you view the result after step 2. By the way, if N is given as number of individuals per two square meters and the cellsize is in square meters, then the result will be twice the number of individuals estimated to be in each region. Keeping careful track of your units of measurement will help you sort this out.
    – whuber
    Jul 16, 2013 at 22:13
  • Thanks I'll try this out at work tomorrow. Curious why you don't think the logistic function would be realistic? I get very few animals in the center of recent clearcuts then more near the edge, even more just into the forest but then it levels out. In reality I will have to do a more complicated version of this with buffers and inside and outside calculations to properly apply the logistic equation.
    – djhocking
    Jul 17, 2013 at 3:44
  • Can you further explain, "This will place NoData values at all cells outside the polygons. Use IsNull and SetNull to place NoData only at cells inside the polygons." If NoData is placed outside the polygons, then I use those commands to place NoData inside the polygons, does that mean that all cells in the whole layer have NoData values? Maybe I'm getting hung up on the "only" in the second sentence. What values go on the outside of the polygons?
    – djhocking
    Jul 17, 2013 at 17:10
  • Any valid value outside the polygon is ok. For more on how this works, please read the help page for Euclidean distance.
    – whuber
    Jul 17, 2013 at 18:52
  • Apparently I'm having trouble with the SetNull function. I converted polygon to raster. Resulting attribute table has values 0 - 8 (representing my 9 polygons) with a count for each. Then in the Python window I import arcpy and set the environment, then run outIsNull = isNull("AreaRas"), which results in a layer with OID (0,1), Value (0,1), and Count (34524, 1460408), where 0 is inside the polygons. Then I try outSetNull = SetNull("outIsNull", "outIsNull", "VALUE = 0"), but I end up with a table of [ID, Value, Count] = [1, 0, 162046] which plots as a strange streaky rectangle. Any ideas?
    – djhocking
    Jul 17, 2013 at 19:55

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.