diff --git a/stat.go b/stat.go
index ee2ad0b6..387d29c8 100644
--- a/stat.go
+++ b/stat.go
@@ -236,10 +236,10 @@ func Covariance(x, y, weights []float64) float64 {
 		w := weights[i]
 		yv := y[i]
 		wxd := w * (xv - xu)
-		wyd := w * (yv - yu)
-		ss += wxd * wyd
+		yd := (yv - yu)
+		ss += wxd * yd
 		xcompensation += wxd
-		ycompensation += wyd
+		ycompensation += w * yd
 		sumWeights += w
 	}
 	// xcompensation and ycompensation are from Chan, et. al.
diff --git a/stat_test.go b/stat_test.go
index a0ee39fa..51e9a54c 100644
--- a/stat_test.go
+++ b/stat_test.go
@@ -128,6 +128,12 @@ func TestCovariance(t *testing.T) {
 			weights: []float64{1, 1.5, 1},
 			ans:     3.2,
 		},
+		{
+			p:       []float64{1, 4, 9},
+			q:       []float64{1, 4, 9},
+			weights: []float64{1, 1.5, 1},
+			ans:     13.142857142857146,
+		},
 	} {
 		c := Covariance(test.p, test.q, test.weights)
 		if math.Abs(c-test.ans) > 1e-14 {