Add CircularMean

Author: rikonor

stat.go

@@ -125,6 +125,31 @@ func ChiSquare(obs, exp []float64) float64 {
 	return result
 }
 
+// CircularMean returns the circular mean of the dataset.
+//	atan2(\sum_i w_i * sin(alpha_i), \sum_i w_i * cos(alpha_i))
+// If weights is nil then all of the weights are 1. If weights is not nil, then
+// len(x) must equal len(weights).
+func CircularMean(x, weights []float64) float64 {
+	if weights != nil && len(x) != len(weights) {
+		panic("stat: slice length mismatch")
+	}
+
+	var aX, aY float64
+	if weights != nil {
+		for i, v := range x {
+			aX += weights[i] * math.Cos(v)
+			aY += weights[i] * math.Sin(v)
+		}
+	} else {
+		for _, v := range x {
+			aX += math.Cos(v)
+			aY += math.Sin(v)
+		}
+	}
+
+	return math.Atan2(aY, aX)
+}
+
 // Correlation returns the weighted correlation between the samples of x and y
 // with the given means.
 //  sum_i {w_i (x_i - meanX) * (y_i - meanY)} / (stdX * stdY)

stat_test.go

@@ -13,6 +13,52 @@ import (
 	"github.com/gonum/floats"
 )
 
+func ExampleCircularMean() {
+	x := []float64{0, 0.25 * math.Pi, 0.75 * math.Pi}
+	weights := []float64{1, 2, 2.5}
+	cmean := CircularMean(x, weights)
+
+	fmt.Printf("The circular mean is %.5f.\n", cmean)
+	// Output:
+	// The circular mean is 1.37037.
+}
+
+func TestCircularMean(t *testing.T) {
+	for i, test := range []struct {
+		x   []float64
+		wts []float64
+		ans float64
+	}{
+		// Values compared against scipy.
+		{
+			x:   []float64{0, 2 * math.Pi},
+			ans: 0,
+		},
+		{
+			x:   []float64{0, 0.5 * math.Pi},
+			ans: 0.78539816339744,
+		},
+		{
+			x:   []float64{-1.5 * math.Pi, 0.5 * math.Pi, 2.5 * math.Pi},
+			wts: []float64{1, 2, 3},
+			ans: 0.5 * math.Pi,
+		},
+		{
+			x:   []float64{0, 0.5 * math.Pi},
+			wts: []float64{1, 2},
+			ans: 1.10714871779409,
+		},
+	} {
+		c := CircularMean(test.x, test.wts)
+		if math.Abs(c-test.ans) > 1e-14 {
+			t.Errorf("Circular mean mismatch case %d: Expected %v, Found %v", i, test.ans, c)
+		}
+	}
+	if !Panics(func() { CircularMean(make([]float64, 3), make([]float64, 2)) }) {
+		t.Errorf("CircularMean did not panic with x, wts length mismatch")
+	}
+}
+
 func ExampleCorrelation() {
 	x := []float64{8, -3, 7, 8, -4}
 	y := []float64{10, 5, 6, 3, -1}