# How to calculate mean aspect in circular statistics?

Measurements of angles (e.g. compass bearings, aspect) are circular in nature (e.g. 0 and 360 degrees coincide) and therefore cannot be meaningfully described by ordinary arithmetic means. I have a shapefile of 22 locations containing aspect value for each point. I want to calculate mean aspect of the locations using circular statistics in R. "circstats" is one of the packages that are suggested to me. Question: Do I need to transform aspect data before feeding into R (aspect data is calculated in GRASS)? How I can calculate the mean direction of aspect in R?

In R, the package CircStats is old and of rather limited scope and has been replaced by the more complete Circular package. There are tutorials and a book, Circular Statistics with R (2013, A. Pewsey, M. Neuhäuser, and G. D. Ruxton, Oxford University Press, 208 pp.) which explains how to use it (The R scripts can be downloaded from the resources site of the book)

1) If the data are directional (= circular with values in the entire 0–360◦ range and 0-360° equivalent), the solution is to use the angular values as vectors an apply the basic vector operations. The mean direction is then found by computing the vector sum of the vectors that represent the various directions in the data.

Since the coordinates of these vectors (normalized) are `Xi =cosθi` and `Yi =sinθi`, the mean direction is given by `atan2(sum(sinθi),sum(cosθi))` ( -> the solution of RadouxJu)

2) If the data are axial (two-headed vectors with a 180◦ ambiguity, as a line) there are others solutions.

Therefore for circular data, you have many solutions with GRASS GIS:

• run R in the GRASS shell and open the GRASS vector with the spgrass6 (GRASS 6.4.x) or rgrass7 (GRASS 7.x) packages (look at spgrass6 and GIS_LOCK, for vectors, the command is `readVECT6`) and use the Circular package

``````library(circular)
angles = c(341.0, 359.0, 334.0, 15.0, 330.0, 301.0, 299.0, 9.0, 7.0, 353.0, 28.0, 25.0, 23.0, 30.0, 350.0, 25.0, 22.0, 8.0, 356.0, 27.0)
anglecir =  circular(angles, type="angles", units="degrees",modulo="2pi", template='geographics')
summary(anglecir)
n     Min.  1st Qu.   Median     Mean  3rd Qu.     Max.      Rho
20.0000  30.0000  23.5000   7.5000   1.0610 -12.2500 -61.0000   0.8962
mean(anglecir)
Circular Data:
Type = angles
Units = degrees
Template = geographics
Modulo = 2pi
Zero = 1.570796
Rotation = clock
[1] 1.060902
``````
• run Python from the GRASS shell (with `grass.script` or `pygrass` modules) or from the Python shell and use the solution of radouxju, very easy to implement in Python.

``````import numpy as np
angles = np.array([341.0, 359.0, 334.0, 15.0, 330.0, 301.0, 299.0, 9.0, 7.0, 353.0, 28.0, 25.0, 23.0, 30.0, 350.0, 25.0, 22.0, 8.0, 356.0, 27.0])
def getCircularMean(angles):
n = len(angles)
vectorMean = np.arctan2(sineMean, cosineMean)
return np.degrees(vectorMean)
print "{:.6f}".format(getCircularMean(angles))
1.060902
``````
• Good references. However, you ought to update your `arctan` solution, which is invalid (as well as being syntactically incorrect in `R`), with a correct one based on `atan2`. May 13, 2015 at 17:34
• done, thanks, but it was pure formalism here because the real trigonometric function is `Arctan`with `atan2` instead of `atan` for computations and correct results.
– gene
May 13, 2015 at 18:03
• gene... and of course this works in ArcMap's python shell as well, for completeness
– user681
May 14, 2015 at 18:23
• ArcGIS has also the command Linear Directional Mean (Spatial Statistics), but for lines only (not points with attributes values)
– gene
May 14, 2015 at 18:55

you can convert your aspects into the sine and cosine, compute the mean of the sine's and the mean of the cosine's, then turn it back to aspect using atan2(sine,cosine). For more details, see Wikipedia

• Do I need also to transform aspect degrees to radius before changing to sin/cos?
– Mah
May 13, 2015 at 14:20
• by using atan(sine/cosine), I have a negative number above 1. Is it suppose to be between -1 and 1?
– Mah
May 13, 2015 at 15:32
• An arctangent is an angle--it can be anything. However, you shouldn't be using `atan` at all--use the `atan2` version (called `atan(x,y)` in the GRASS manual; its value is in degrees). Otherwise you will only get a 180 degree range of angles rather than the full 360 degrees. May 13, 2015 at 15:42
• @whuber Thank you, with atan2(x,y) you don't need to test the sign of sine and cos. I've updated my answer. May 13, 2015 at 17:30
• @Mah yes, R works with angles in radian, so you need to transform from degree to radian before sin and cos, then from radian to degree after atan2 May 13, 2015 at 17:31

Most textbooks suggest using atan2(Sigma(sin(x)), Sigma(cos(x))), however this is not always the right thing to do. For example, the average of 0, 0 and 90 degrees is atan( (sin(0)+sin(0)+sin(90)) / (cos(0)+cos(0)+cos(90)) ) = atan(1/2)= 26.56 deg, and not 30 deg as one may expect.

Take a look at my article on CodeProject "Circular Values Math and Statistics with C++11", especially section 26 - Averaging n circular values, and section 29 - Circular parameter estimation based on noisy measurements, were I suggest a method that is consistent with the general definition on an average, and also more accurate when the data has a wrapped normal distribution.

An equivalent method is suggested in Edwin Olson's paper "On computing the average orientation of vectors and lines" (Algorithm 1).

• This is an interesting and valuable observation (+1). Ultimately, the choice of how to average directional values depends on the purpose of the analysis as well as the nature of the random errors (the "noise"). Aspects are computed from DEMs, typically, and as such tend to be functions of least squares estimates of locally fitted surfaces. Although the distributions of those parameter estimates tend to be approximately Normal, the aspects (as angles) will not be. In effect, the atan2(sin,cos) procedure appears to be a more appropriate mean than Olson's method. May 14, 2015 at 19:26