4

I want to create a Stratified Random Sample with a land cover classification map. The purpose for the Stratrified Random Sample is to produce an Accuracy Assessment. I am using ArcGIS 10.2 with access to Spatial Analyst and all other tools.

I have six land cover types and I want 300 random points for three of the classes I am interested in.My first step was converting my raster to polygon. The result gives me over 3 million polygons to work with.

I have tried the Dissolve tool and making multi-part polygons, but some of my classes are too large, which results in an error message "Topology Invalid (Out of Memory).

I have tried using the Merge in the Editor toolbar. Again, due to my large number of polygons in some of my classes, it freezes and crashes ArcMap.

I would like to do the stratified random sample on a raster, rather than converting it to polygon form. Is this possible?

Are there any other tools/programs that I should use that may help combat my large dataset issue? I have consistently run into problems because I simply have too many polygons per land cover type.

3
  • Can you explain why you need the dissolve when you can summarize on an attribute. The only thing that I see the dissolve required for is to make the vector table look like the raster table (ie. zones to multiparts)
    – user681
    Commented May 14, 2015 at 18:29
  • @DanPatterson my thinking was that dissolve will create a multi-part polygon (combines all the individual polygons into one feature). This will allow me to use the Create Random Sample tool and create 300 points for the classes I am interested in. Otherwise, when I select all the individual polygons for a specific class, the Create Random Sample tool creates 300 points in each individual polygon. It at least tries to do that, but always ends with Error (Out of Memory) or crashes ArcMap. Commented May 14, 2015 at 18:39
  • strange...I always found that to be a problem converting raster to vector is that you end up with multipart polygons. Do you have a sample of your raster table? Did you produce the polygons based upon the class values? You could also remove unnecessary classes first before the conversion to polygon perhaps doing it class by class, then dissolving, then sampling, then merging them back together
    – user681
    Commented May 14, 2015 at 19:25

2 Answers 2

4

Stay in raster format.

If ArcGIS would provide a zonal percentile function you could do this in a couple easy steps. Since it doesn't, you have to work around that limitation. There are many ways to proceed. The key ideas are:

  1. Determine the proportion of cells in each class represented by the desired sample of 300.

  2. By thresholding a grid of uniform random values, you can identify slightly more than the number of points needed in each class.

  3. Extract these points, then post-process them by retaining only those having the lowest random values within each class.

This procedure obtains a simple random sample without replacement, and independently, of each land class.


Here are some details. Let there be N classes (n = 3 in the question), each with n(i) cells (i = 1, 2, ..., N). Suppose you intend to select k(i) randomly from each group (k(i) = 300 for the class indexes i = 1, 2, 3 in the question). You will do that by generating a uniform random grid and comparing it to a threshold grid [T]. It has the constant value p(i) at each cell in class i, with p(i) somewhere (to be chosen appropriately) between 0 and 1.

The result of this comparison in class i is a random count X(i) having a Binomial(n(i), p(i)) distribution. Its expectation therefore is m(i) = n(i)*p(i), which needs to be close to k(i). However, X(i) can randomly be less than that. Its variance is v(i) = n(i)*p(i)*(1-p(i)). To be really sure that X(i) will be at least k(i), we need a very small lower percentile of the distribution of X(i) to exceed k(i). You can work out what the resulting value of p(i) has to be by means of a statistical calculator, but an effective rule of thumb (for moderately small numbers of classes N) is to ensure that k(i) is less than several standard deviations less than m(i). Letting a reasonable starting value for "several" be L, we obtain the quadratic inequality

k(i) <= m(i) - L*sqrt(n(i)*p(i)*(1-p(i)))

whose approximate solution is

p(i) = (k(i)/n(i)) * (1 + z(i) + z(i)^2/2)

where

z(i) = L/sqrt(k(i)).

In the question, k(i) = 300 and it looks like the n(i) are in the millions. Taking L = 3, say, gives z(i) = 3 / sqrt(300), which is about 1/6 (and z(i)^2/2 is small enough to neglect) This tells us to ask for 1/6 more points than the desired proportion k(i)/n(i).

With this (simple) calculation out of the way,

  • Reclassify the land classes i into the values of the p(i) you computed. Call the resulting grid [P].

  • Create a uniform random grid [R].

  • Set the values of [R] to null wherever [R] > [P].

This last grid will have slightly more than X = k(1) + k(2) + ... + k(N) (= 900 in the question) cells that are not null. At this juncture use your favorite method (there are many that will work) to

  1. combine it with the land class grid and

  2. output all tuples (i, R(i), x coordinate, y coordinate) as a point shapefile.

Post-process the shapefile by retaining only the lowest k(i) values for each class i: that's the desired sample. If you are so extremely unlucky as to not obtain enough points in some class, just repeat the whole process from the beginning. (Consider increasing the value of L when you do so.)


Alternatively, doing the work in R would be simple, because it offers functions to find the lowest k(i) random values in each class directly, rather than having to take several steps and convert things into a point shapefile format. The same conceptual workflow is involved.

1

It seems like all you really need is a random selection of cells (for the raster) or features (for vector), for each class. When I had to do this I broke my classification into separate rasters for each class, then converted to points, then used a random number generator and a definition query to keep only those features with matching FIDs. You could also use something like Hawth's Tools to create a random selection on your point file. In either case just export the selection as your accuracy assessment locations. I did this exact thing for a very large land cover classification for my thesis and it worked well.

3
  • Such methods often do work, although they can be inefficient due to the need to convert from grid to vector formats. A truly "very large land cover classification" would likely produce unwieldy datasets; for instance, a giga-cell grid would expand from perhaps a few hundred MB on disk to a vector dataset of at least 24 GB (and likely larger), so I am curious how large your dataset was. Did you notice that the question asks explicitly for a solution that works directly on a raster?
    – whuber
    Commented Dec 19, 2015 at 4:05
  • I did see that he asked for a solution that works directly on a raster. I just figured this solution would work too, since it worked for me. My land cover classification covered several thousand square miles, and yea, my data was unwieldy - I was dealing with millions of points. It wasn't quick, but it did work. I'm not sure a very efficient way exists to do this in ArcGis, but I'm guessing the conversion process could be automated.
    – CSB
    Commented Dec 19, 2015 at 4:32
  • The procedure I described works in ArcGIS and is efficient. The main drawback (for non-students) is that it requires the Spatial Analyst extension. It's a powerful tool to have and to know how to use, because as a general principle calculations on raster data proceed better when carried out directly on those data instead of converting them to a vector format.
    – whuber
    Commented Dec 19, 2015 at 4:53

Your Answer

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

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