Revert time averaged values to instant values with numpy

user_4390189 0
Viewed 46

I am having a 3D field of an averaged quantity. Now I am looking for a pythonic way to revert the averaging to get instantaneous values for each time stamp. Note: The averaging is made from the beginning of the whole period. It is like rolling mean with the window size adapted to the index of the value to be averaged.

For better clarification I give an 1D example:

import numpy as np

input_array = np.array([
       [0.      ],
       [0.5     ],
       [1.      ],
       [1.5     ],
       [2.      ],
       [2.5     ],
       [3.      ],
       [3.25    ],
       [3.333333],
       [3.3     ],
       [3.181818],
       [3.      ],
       [2.769231]
])

exp_result = de_average(input_array)

The expected result exp_result should look like this:

exp_result= np.array([
       [0],
       [1],
       [2],
       [3],
       [4],
       [5],
       [6],
       [5],
       [4],
       [3],
       [2],
       [1],
       [0]])
grib python-3.x numpy python
2 Answers

Something like this should work. You're problem is very related to the reverse cumsum problem, so I'm using part of the solution given here.

from itertools import tee 

def pairwise(iterable):
    "s -> (s0,s1), (s1,s2), (s2, s3), ..."
    a, b = tee(iterable)
    next(b, None)
    return zip(a, b)

def inverse_cumsum(cumulative):
    """ Inverse an array obtained by cumulative sum
       [1, 3, 6, 10, 15, 21, 28, 36, 45] -> [1, 2, 3, 4, 5, 6, 7, 8, 9] 
    """
    yield cumulative[0]
    for a, b in pairwise(cumulative):
        yield b-a

def inverse_average(a, decimals=1):
   """ Inverse an array averaged (where each entry i is the average up to i) """ 
   deav = inverse_cumsum([a * (i + 1) for i, a in enumerate(a)])
   return np.array(list(deav)).round(decimals)

inverse_average(input_array)
user_74417570

For 1-D numpy arrays:

import numpy as np

def de_average(input_array):
     timesteps = np.arange(1, input_array.shape[0] + 1).resize(-1, 1)
     _sum = timesteps.dot(input_array.T).diagonal()
     original = np.empty(input_array.shape, dtype=int)
     original[0,0] = _sum[0]
     original[1:,0] = _sum[1:] - _sum[:-1]
     return original

As mentioned in the comments, this will not work for higher dimensional arrays.

For 2-D and 3-D numpy arrays, try:

import numpy as np

# de-averages along axis 0
def de_average(input_array):
     _sum = np.apply_along_axis(
            lambda x:x*range(1, input_array.shape[0]+1), 
            0, 
            input_array
     )
     original = np.array(input_array)
     original[1:] = _sum[1:] - _sum[:-1]
     return original.astype(int)
user_143936890
Related