Htykuut

I have leveraged the rolling window examples using as_strided to create various sliding versions of numpy functions.

def std(a,window):

    shape = a.shape[:-1] + (a.shape[-1] - window + 1, window)

    strides = a.strides + (a.strides[-1],)

    return np.std(np.lib.stride_tricks.as_strided(a, shape=shape,strides=strides),axis=1)

Now I'm attempt to leverage the same as_strided method on a linear regression function. y = a + bx

def linear_regress(x,y):

  sum_x = np.sum(x)

  sum_y = np.sum(y)

  sum_xy = np.sum(np.multiply(x,y))

  sum_xx = np.sum(np.multiply(x,x))

  sum_yy = np.sum(np.multiply(y,y))

  number_of_records = len(x)

  A = (sum_y*sum_xx - sum_x*sum_xy)/(number_of_records*sum_xx - sum_x*sum_x)

  B = (number_of_records*sum_xy - sum_x*sum_y)/(number_of_records*sum_xx - sum_x*sum_x)

  return A + B*number_of_records

You can get the same result using stats

from scipy import stats

slope, intercept, r_value, p_value, std_err = stats.linregress(x,p)

xValue = len(x)

y = (slope + intercept*xValue) 

I am not sure how to fit the above functions into the as_strided method when two arrays are passed. I guess I would have to create two shapes and pass both through?

def rolling_lr(x,y,window):

   shape = y.shape[:-1] + (y.shape[-1] - window + 1, window)

   strides = y.strides + (y.strides[-1],)

   return linear_regress(np.lib.stride_tricks.as_strided(x,y shape=shape, strides=strides))

Any help is appreciated.

Here is what i came up with. Not the prettiest but works.

def linear_regress(x,y,window):

    x_shape = x.shape[:-1] + (x.shape[-1] - window + 1, window)

    x_strides = x.strides + (x.strides[-1],)

    y_shape = y.shape[:-1] + (y.shape[-1] - window + 1, window)

    y_strides = y.strides + (y.strides[-1],)

    sx = np.lib.stride_tricks.as_strided(x,shape=x_shape, strides=x_strides)

    sy = np.lib.stride_tricks.as_strided(y,shape=y_shape, strides=y_strides)

    sum_x = np.sum(sx,axis=1)

    sum_y = np.sum(sy,axis=1)

    sum_xy = np.sum(np.multiply(sx,sy),axis=1)

    sum_xx = np.sum(np.multiply(sx,sx),axis=1)

    sum_yy = np.sum(np.multiply(sy,sy),axis=1)

    m = (sum_y*sum_xx - sum_x*sum_xy)/(window*sum_xx - sum_x*sum_x)

    b = (window*sum_xy - sum_x*sum_y)/(window*sum_xx - sum_x*sum_x)      

Interesting, I have never seen the stride function. The docs do warn about this method however. An alternative method would be using linear algebra to do the regression on the windows. Here is a trivial example:

from numpy import *



# generate points

N            = 30

x            = linspace(0, 10, N)[:, None]

X            = ones((N, 1)) * x

Y            = X * array([1, 30])  + random.randn(*X.shape)*1e-1

XX           = concatenate((ones((N,1)), X), axis = 1)



# window data

windowLength = 10

windows = array([roll(

          XX, -i * windowLength, axis = 0)[:windowLength, :]

           for i in range(len(XX) - windowLength)])

windowsY = array([roll(Y, -i * windowLength, axis = 0)[:windowLength, :]

           for i in range(len(Y) - windowLength)])





# linear regression on windows

reg = array([

              ((linalg.pinv(wx.T.dot(wx))).dot(wx.T)).dot(wy) for 

                wx, wy in zip(windows, windowsY)])



# plot regression on windows

from matplotlib import style

style.use('seaborn-poster')

from matplotlib.pyplot import subplots, cm

fig, ax = subplots()



colors = cm.tab20(linspace(0, 1, len(windows)))

for win, color, coeffs, yi in zip(windows, colors, reg, windowsY):



    ax.plot(win, yi,'.', alpha = .5, color = color)

    ax.plot(win[:, 1], win.dot(coeffs), alpha = .5, color = color)

    x += 1

ax.set(**dict(xlabel = 'x', ylabel = 'y'))

Which produces: