Consider a small 2d array:

In [180]: A=np.arange(12).reshape(3,4)
In [181]: A
Out[181]: 
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])

Sum across rows; the result is a (3,) array

In [182]: A.sum(axis=1)
Out[182]: array([ 6, 22, 38])

But to sum (or divide) A by the sum requires reshaping

In [183]: A-A.sum(axis=1)
...
ValueError: operands could not be broadcast together with shapes (3,4) (3,) 
In [184]: A-A.sum(axis=1)[:,None]   # turn sum into (3,1)
Out[184]: 
array([[ -6,  -5,  -4,  -3],
       [-18, -17, -16, -15],
       [-30, -29, -28, -27]])

If I use keepdims, "the result will broadcast correctly against" A.

In [185]: A.sum(axis=1, keepdims=True)   # (3,1) array
Out[185]: 
array([[ 6],
       [22],
       [38]])
In [186]: A-A.sum(axis=1, keepdims=True)
Out[186]: 
array([[ -6,  -5,  -4,  -3],
       [-18, -17, -16, -15],
       [-30, -29, -28, -27]])

If I sum the other way, I don't need the keepdims. Broadcasting this sum is automatic: A.sum(axis=0)[None,:]. But there's no harm in using keepdims.

In [190]: A.sum(axis=0)
Out[190]: array([12, 15, 18, 21])    # (4,)
In [191]: A-A.sum(axis=0)
Out[191]: 
array([[-12, -14, -16, -18],
       [ -8, -10, -12, -14],
       [ -4,  -6,  -8, -10]])

If you prefer, these actions might make more sense with np.mean, normalizing the array over columns or rows. In any case it can simplify further math between the original array and the sum/mean.

@Ney @hpaulj is correct, you need to experiment, but I suspect you don't realize that summation for some arrays can occur along axes. Observe the following which reading the documentation

>>> a
array([[0, 0, 0],
       [0, 1, 0],
       [0, 2, 0],
       [1, 0, 0],
       [1, 1, 0]])
>>> np.sum(a, keepdims=True)
array([[6]])
>>> np.sum(a, keepdims=False)
6
>>> np.sum(a, axis=1, keepdims=True)
array([[0],
       [1],
       [2],
       [1],
       [2]])
>>> np.sum(a, axis=1, keepdims=False)
array([0, 1, 2, 1, 2])
>>> np.sum(a, axis=0, keepdims=True)
array([[2, 4, 0]])
>>> np.sum(a, axis=0, keepdims=False)
array([2, 4, 0])

You will notice that if you don't specify an axis (1st two examples), the numerical result is the same, but the keepdims = True returned a 2D array with the number 6, whereas, the second incarnation returned a scalar. Similarly, when summing along axis 1 (across rows), a 2D array is returned again when keepdims = True. The last example, along axis 0 (down columns), shows a similar characteristic... dimensions are kept when keepdims = True.
Studying axes and their properties is critical to a full understanding of the power of NumPy when dealing with multidimensional data.

Ok, CAUTION I don't have my masters in numpyology yet, but just playing around, I came up with:

>>> np.average(a,axis=-1).repeat(a.shape[-1]).reshape(a.shape)
array([[[ 0.2       ,  0.2       ,  0.2       ],
        [ 0.29999998,  0.29999998,  0.29999998]],

       [[ 0.40000001,  0.40000001,  0.40000001],
        [ 0.69999999,  0.69999999,  0.69999999]]], dtype=float32)