ENH: recurrence relation calculations

[jreback]

I think an implementation of a general recurrent relation might be interesting. These
are not able to evaluated in numpy/numexpr directly

see here:
http://stackoverflow.com/questions/9984373/calculate-two-coupled-equation-with-numpy
http://stackoverflow.com/questions/21336794/python-recursive-vectorization-with-timeseries/21338665
http://stackoverflow.com/questions/23996260/python-pandas-using-the-previous-value-in-dataframe/23997030#23997030

R[i<=k] = start(i) 
R[i] = R[i-k]*a(i) + b(i) 

in the case where a(i) and b(i) are constants then this can be completely cythonized
if a and/or b are passed functions then the shell can still be cythonized (with the function called at each iteration)

Example

values = df['input']

In [62]: result = np.empty(len(values))

In [65]: result[0] = values[0]

In [66]: for i in xrange(len(result)):
   ....:     if i > 0:
   ....:         result[i] = result[i-1]*.8
   ....:         

In [67]: result
Out[67]: array([ 5.    ,  4.    ,  3.2   ,  2.56  ,  2.048 ,  1.6384])

I propose that pandas have a recurrence wrapper (that can then call a pure constant cython version and a cython version that takes functions)
Probably just for a 1-dim array atm

In theory this can also be done by a rolling_apply that has memory (or I guess that is sort of the same thing)

related #4546
https://groups.google.com/forum/#!topic/pydata/0MCWhwurOWs

[timmie]

Here a python function solving this partly:

def apply_sequencial(series, func, dtype_internal=float, dtype_res=dtype, precis=None):
    result_ser = series.copy()
    #result_ser = result_ser.astype(dtype_internal)
    result_ser .astype(dtype_internal)
    try:
        for i in xrange(len(result_ser)):
            if i > 0:
                input = (result_ser[i-1])
                print input, result_ser[i]
                result_ser[i] = func(input)
    except OverflowError:
        result_ser = result_ser.astype(dtype_internal)
        for i in xrange(len(result_ser)):
            if i > 0:
                input = (result_ser[i-1])
                result_ser[i] = func(input)
    if dtype_res:       
        result_ser = result_ser.astype(dtype_res)
    elif result_ser.dtype == dtype_internal:
        result_ser = result_ser.astype(series.dtype)
    else:
        result_ser = result_ser.astype(series.dtype)

    if precis:
        from numpy import round
        result_ser = round(result_ser, precis)

    return result_ser

and can be called like:

import pandas as pd
import numpy as np
ts = pd.Series(range(1, 7), index=pd.date_range('1/1/2000', periods=6))
data = (range(1, 13))
data = np.array(data)
df = pd.DataFrame(data, index=ts.index, columns=['A', 'B'])
df['C'] = apply_sequencial(df.A, func=(lambda x: (x * (0.8))), dtype_res=float, precis=4)

Maybe a hack from performance point but works ;-)

[dirkbike]

I think it's important to make a recurrent relation apply to an entire data frame and not just a single series, otherwise this enhancement becomes equivalent to the rolling_apply function as mentioned earlier.

[jreback]

of course, that is just an example of the impl; should for sure handle 2-d

[dirkbike]

Or even n-dimensions. The recurrence relation could take (n-1)-dimensional slices of the data frame along the remaining dimension and provide the designated apply function with the previous slice. That apply function would then return a new (n-1)-dimensional data frame representing the new current slice.

[jreback]

you could do that with a custom function, when I mean 2-d I mean it basically respects the axis parameter, like apply

[cpcloud]

Should get a working 1d version before we start generalizing

[MichaelWS]

This would be a great piece. I use talib to do tons of financial calculations that are recurrence relations

[jbrockmendel]

this could be neat, but probably should live outside of pandas

[jreback]

these are actually trivially implemented with numba :->

[glopesdev]

@jreback would you elaborate how these are trivially implemented with numba?

[jreback]

look in the window functions we already have ewm implemented this way