Added dyn_stdev

gputils.py

@@ -1,3 +1,5 @@
+import math
+
 def dyn_mean(val, prev_mean, n):
     """Dynamic mean: computes the mean based on a previous mean plus a new value. Useful when mean is built
     incrementally, it saves the usage of huge arrays.
@@ -11,4 +13,21 @@ def dyn_mean(val, prev_mean, n):
         raise ValueError("n < 1, mean only defined for a positive number of elements")
     if n == 1:
         return val
-    return (prev_mean*(n-1)+val) / n
+    return (prev_mean*(n-1)+val) / n
+
+
+def dyn_stdev(val, prev_stdev, prev_mean, n):
+    """Dynamic stdev: computes the standard deviation based on a previous stdev plus a new value. Useful when stdev
+     is built incrementally, it saves the usage of huge arrays.
+
+    Keyword arguments:
+    val -- new val to add to the mean
+    prev_stdev -- previous stdev value
+    n -- number of total elements in the mean including the new val
+    """
+    if n < 1:
+        raise ValueError("n < 1, stdev only defined for a positive number of elements")
+    if n == 1:
+        return 0
+    curr_mean = dyn_mean(val, prev_mean, n)
+    return math.sqrt(((n-1)*prev_stdev*prev_stdev + (val - prev_mean)*(val - curr_mean)) / float(n))

tests/tests.py

@@ -1,6 +1,6 @@
-import unittest
+import unittest, time
 import numpy as np
-from gputils import dyn_mean
+from gputils import dyn_mean, dyn_stdev
 
 
 class TestMethods(unittest.TestCase):
@@ -20,6 +20,27 @@ def test_dyn_mean(self):
             test_mean = dyn_mean(values[-1], prev_mean, len(values))
             self.assertAlmostEqual(curr_mean, test_mean, 10)
 
+    def test_dyn_stdev(self):
+        # Special cases
+        with self.assertRaises(ValueError):
+            dyn_stdev(69, 0, 0, 0)
+        self.assertEqual(np.std(np.array(12)),dyn_stdev(12, 0, 0, 1))
+
+        # Normal cases
+        for i in range(10):
+            size = np.random.randint(10000, 100000)
+            values = np.random.rand(size) * 100
+            prev_mean = np.mean(values[:-1])
+            prev_std = np.std(values[:-1])
+            t0_trusted = time.time()
+            curr_std = np.std(values)
+            t1_trusted = time.time()
+            t0_test = time.time()
+            test_stdev = dyn_stdev(values[-1], prev_std, prev_mean, len(values))
+            t1_test = time.time()
+            self.assertAlmostEqual(curr_std, test_stdev, 10)
+            self.assertLessEqual(t1_test - t0_test, t1_trusted - t0_trusted) # ensuring we are faster
+
 
 if __name__ == '__main__':
     unittest.main()