Add axis arguments to various functions so as not to rely on the defaults.

16 files changed, 72 insertions, 83 deletions

diff --git a/numpy/add_newdocs.py b/numpy/add_newdocs.py
index ea8915af3..c91e8482c 100644
--- a/numpy/add_newdocs.py
+++ b/numpy/add_newdocs.py
@@ -657,7 +657,7 @@ add_newdoc('numpy.core.multiarray', 'ndarray', ('mean',
   is given, this equals:
       a.sum(axis, dtype) * 1.0 / len(a).
   If axis is None, this equals:
-      a.sum(axis, dtype) * 1.0 / product(a.shape)
+      a.sum(axis, dtype) * 1.0 / product(a.shape,axis=0)
   The optional dtype argument is the data type for intermediate
   calculations in the sum.;
 
@@ -834,7 +834,7 @@ distribution.
 
 The standard deviation is the square root of the average of the
 squared deviations from the mean, i.e.
-    std = sqrt(mean((x - x.mean())**2)).
+    std = sqrt(mean((x - x.mean())**2,axis=0)).
 
 For multidimensional arrays, std is computed by default along the
 first axis.
@@ -853,7 +853,7 @@ value and intermediate calculations.  The default is to upcast
 (promote) smaller integer types to the platform-dependent int.
 For example, on 32-bit platforms:
 
-  a.dtype                         default sum() dtype
+  a.dtype                         default sum dtype
   ---------------------------------------------------
   bool, int8, int16, int32        int32
 
diff --git a/numpy/core/ma.py b/numpy/core/ma.py
index 2e2d3f8ec..67335019f 100644
--- a/numpy/core/ma.py
+++ b/numpy/core/ma.py
@@ -1613,7 +1613,7 @@ def average (a, axis=None, weights=None, returned = 0):
        If axis is None, average over the entire array
        Inputs can be integer or floating types; result is of type float.
 
-       If weights are given, result is sum(a*weights)/(sum(weights)*1.0)
+       If weights are given, result is sum(a*weights,axis=0)/(sum(weights,axis=0)*1.0)
        weights must have a's shape or be the 1-d with length the size
        of a in the given axis.
 
diff --git a/numpy/core/tests/test_ma.py b/numpy/core/tests/test_ma.py
index 3607b965b..d7380a127 100644
--- a/numpy/core/tests/test_ma.py
+++ b/numpy/core/tests/test_ma.py
@@ -167,14 +167,14 @@ class test_ma(NumpyTestCase):
         (x, y, a10, m1, m2, xm, ym, z, zm, xf, s) = self.d
         self.failUnless (eq(numpy.add.reduce(x), add.reduce(x)))
         self.failUnless (eq(numpy.add.accumulate(x), add.accumulate(x)))
-        self.failUnless (eq(4, sum(array(4))))
+        self.failUnless (eq(4, sum(array(4),axis=0)))
         self.failUnless (eq(4, sum(array(4), axis=0)))
-        self.failUnless (eq(numpy.sum(x), sum(x)))
-        self.failUnless (eq(numpy.sum(filled(xm,0)), sum(xm)))
+        self.failUnless (eq(numpy.sum(x,axis=0), sum(x,axis=0)))
+        self.failUnless (eq(numpy.sum(filled(xm,0),axis=0), sum(xm,axis=0)))
         self.failUnless (eq(numpy.sum(x,0), sum(x,0)))
-        self.failUnless (eq(numpy.product(x), product(x)))
+        self.failUnless (eq(numpy.product(x,axis=0), product(x,axis=0)))
         self.failUnless (eq(numpy.product(x,0), product(x,0)))
-        self.failUnless (eq(numpy.product(filled(xm,1)), product(xm)))
+        self.failUnless (eq(numpy.product(filled(xm,1),axis=0), product(xm,axis=0)))
         if len(s) > 1:
             self.failUnless (eq(numpy.concatenate((x,y),1), concatenate((xm,ym),1)))
             self.failUnless (eq(numpy.add.reduce(x,1), add.reduce(x,1)))
@@ -308,9 +308,9 @@ class test_ma(NumpyTestCase):
         m = [1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 1]
         i = numpy.nonzero(m)[0]
         putmask(xm, m, z)
-        assert take(xm, i) == z
+        assert take(xm, i,axis=0) == z
         put(ym, i, zm)
-        assert take(ym, i) == zm
+        assert take(ym, i,axis=0) == zm
 
     def check_testOddFeatures(self):
         "Test of other odd features"
@@ -513,13 +513,13 @@ class test_ma(NumpyTestCase):
     def check_testAverage1(self):
         "Test of average."
         ott = array([0.,1.,2.,3.], mask=[1,0,0,0])
-        self.failUnless(eq(2.0, average(ott)))
+        self.failUnless(eq(2.0, average(ott,axis=0)))
         self.failUnless(eq(2.0, average(ott, weights=[1., 1., 2., 1.])))
         result, wts = average(ott, weights=[1.,1.,2.,1.], returned=1)
         self.failUnless(eq(2.0, result))
         self.failUnless(wts == 4.0)
         ott[:] = masked
-        self.failUnless(average(ott) is masked)
+        self.failUnless(average(ott,axis=0) is masked)
         ott = array([0.,1.,2.,3.], mask=[1,0,0,0])
         ott=ott.reshape(2,2)
         ott[:,1] = masked
@@ -539,20 +539,20 @@ class test_ma(NumpyTestCase):
         y=array([arange(6), 2.0*arange(6)])
         self.failUnless(allclose(average(y, None), numpy.add.reduce(numpy.arange(6))*3./12.))
         self.failUnless(allclose(average(y, axis=0), numpy.arange(6) * 3./2.))
-        self.failUnless(allclose(average(y, axis=1), [average(x), average(x) * 2.0]))
+        self.failUnless(allclose(average(y, axis=1), [average(x,axis=0), average(x,axis=0) * 2.0]))
         self.failUnless(allclose(average(y, None, weights=w2), 20./6.))
         self.failUnless(allclose(average(y, axis=0, weights=w2), [0.,1.,2.,3.,4.,10.]))
-        self.failUnless(allclose(average(y, axis=1), [average(x), average(x) * 2.0]))
+        self.failUnless(allclose(average(y, axis=1), [average(x,axis=0), average(x,axis=0) * 2.0]))
         m1 = zeros(6)
         m2 = [0,0,1,1,0,0]
         m3 = [[0,0,1,1,0,0],[0,1,1,1,1,0]]
         m4 = ones(6)
         m5 = [0, 1, 1, 1, 1, 1]
-        self.failUnless(allclose(average(masked_array(x, m1)), 2.5))
-        self.failUnless(allclose(average(masked_array(x, m2)), 2.5))
-        self.failUnless(average(masked_array(x, m4)) is masked)
-        self.assertEqual(average(masked_array(x, m5)), 0.0)
-        self.assertEqual(count(average(masked_array(x, m4))), 0)
+        self.failUnless(allclose(average(masked_array(x, m1),axis=0), 2.5))
+        self.failUnless(allclose(average(masked_array(x, m2),axis=0), 2.5))
+        self.failUnless(average(masked_array(x, m4),axis=0) is masked)
+        self.assertEqual(average(masked_array(x, m5),axis=0), 0.0)
+        self.assertEqual(count(average(masked_array(x, m4),axis=0)), 0)
         z = masked_array(y, m3)
         self.failUnless(allclose(average(z, None), 20./6.))
         self.failUnless(allclose(average(z, axis=0), [0.,1.,99.,99.,4.0, 7.5]))
@@ -678,10 +678,10 @@ class test_ufuncs(NumpyTestCase):
 
     def test_reduce(self):
         a = self.d[0]
-        self.failIf(alltrue(a))
-        self.failUnless(sometrue(a))
-        self.failUnlessEqual(sum(a[:3]), 0)
-        self.failUnlessEqual(product(a), 0)
+        self.failIf(alltrue(a,axis=0))
+        self.failUnless(sometrue(a,axis=0))
+        self.failUnlessEqual(sum(a[:3],axis=0), 0)
+        self.failUnlessEqual(product(a,axis=0), 0)
 
     def test_minmax(self):
         a = arange(1,13).reshape(3,4)
@@ -738,7 +738,7 @@ class test_array_methods(NumpyTestCase):
         mXdiag = mX.diagonal()
         self.assertEqual(mX.trace(), mX.diagonal().compressed().sum())
         self.failUnless(eq(mX.trace(),
-                           X.trace() - sum(mXdiag.mask*X.diagonal())))
+                           X.trace() - sum(mXdiag.mask*X.diagonal(),axis=0)))
 
     def test_clip(self):
         (x,X,XX,m,mx,mX,mXX,) = self.d
