This posts looks at being mindful of creating ambiguous results in your spatial analysis, in a raster context

The GIS Team recently did some work with Sonya Olykan and Derrick Moot on a North Island dryland farm.  The aim was to identify where subterranean clover might be best  grown on the Tokaroa Farm, near Martinborough.  It’s a 608 hectare farm running sheep and beef (and do I sound like I know what I’m talking about here?).

The aim of this analysis was to find the best places to grow subterranean clover on the farm.  The criteria were pretty straightforward and only considered slope and aspect over two parts of the farm: the flats and the ranges.

For slope there were three classes:

  • Flat: (0 – 7º)
  • Drillable: (8 – 20º)
  • Steep: (> 20º)

Flat and drillable are the most preferred areas for slope.  The slope layer was reclassified so that:

  • Flat = 1
  • Drillable = 2
  • Steep = 3

(NB: there is a pretty big issue hiding in these classes, but I’ll overlook it for now.)

For Aspect, two classes with their reclassified values:

  • Sunny (northerly):  W + NW + N + NE = 1
  • Shady (southerly):  E + SE + S + SW = 0

As you might imagine, the north facing slopes are preferable.

I now need to find a way to combine those two layers.

So not too challenging to get these layers – slope and aspect are easy enough to derive once you’ve got a DEM, and then I can reclassify those layers to fit these criteria.

The bigger challenge (and the reason why this is a blog post) is the next step of bringing these layers together in a meaningful way.  Part of that means not creating a ambiguous results, or results that are difficult to interpret.  GIS can be very much a GIGO kind of thing – good analysts will think ahead to any problems that might arise in a workflow and do their best to sidestep it ahead of time.  Let’s start with our layers to illustrate this.  Here’s my reclassified slope layer using the classes from above:

And here’s my aspect, simpler with just 0 (sunny) and 1 (shady).

I next want to combine these layers so I can clearly know where the different conditions are, e.g flat/sunnyor drillable/shady.  The key to do that is to think about how to combine these layers in a raster calculation to avoid anything ambiguous.  As it is I’m combining 0s and 1s with 1s, 2s and 3s.  I’ve got two simple options here, for the sake of argument: I could either add the layers together or multiply them.  I’ll use a table here to illustrate the possible outputs:

S + A 1 2 3
0 1 2 3
1 2 3 4

If I do this, I’ll get some ambiguous results – I won’t be able to separate the flat/sunny slopes from the drillable/shady ones.  Same for drillable/sunny and steep/sunny.

What if I multiple them together?

S * A 1 2 3
0 0 0 0
1 1 2 3

Clearly that’s not going to get me any closer.  So I’ll do a bit of hybrid.  I’ll multiply my slope classes by 10 and then add them to my aspect classes.  This is just a raster calculation which I’ll do in one fell swoop:

What do I get from this?

10*S + A 10 20 30
0 10 20 30
1 11 21 31

That’s cleared up any ambiguity for me now.  I can clearly separate out the classes and know what’s what with certainty.

  • 10 = shady and flat
  • 11 = sunny and flat
  • 20 = shady and drillable
  • 21 = sunny and drillable
  • 30 = shady and steep
  • 31 = sunny and steep

Here’s a subtle point that goes along with what were doing here.  I started with Slope and Aspect grids – the values in those grids had a physical meaning.  As soon as I reclassified them into their different classes, the numbers ceased being used as numbers.  Now they’re being used like names; their values signify classes, or better still, categories.  I have to do this with numbers because when we do raster calculations, it’s only a grid’s values that is worked with – none of the other attributes can really be used.  Here I’m using numbers as names (in the parlance of the data analysis trade, they are now nominal data).  A value of 21 in my table above has no physical meaning whatsoever, but it points to a category of cells that are drillable/northerly slopes.

For the final mapping, I need to be able to separate these out as separate areas on the map, so removing any ambiguity is pretty critical.  Here’s a final map of these classes, highlighting the most and least preferable areas:

(Hmmmm…really debating if that hillshade underneath is helping or hurting.  What do you think?)

Behind the scenes you can see how I’ve grouped the classes together for the symbology – note the raster values that tell me what’s what:

Of course, you could easily make the argument that I should have thought this through when I assigned the new values in my reclassified grids – and right you would be – but then what would I have to blather on about?