Merge pull request #961 from Madhav2310/stats4

Covariance and Correlation functions implemented

Author: certik

integration_tests/test_statistics.py

@@ -1,5 +1,5 @@
 from statistics import (mean, fmean, geometric_mean, harmonic_mean,
-                        variance, stdev)
+                        variance, stdev, covariance, correlation)
 from ltypes import i32, f64, i64
 
 eps: f64
@@ -66,6 +66,46 @@ def test_variance():
     k = variance(b)
     assert abs(k - 0.40924) < eps
 
+def test_covariance():
+    a: list[i32]
+    a = [1, 2, 3, 4, 5, 6, 7, 8, 9]
+    b: list[i32]
+    b = [1, 2, 3, 1, 2, 3, 1, 2, 3]
+    j: f64
+    j = covariance(a,b)
+    assert abs(j - 0.75) < eps
+
+    c: list[f64]
+    c = [2.74, 1.23, 2.63, 2.22, 3.0, 1.98]
+    d: list[f64]
+    d = [9.4, 1.23, 2.63, 22.4, 1.9, 13.98]
+    k: f64
+    k = covariance(c,d)
+    assert abs(k + 0.24955999999999934) < eps
+
+def test_correlation():
+    a: list[i32]
+    a = [11, 2, 7, 4, 15, 6, 10, 8, 9, 1, 11, 5, 13, 6, 15]
+    b: list[i32]
+    b = [2, 5, 17, 6, 10, 8, 13, 4, 6, 9, 11, 2, 5, 4, 7]
+
+    j: f64
+    j = correlation(a,b)
+    assert abs(j - 0.11521487988958108) < eps
+
+    c: list[i32]
+    c = [1, 2, 3, 4, 5, 6, 7, 8, 9]
+    d: list[i32]
+    d = [9, 8, 7, 6, 5, 4, 3, 2, 1]
+
+    k: f64
+    k = correlation(c,c)
+    assert k == 1.0
+
+    l: f64
+    l = correlation(c,d)
+    assert l == -1.0
+
 
 def test_stdev():
     a: list[i32]
@@ -88,5 +128,7 @@ def check():
     test_fmean()
     test_variance()
     test_stdev()
+    test_covariance()
+    test_correlation()
 
 check()

src/runtime/statistics.py

@@ -1,5 +1,5 @@
 from ltypes import i32, f64, overload
-
+from math import sqrt
 
 
 @overload
@@ -117,7 +117,7 @@ def geometric_mean(x: list[i32]) -> f64:
     for i in range(k):
         product *= float(x[i])
 
-    return product**(1/k)
+    return product ** (1 / k)
 
 
 def harmonic_mean(x: list[i32]) -> f64:
@@ -136,7 +136,7 @@ def harmonic_mean(x: list[i32]) -> f64:
             return 0.0
         sum += 1 / x[i]
 
-    return k/sum
+    return k / sum
 
 @overload
 def variance(x: list[f64]) -> f64:
@@ -153,8 +153,8 @@ def variance(x: list[f64]) -> f64:
     num = 0.0
     i: i32
     for i in range(n):
-        num += (x[i]-xmean)**2
-    return num/(n-1)
+        num += (x[i] - xmean)**2
+    return num / (n-1)
 
 @overload
 def variance(x: list[i32]) -> f64:
@@ -171,8 +171,8 @@ def variance(x: list[i32]) -> f64:
     num = 0.0
     i: i32
     for i in range(n):
-        num += (x[i]-xmean)**2
-    return num/(n-1)
+        num += (x[i] - xmean)**2
+    return num / (n-1)
 
 
 @overload
@@ -188,3 +188,104 @@ def stdev(x: list[i32]) -> f64:
     Returns the standard deviation of a data sequence of numbers
     """
     return variance(x)**0.5
+
+@overload
+def covariance(x: list[i32], y: list[i32]) -> f64:
+    """
+    Returns the covariance of a data sequence of numbers
+    """
+    n: i32 = len(x)
+    m: i32 = len(y)
+    if (n < 2 or m < 2) or n != m:
+        raise Exception("Both inputs must be of the same length (no less than two)")
+    xmean: f64 = mean(x)
+    ymean: f64 = mean(y)
+    num: f64
+    num = 0.0
+    i: i32
+    for i in range(n):
+        num += (x[i] - xmean) * (y[i] - ymean)
+    return num / (n-1)
+
+@overload
+def covariance(x: list[f64], y: list[f64]) -> f64:
+    """
+    Returns the covariance of a data sequence of numbers
+    """
+    n: i32 = len(x)
+    m: i32 = len(y)
+    if (n < 2 or m < 2) or n != m:
+        raise Exception("Both inputs must be of the same length (no less than two)")
+    xmean: f64 = mean(x)
+    ymean: f64 = mean(y)
+    num: f64
+    num = 0.0
+    i: i32
+    for i in range(n):
+        num += (x[i] - xmean) * (y[i] - ymean)
+    return num / (n-1)
+
+@overload
+def correlation(x: list[i32], y: list[i32]) -> f64:
+    """
+    Return the Pearson's correlation coefficient for two inputs.
+    """
+    n: i32 = len(x)
+    m: i32 = len(y)
+    if n != m:
+        raise Exception("correlation requires that both inputs have same number of data points")
+    if n < 2:
+        raise Exception("correlation requires at least two data points")
+    xmean: f64 = mean(x)
+    ymean: f64 = mean(y)
+
+    sxy: f64 = 0.0
+    i: i32
+    for i in range(n):
+        sxy += (x[i] - xmean) * (y[i] - ymean)
+
+    sxx: f64 = 0.0
+    j: i32
+    for j in range(n):
+        sxx += (x[j] - xmean) ** 2
+
+    syy: f64 = 0.0
+    k: i32
+    for k in range(n):
+        syy += (y[k] - ymean) ** 2
+    if sqrt(sxx * syy) == 0:
+        raise Exception('at least one of the inputs is constant')
+    return sxy / sqrt(sxx * syy)
+
+@overload
+def correlation(x: list[f64], y: list[f64]) -> f64:
+    """
+    Return the Pearson's correlation coefficient for two inputs.
+    """
+    n: i32 = len(x)
+    m: i32 = len(y)
+    if n != m:
+        raise Exception("correlation requires that both inputs have same number of data points")
+    if n < 2:
+        raise Exception("correlation requires at least two data points")
+    xmean: f64 = mean(x)
+    ymean: f64 = mean(y)
+
+    sxy: f64 = 0.0
+    i: i32
+    for i in range(n):
+        sxy += (x[i] - xmean) * (y[i] - ymean)
+
+    sxx: f64 = 0.0
+    j: i32
+    for j in range(n):
+        sxx += (x[j] - xmean) ** 2
+
+    syy: f64 = 0.0
+    k: i32
+    for k in range(n):
+        syy += (y[k] - ymean) ** 2
+    if sqrt(sxx * syy) == 0:
+        raise Exception('at least one of the inputs is constant')
+    return sxy / sqrt(sxx * syy)
+