diff --git a/numpy/core/tests/test_numeric.py b/numpy/core/tests/test_numeric.py
index f6354b70b..44a69f8e8 100644
--- a/numpy/core/tests/test_numeric.py
+++ b/numpy/core/tests/test_numeric.py
@@ -228,8 +228,8 @@ class test_fromiter(NumpyTestCase):
         expected = array(list(self.makegen()))
         a = fromiter(self.makegen(), int)
         a20 = fromiter(self.makegen(), int, 20)
-        self.failUnless(alltrue(a == expected))
-        self.failUnless(alltrue(a20 == expected[:20]))
+        self.failUnless(alltrue(a == expected,axis=0))
+        self.failUnless(alltrue(a20 == expected[:20],axis=0))
 
 class test_index(NumpyTestCase):
     def test_boolean(self):
diff --git a/numpy/f2py/tests/array_from_pyobj/tests/test_array_from_pyobj.py b/numpy/f2py/tests/array_from_pyobj/tests/test_array_from_pyobj.py
index e1d4a47a6..ab8b62f4b 100644
--- a/numpy/f2py/tests/array_from_pyobj/tests/test_array_from_pyobj.py
+++ b/numpy/f2py/tests/array_from_pyobj/tests/test_array_from_pyobj.py
@@ -222,7 +222,7 @@ class Array:
         if arr1.shape != arr2.shape:
             return False
         s = arr1==arr2
-        return alltrue(s.flatten())
+        return alltrue(s.flatten(),axis=0)
 
     def __str__(self):
         return str(self.arr)
diff --git a/numpy/fft/fftpack.py b/numpy/fft/fftpack.py
index 1cb24f2b7..ffa6ac18d 100644
--- a/numpy/fft/fftpack.py
+++ b/numpy/fft/fftpack.py
@@ -200,7 +200,7 @@ def _cook_nd_args(a, s=None, axes=None, invreal=0):
         if axes == None:
             s = list(a.shape)
         else:
-            s = take(a.shape, axes)
+            s = take(a.shape, axes,axis=0)
     else:
         shapeless = 0
     s = list(s)
diff --git a/numpy/lib/arraysetops.py b/numpy/lib/arraysetops.py
index 7bd666029..b98517f3d 100644
--- a/numpy/lib/arraysetops.py
+++ b/numpy/lib/arraysetops.py
@@ -179,7 +179,7 @@ def test_unique1d_speed( plot_results = False ):
         dt1s.append( dt1 )
         dt2s.append( dt2 )
 
-        assert numpy.alltrue( b == c )
+        assert numpy.alltrue( b == c)
 
 
     print nItems
diff --git a/numpy/lib/function_base.py b/numpy/lib/function_base.py
index d8cd30f2c..a202f67cb 100644
--- a/numpy/lib/function_base.py
+++ b/numpy/lib/function_base.py
@@ -108,17 +108,16 @@ def average(a, axis=None, weights=None, returned=False):
     """average(a, axis=None weights=None, returned=False)
 
     Average the array over the given axis.  If the axis is None, average
-    over all dimensions of the array.  Equivalent to a.mean(axis), but
-    with a default axis of 0 instead of None.
+    over all dimensions of the array.  Equivalent to a.mean(axis)
 
     If an integer axis is given, this equals:
         a.sum(axis) * 1.0 / len(a)
 
     If axis is None, this equals:
-        a.sum(axis) * 1.0 / product(a.shape)
+        a.sum(axis) * 1.0 / product(a.shape,axis=0)
 
     If weights are given, result is:
-        sum(a * weights) / sum(weights),
+        sum(a * weights,axis) / sum(weights,axis),
     where the weights must have a's shape or be 1D with length the
     size of a in the given axis. Integer weights are converted to
     Float.  Not specifying weights is equivalent to specifying
@@ -541,9 +540,9 @@ def extract(condition, arr):
     """Return the elements of ravel(arr) where ravel(condition) is True
     (in 1D).
 
-    Equivalent to compress(ravel(condition), ravel(arr)).
+    Equivalent to compress(ravel(condition), ravel(arr),0).
     """
-    return _nx.take(ravel(arr), nonzero(ravel(condition))[0])
+    return _nx.take(ravel(arr), nonzero(ravel(condition))[0],axis=0)
 
 def place(arr, mask, vals):
     """Similar to putmask arr[mask] = vals but the 1D array vals has the
