CODEWITHSUNDEEP
pythonnumpystatisticslinear-algebraarray-broadcasting
User1291
User1291

Reputation: 8182

numpy - do operation along specified axis

So I want to implement a matrix standardisation method. To do that, I've been told to

subtract the mean and divide by the standard deviation for each dimension

And to verify:

after this processing, each dimension has zero mean and unit variance.

That sounds simple enough ...

import numpy as np
def standardize(X : np.ndarray,inplace=True,verbose=False,check=False):

    ret = X
    if not inplace:
        ret = X.copy()

    ndim = np.ndim(X)

    for d in range(ndim):
        m = np.mean(ret,axis=d)
        s = np.std(ret,axis=d)

        if verbose:
            print(f"m{d} =",m)
            print(f"s{d} =",s)

        # TODO: handle zero s
        # TODO: subtract m along the correct axis
        # TODO: divide by s along the correct axis

    if check:    
        means = [np.mean(X,axis=d) for d in range(ndim)]
        stds  = [np.std(X,axis=d)  for d in range(ndim)]
        if verbose:
            print("means=\n",means)
            print("stds=\n",stds)

        assert all(all(m < 1e-15 for m in mm) for mm in means)
        assert all(all(s == 1.0 for s in ss) for ss in stds)

    return ret

e.g. for ndim == 2, we could get something like

A=
 [[ 0.40923704  0.91397416  0.62257397]
  [ 0.15614258  0.56720836  0.80624135]]
m0 = [ 0.28268981  0.74059126  0.71440766]  # can broadcast with ret -= m0
s0 = [ 0.12654723  0.1733829   0.09183369]  # can broadcast with ret /= s0
m1 = [ 0.33333333 -0.33333333]  # ???
s1 = [ 0.94280904  0.94280904]  # ???

How do I do that?

Judging by Broadcast an operation along specific axis in python , I thought I may be looking for a way to create

m[None, None, None, .., None, : , None, None, .., None]

Where there is exactly one : at index d.

But even if I knew how to do that, I'm not sure it'd work.

Upvotes: 0

Views: 383

Answers (1)

Stefan T. Huber
Stefan T. Huber

Reputation: 51

You can swap your axes such that the first axes is the one you want to normalize. This should also work inplace, since swapaxes just returns a view on your data.

Using the numpy command swapaxes:

for d in range(ndim):

    m = np.mean(ret,axis=d)
    s = np.std(ret,axis=d)

    ret = np.swapaxes(ret, 0, d)

    # Perform Normalisation of Axis
    ret -= m
    ret /= s

    ret = np.swapaxes(ret, 0, d)

Upvotes: 1

Related Questions