window_regres.py

"""
x[t] = fp(u)[t] = tanh(v[t])
v[t] = mean( u[t-T+1:t+1] ) = ( u[t-T+1] + u[t-T+2] + ... + u[t-1] + u[t] )/T

E = sum( 0.5 * (x - y)**2 )
dE/dv[t] = (x[t] - y[t]) * (1 - x[t]**2)
dE/du[t] = mean( dE/dv[t:t+T] )
         = ( dE/dv[t+T-1] + dE/dv[t+T-2] + ... + dE/dv[t+1] + dE/dv[t] )/T

See result here: http://www.flickr.com/photos/arataka/tags/up20091213135244/show/
"""

import pylab
import numpy
tanh = numpy.tanh

class WindowRegres(object):

    def __init__(self, steps, wlen, _learn_rate=5, title=''):
        self.t          = numpy.arange(steps+wlen)
        self.x          = numpy.zeros((steps,), dtype=numpy.float)
        self.x          = numpy.zeros((steps,), dtype=numpy.float)
        self.v          = numpy.zeros((steps,), dtype=numpy.float)
        self.u          = numpy.zeros((steps+wlen,), dtype=numpy.float)
        self.dEdu       = numpy.zeros((steps+wlen,), dtype=numpy.float)
        self.dlist      = range(wlen+1)
        self.wlen       = wlen
        self.learn_rate = float(_learn_rate)/steps
        self.title = title

    def fp(self):
        u = self.u
        v = self.v
        v.fill(0)
        steps = self.v.shape[0]
        wlen = self.wlen
        for d in self.dlist:
            v += u[d:steps+d]
        v /= (wlen+1)
        self.x = tanh(v)

    def bp(self):
        y = self.y
        x = self.x
        steps = self.v.shape[0]
        wlen = self.wlen
        self.E = numpy.sum( 0.5*(x-y)**2 )/steps
        dEdv = (x-y) * (1-x**2)
        dEdu = self.dEdu
        dEdu.fill(0)
        for d in self.dlist:
            dEdu[d:steps+d] += dEdv
        #dEdu /= (wlen+1)

        self.u -= self.learn_rate * dEdu

    def learn(self, epochs, draw_step=None):
        self.rec = numpy.zeros((epochs,1), dtype=numpy.float)
        for ep in xrange(epochs):
            self.fp()
            self.bp()
            self.rec[ep,0] = self.E
            if draw_step != None and ep % draw_step == 0:
                self.plot()
        self.plot()

    def plot(self):
        wlen = self.wlen
        pylab.clf()
        
        pylab.subplot(311)
        pylab.plot(self.rec[:,0])
        pylab.ylim(0,None)
        pylab.title(self.title)
        pylab.ylabel('Error')
        
        pylab.subplot(312)
        pylab.plot(self.t[wlen:], self.y, linewidth=0.5, label='y')
        pylab.plot(self.t[wlen:], self.x, linewidth=2.0, label='x')
        pylab.ylim(-1.05,1.05)
        pylab.legend()
        pylab.ylabel('x, y')
        
        pylab.subplot(313)
        pylab.plot(self.t, self.u/(wlen+1), label='u/T')
        pylab.plot(self.t[wlen:], self.v, label='v')
        pylab.legend()
        pylab.ylabel('u/T, v')
        
        pylab.draw()
        pylab.show()


tch_len = 500
for wlen in range(0,60,10):
    wr = WindowRegres(tch_len,wlen,10, title='T=%d'%wlen)
    t = 2 * numpy.pi * numpy.arange(tch_len) / 50
    wr.y = numpy.sin(t) * 0.8 + numpy.random.randn(tch_len) * 0.5
    wr.learn(1000)
    pylab.savefig('T%03d.png'%wlen)