diff --git a/numpy/lib/shape_base.py b/numpy/lib/shape_base.py
index df0e9876d..03db2570a 100644
--- a/numpy/lib/shape_base.py
+++ b/numpy/lib/shape_base.py
@@ -32,7 +32,7 @@ def apply_along_axis(func1d,axis,arr,*args):
     if isscalar(res):
         outarr = zeros(outshape,asarray(res).dtype)
         outarr[ind] = res
-        Ntot = product(outshape)
+        Ntot = product(outshape,axis=0)
         k = 1
         while k < Ntot:
             # increment the index
@@ -48,7 +48,7 @@ def apply_along_axis(func1d,axis,arr,*args):
             k += 1
         return outarr
     else:
-        Ntot = product(outshape)
+        Ntot = product(outshape,axis=0)
         holdshape = outshape
         outshape = list(arr.shape)
         outshape[axis] = len(res)
@@ -326,12 +326,7 @@ def array_split(ary,indices_or_sections,axis = 0):
         Caveats:
            Currently, the default for axis is 0.  This
            means a 2D array is divided into multiple groups
-           of rows.  This seems like the appropriate default, but
-           we've agreed most other functions should default to
-           axis=-1.  Perhaps we should use axis=-1 for consistency.
-           However, we could also make the argument that NumPy
-           works on "rows" by default.  sum() sums up rows of
-           values.  split() will split data into rows.  Opinions?
+           of rows.  This seems like the appropriate default,
     """
     try:
         Ntotal = ary.shape[axis]
@@ -391,12 +386,7 @@ def split(ary,indices_or_sections,axis=0):
         Caveats:
            Currently, the default for axis is 0.  This
            means a 2D array is divided into multiple groups
-           of rows.  This seems like the appropriate default, but
-           we've agreed most other functions should default to
-           axis=-1.  Perhaps we should use axis=-1 for consistency.
-           However, we could also make the argument that NumPy
-           works on "rows" by default.  sum() sums up rows of
-           values.  split() will split data into rows.  Opinions?
+           of rows.  This seems like the appropriate default
     """
     try: len(indices_or_sections)
     except TypeError:
diff --git a/numpy/lib/tests/test_function_base.py b/numpy/lib/tests/test_function_base.py
index fdb2f270f..b548ce386 100644
--- a/numpy/lib/tests/test_function_base.py
+++ b/numpy/lib/tests/test_function_base.py
@@ -42,11 +42,11 @@ class test_all(NumpyTestCase):
 class test_average(NumpyTestCase):
     def check_basic(self):
         y1 = array([1,2,3])
-        assert(average(y1) == 2.)
+        assert(average(y1,axis=0) == 2.)
         y2 = array([1.,2.,3.])
-        assert(average(y2) == 2.)
+        assert(average(y2,axis=0) == 2.)
         y3 = [0.,0.,0.]
-        assert(average(y3) == 0.)
+        assert(average(y3,axis=0) == 0.)
 
         y4 = ones((4,4))
         y4[0,1] = 0
@@ -117,7 +117,7 @@ class test_amin(NumpyTestCase):
 class test_ptp(NumpyTestCase):
     def check_basic(self):
         a = [3,4,5,10,-3,-5,6.0]
-        assert_equal(ptp(a),15.0)
+        assert_equal(ptp(a,axis=0),15.0)
         b = [[3,6.0, 9.0],
              [4,10.0,5.0],
              [8,3.0,2.0]]
@@ -132,7 +132,7 @@ class test_cumsum(NumpyTestCase):
                       float32,float64,complex64,complex128]:
             a = array(ba,ctype)
             a2 = array(ba2,ctype)
-            assert_array_equal(cumsum(a), array([1,3,13,24,30,35,39],ctype))
+            assert_array_equal(cumsum(a,axis=0), array([1,3,13,24,30,35,39],ctype))
             assert_array_equal(cumsum(a2,axis=0), array([[1,2,3,4],[6,8,10,13],
                                                          [16,11,14,18]],ctype))
             assert_array_equal(cumsum(a2,axis=1),
@@ -153,7 +153,7 @@ class test_prod(NumpyTestCase):
                 self.failUnlessRaises(ArithmeticError, prod, a2, 1)
                 self.failUnlessRaises(ArithmeticError, prod, a)
             else:
-                assert_equal(prod(a),26400)
+                assert_equal(prod(a,axis=0),26400)
                 assert_array_equal(prod(a2,axis=0),
                                    array([50,36,84,180],ctype))
                 assert_array_equal(prod(a2,axis=-1),array([24, 1890, 600],ctype))
@@ -305,35 +305,35 @@ class test_filterwindows(NumpyTestCase):
         w=hanning(10)
         assert_array_almost_equal(w,flipud(w),7)
         #check known value
-        assert_almost_equal(sum(w),4.500,4)
+        assert_almost_equal(sum(w,axis=0),4.500,4)
 
     def check_hamming(self):
         #check symmetry
         w=hamming(10)
         assert_array_almost_equal(w,flipud(w),7)
         #check known value
-        assert_almost_equal(sum(w),4.9400,4)
+        assert_almost_equal(sum(w,axis=0),4.9400,4)
 
     def check_bartlett(self):
         #check symmetry
         w=bartlett(10)
         assert_array_almost_equal(w,flipud(w),7)
         #check known value
-        assert_almost_equal(sum(w),4.4444,4)
+        assert_almost_equal(sum(w,axis=0),4.4444,4)
 
     def check_blackman(self):
         #check symmetry
         w=blackman(10)
         assert_array_almost_equal(w,flipud(w),7)
         #check known value
-        assert_almost_equal(sum(w),3.7800,4)
+        assert_almost_equal(sum(w,axis=0),3.7800,4)
 
 
 class test_trapz(NumpyTestCase):
     def check_simple(self):
         r=trapz(exp(-1.0/2*(arange(-10,10,.1))**2)/sqrt(2*pi),dx=0.1)
         #check integral of normal equals 1
-        assert_almost_equal(sum(r),1,7)
+        assert_almost_equal(sum(r,axis=0),1,7)
 
 class test_sinc(NumpyTestCase):
     def check_simple(self):
@@ -348,7 +348,7 @@ class test_histogram(NumpyTestCase):
         v=rand(n)
         (a,b)=histogram(v)
         #check if the sum of the bins equals the number of samples
-        assert(sum(a)==n)
+        assert(sum(a,axis=0)==n)
         #check that the bin counts are evenly spaced when the data is from a linear function
         (a,b)=histogram(linspace(0,10,100))
         assert(all(a==10))
diff --git a/numpy/lib/tests/test_type_check.py b/numpy/lib/tests/test_type_check.py
index b95760c8b..fba89e4c1 100644
--- a/numpy/lib/tests/test_type_check.py
+++ b/numpy/lib/tests/test_type_check.py
@@ -91,7 +91,7 @@ class test_iscomplex(NumpyTestCase):
     def check_fail(self):
         z = array([-1,0,1])
         res = iscomplex(z)
-        assert(not sometrue(res))
+        assert(not sometrue(res,axis=0))
     def check_pass(self):
         z = array([-1j,1,0])
         res = iscomplex(z)
@@ -125,7 +125,7 @@ class test_isnan(NumpyTestCase):
     def check_goodvalues(self):
         z = array((-1.,0.,1.))
         res = isnan(z) == 0
-        assert_all(alltrue(res))
+        assert_all(alltrue(res,axis=0))
     def check_posinf(self):
         assert_all(isnan(array((1.,))/0.) == 0)
     def check_neginf(self):
@@ -145,7 +145,7 @@ class test_isfinite(NumpyTestCase):
     def check_goodvalues(self):
         z = array((-1.,0.,1.))
         res = isfinite(z) == 1
-        assert_all(alltrue(res))
+        assert_all(alltrue(res,axis=0))
     def check_posinf(self):
         assert_all(isfinite(array((1.,))/0.) == 0)
     def check_neginf(self):
@@ -165,7 +165,7 @@ class test_isinf(NumpyTestCase):
     def check_goodvalues(self):
         z = array((-1.,0.,1.))
         res = isinf(z) == 0
-        assert_all(alltrue(res))
+        assert_all(alltrue(res,axis=0))
     def check_posinf(self):
         assert_all(isinf(array((1.,))/0.) == 1)
     def check_posinf_scalar(self):
diff --git a/numpy/lib/utils.py b/numpy/lib/utils.py
index db7c00db6..deca8fa06 100644
--- a/numpy/lib/utils.py
+++ b/numpy/lib/utils.py
@@ -126,7 +126,7 @@ def who(vardict=None):
                 namestr = name
                 original=1
             shapestr = " x ".join(map(str, var.shape))
-            bytestr = str(var.itemsize*product(var.shape))
+            bytestr = str(var.itemsize*product(var.shape,axis=0))
             sta.append([namestr, shapestr, bytestr, var.dtype.name,
                         original])
 
diff --git a/numpy/linalg/linalg.py b/numpy/linalg/linalg.py
index 8dec67318..e0b3268a6 100644
--- a/numpy/linalg/linalg.py
+++ b/numpy/linalg/linalg.py
@@ -661,11 +661,11 @@ Singular values less than s[0]*rcond are treated as zero.
     if one_eq:
         x = array(ravel(bstar)[:n], dtype=result_t, copy=True)
         if results['rank']==n and m>n:
-            resids = array([sum((ravel(bstar)[n:])**2)], dtype=result_t)
+            resids = array([sum((ravel(bstar)[n:])**2,axis=0)], dtype=result_t)
     else:
         x = array(transpose(bstar)[:n,:], dtype=result_t, copy=True)
         if results['rank']==n and m>n:
-            resids = sum((transpose(bstar)[n:,:])**2).astype(result_t)
+            resids = sum((transpose(bstar)[n:,:])**2,axis=0).astype(result_t)
     st = s[:min(n,m)].copy().astype(_realType(result_t))
     return wrap(x), resids, results['rank'], st
 
@@ -687,7 +687,7 @@ def norm(x, ord=None):
        For vectors ord can be any real number including Inf or -Inf.
          ord = Inf, computes the maximum of the magnitudes
          ord = -Inf, computes minimum of the magnitudes
-         ord is finite, computes sum(abs(x)**ord)**(1.0/ord)
+         ord is finite, computes sum(abs(x)**ord,axis=0)**(1.0/ord)
 
        For matrices ord can only be one of the following values:
          ord = 2 computes the largest singular value
@@ -696,7 +696,7 @@ def norm(x, ord=None):
          ord = -1 computes the smallest column sum of absolute values
          ord = Inf computes the largest row sum of absolute values
          ord = -Inf computes the smallest row sum of absolute values
-         ord = 'fro' computes the frobenius norm sqrt(sum(diag(X.H * X)))
+         ord = 'fro' computes the frobenius norm sqrt(sum(diag(X.H * X),axis=0))
 
        For values ord < 0, the result is, strictly speaking, not a
        mathematical 'norm', but it may still be useful for numerical purposes.
diff --git a/numpy/numarray/functions.py b/numpy/numarray/functions.py
index 55922e9fe..95738d220 100644
--- a/numpy/numarray/functions.py
+++ b/numpy/numarray/functions.py
@@ -206,7 +206,7 @@ def fromfile(infile, type=None, shape=None, sizing=STRICT,
     ##file whose size may be determined before allocation, should be
     ##quick -- only one allocation will be needed.
     
-    recsize = dtype.itemsize * N.product([i for i in shape if i != -1])
+    recsize = dtype.itemsize * N.product([i for i in shape if i != -1],axis=0)
     blocksize = max(_BLOCKSIZE/recsize, 1)*recsize
 
     ##try to estimate file size
@@ -268,7 +268,7 @@ def fromstring(datastring, type=None, shape=None, typecode=None, dtype=None):
     if shape is None:
         count = -1
     else:
-        count = N.product(shape)*dtype.itemsize
+        count = N.product(shape,axis=0)*dtype.itemsize
     res = N.fromstring(datastring, count=count)
     if shape is not None:
         res.shape = shape
diff --git a/numpy/oldnumeric/random_array.py b/numpy/oldnumeric/random_array.py
index c5d163196..550af720c 100644
--- a/numpy/oldnumeric/random_array.py
+++ b/numpy/oldnumeric/random_array.py
@@ -166,7 +166,7 @@ def multinomial(trials, probs, shape=[]):
            trials is the number of trials in each multinomial distribution.
            probs is a one dimensional array. There are len(prob)+1 events.
            prob[i] is the probability of the i-th event, 0<=i<len(prob).
-           The probability of event len(prob) is 1.-Numeric.sum(prob).
+           The probability of event len(prob) is 1.-Numeric.sum(prob,axis=0).
 
        The first form returns a single 1-D array containing one multinomially
            distributed vector.
@@ -188,14 +188,14 @@ def poisson(mean, shape=[]):
 
 def mean_var_test(x, type, mean, var, skew=[]):
     n = len(x) * 1.0
-    x_mean = Numeric.sum(x)/n
+    x_mean = Numeric.sum(x,axis=0)/n
     x_minus_mean = x - x_mean
-    x_var = Numeric.sum(x_minus_mean*x_minus_mean)/(n-1.0)
+    x_var = Numeric.sum(x_minus_mean*x_minus_mean,axis=0)/(n-1.0)
     print "\nAverage of ", len(x), type
     print "(should be about ", mean, "):", x_mean
     print "Variance of those random numbers (should be about ", var, "):", x_var
     if skew != []:
-       x_skew = (Numeric.sum(x_minus_mean*x_minus_mean*x_minus_mean)/9998.)/x_var**(3./2.)
+       x_skew = (Numeric.sum(x_minus_mean*x_minus_mean*x_minus_mean,axis=0)/9998.)/x_var**(3./2.)
        print "Skewness of those random numbers (should be about ", skew, "):", x_skew
 
 def test():
@@ -205,12 +205,12 @@ def test():
     if (obj2[1] - obj[1]).any():
         raise SystemExit, "Failed seed test."
     print "First random number is", random()
-    print "Average of 10000 random numbers is", Numeric.sum(random(10000))/10000.
+    print "Average of 10000 random numbers is", Numeric.sum(random(10000),axis=0)/10000.
     x = random([10,1000])
     if len(x.shape) != 2 or x.shape[0] != 10 or x.shape[1] != 1000:
         raise SystemExit, "random returned wrong shape"
     x.shape = (10000,)
-    print "Average of 100 by 100 random numbers is", Numeric.sum(x)/10000.
+    print "Average of 100 by 100 random numbers is", Numeric.sum(x,axis=0)/10000.
     y = uniform(0.5,0.6, (1000,10))
     if len(y.shape) !=2 or y.shape[0] != 1000 or y.shape[1] != 10:
         raise SystemExit, "uniform returned wrong shape"
@@ -239,7 +239,7 @@ def test():
     print x
     if x.shape != (4,3,2): raise SystemExit, "multivariate_normal returned wrong shape"
     x = multivariate_normal(Numeric.array([-100,0,100]), Numeric.array([[3,2,1],[2,2,1],[1,1,1]]), 10000)
-    x_mean = Numeric.sum(x)/10000.
+    x_mean = Numeric.sum(x,axis=0)/10000.
     print "Average of 10000 multivariate normals with mean [-100,0,100]"
     print x_mean
     x_minus_mean = x - x_mean
@@ -262,7 +262,7 @@ def test():
     print "\nEach row is the result of 16 multinomial trials with probabilities [0.1, 0.5, 0.1 0.3]:"
     x = multinomial(16, [0.1, 0.5, 0.1], 8)
     print x
-    print "Mean = ", Numeric.sum(x)/8.
+    print "Mean = ", Numeric.sum(x,axis=0)/8.
 
 if __name__ == '__main__':
     test()
diff --git a/numpy/oldnumeric/rng_stats.py b/numpy/oldnumeric/rng_stats.py
index 60248247e..8c7fec433 100644
--- a/numpy/oldnumeric/rng_stats.py
+++ b/numpy/oldnumeric/rng_stats.py
@@ -9,7 +9,7 @@ def average(data):
 
 def variance(data):
     data = Numeric.array(data)
-    return Numeric.add.reduce((data-average(data))**2)/(len(data)-1)
+    return Numeric.add.reduce((data-average(data,axis=0))**2)/(len(data)-1)
 
 def standardDeviation(data):
     data = Numeric.array(data)
@@ -25,7 +25,7 @@ def histogram(data, nbins, range = None):
         data = Numeric.repeat(data,
                               Numeric.logical_and(Numeric.less_equal(data, max),
                                                   Numeric.greater_equal(data,
-                                                                        min)))
+                                                                        min)),axis=0)
     bin_width = (max-min)/nbins
     data = Numeric.floor((data - min)/bin_width).astype(Numeric.Int)
     histo = Numeric.add.reduce(Numeric.equal(