tests.py

The library is a C++/Python implementation of the variational building block framework introduced in

        m2 = nf.GetGaussian("m2", c0, c_5)
        v2 = nf.GetGaussian("v2", c0, c_5)

        nf.EvidenceNode(m1, mean=-1.0)
        nf.EvidenceNode(m2, mean=1.0)

        vx = nf.GetConstant("vx", 9.21)

        d = []
        s = []
        x = []
        for i in range(dim):
            d.append(nf.GetDiscreteDirichlet(Label("d", i), c))
            s.append(nf.GetMoG(Label("s", i), d[i]))
            s[i].AddComponent(m1, v1)
            s[i].AddComponent(m2, v2)
            x.append(nf.GetGaussian(Label("x", i), s[i], vx))
            x[i].Clamp(self.data[i])

        return net

    def costdump(self, vector, scalar):
        variables = ["c", "m1", "v1", "m2", "v2"]

        for l in variables:
            n1 = vector.GetVariable(l)
            n2 = scalar.GetVariable(l)
            print "%s : %f : %f" % (l, n1.Cost(), n2.Cost())

        variables = ["d", "s", "x"]

        for l in variables:
            c1 = vector.GetVariable(l).Cost()
            c2 = MLab.sum([n.Cost() for n in scalar.GetVariableArray(l)])
            print "%s : %f : %f" % (l, c1, c2)

    def comparenets(self, vector, scalar):
        dim = self.data.shape[0]

        #self.costdump(vector, scalar)

        # this alone is quite reliable indicator
        costdiff = abs(vector.Cost() - scalar.Cost())
        self.assert_(costdiff < dim * 1e-8, "costdiff = %f" % costdiff)

        s1 = GetVarV(vector.GetVariable("s"))
        s2 = GetVar(scalar.GetVariableArray("s"))
##
##        import plot
##        pl = plot.Plotter()
##        pl.multiplot(Numeric.log10([s1,s2]))

        self.assert_(MLab.max(Numeric.absolute(s1 - s2)) < 1e-8)

        for l in ["m1", "v1", "m2", "v2"]:
            n1 = GetMean(vector.GetVariable(l))
            n2 = GetMean(scalar.GetVariable(l))
            #print "%s : %f : %f" % (l, n1, n2)
            self.assert_(abs(n1 - n2) < 1e-8)

        d1 = GetDiscreteV(vector.GetVariable("d"))
        d2 = Numeric.transpose(GetDiscrete(scalar.GetVariableArray("d")))

        #print d1
        #print d2

        e = MLab.max(Numeric.reshape(
            Numeric.absolute(d1 - d2), (d1.shape[0]*d1.shape[1],)))

        self.assert_(e < 1e-8, "max e = %f" % e)

        c1 = vector.GetVariable("c")
        c2 = scalar.GetVariable("c")

        #print GetMeanV(c1)

        self.assert_(
            MLab.mean(Numeric.absolute(GetMeanV(c1) - GetMeanV(c2))) < 1e-8)
        self.assert_(
            MLab.mean(Numeric.absolute(GetVarV(c1) - GetVarV(c2))) < 1e-8)
        
    def test1(self):
        """Sameness of scalar and vector implementation of MoGs."""
        net1 = self.buildnet1()
        net2 = self.buildnet2()
        
        for i in range(50):
            net1.UpdateAll()
            net2.UpdateAll()

            #print "%d : %f : %f" % (i, net1.Cost(), net2.Cost())

        self.comparenets(net1, net2)

    def test2(self):
        """Load/save interface."""
        net1 = self.buildnet1()
        net2 = self.buildnet2()

        for i in range(50):
            net1.UpdateAll()
            net2.UpdateAll()

        self.comparenets(net1, net2)

        from PickleHelpers import load, save

        save(net2, "/tmp/bblocks_tmpnet")
        net2alt = load("/tmp/bblocks_tmpnet")

        self.comparenets(net1, net2alt)

class GaussRectTestCase(unittest.TestCase):
    def setUp(self):
        t = Numeric.arange(20) / float(20)
        s = Numeric.sin(4*Numeric.pi*t) + 0.5

        # N(0, 0.01^2)
        noise = Numeric.array(
            [-0.00020618, -0.00214195, -0.00860223,  0.0083709 ,
             0.00493691, -0.00303564,  0.00467232,  0.00875158,
             0.01462178,  0.00804469,  0.00791851,  0.00444479,
             -0.00398483, -0.00677404, -0.02249393,  0.00172111,
            -0.00420562,  0.00287512,  0.02004725,  0.00479894],
            Numeric.Float)

        self.data = s * (s > 0) + noise

    def buildstatnet1(self):
        net = PyNet.PyNet(self.data.shape[0])
        nf = PyNet.PyNodeFactory(net)

        c0 = nf.GetConstant("c0", 0.0)
        c_5 = nf.GetConstant("c_5", -5.0)

        ms = nf.GetGaussian("ms", c0, c_5)
        vs = nf.GetGaussian("vs", c0, c_5)

        nf.EvidenceNode(ms, mean=1.0)
        nf.EvidenceNode(vs, mean=0.0)

        s = nf.GetGaussRectV("s", ms, vs)
        r = nf.GetRectificationV("r", s)

        vx = nf.GetConstant("vx", -2*Numeric.log(0.01))
        #vx = nf.GetGaussian("vx", c0, c_5)
        x = nf.GetGaussianV("x", r, vx)

        x.Clamp(Helpers.Array2DV(self.data))

        return net

    def buildstatnet2(self):
        net = PyNet.PyNet(1)
        nf = PyNet.PyNodeFactory(net)

        c0 = nf.GetConstant("c0", 0.0)
        c_5 = nf.GetConstant("c_5", -5.0)

        ms = nf.GetGaussian("ms", c0, c_5)
        vs = nf.GetGaussian("vs", c0, c_5)

        nf.EvidenceNode(ms, mean=1.0)
        nf.EvidenceNode(vs, mean=0.0)

        vx = nf.GetConstant("vx", -2*Numeric.log(0.01))

        for i in range(self.data.shape[0]):
            s = nf.GetGaussRect(Label("s", i), ms, vs)
            r = nf.GetRectification(Label("r", i), s)
            x = nf.GetGaussian(Label("x", i), r, vx)
            x.Clamp(self.data[i])

        return net

    def comparenets(self, vector, scalar):
        dim = self.data.shape[0]

        costdiff = abs(vector.Cost() - scalar.Cost())
        self.assert_(costdiff < dim * 1e-8, "costdiff = %f" % costdiff)

        s1 = GetMeanV(vector.GetVariable("s"))
        s2 = GetMean(scalar.GetVariableArray("s"))

        e = MLab.max(Numeric.absolute(s1 - s2))
        self.assert_(e < 1e-8)
        
        r1 = GetMeanV(vector.GetNode("r"))
        r2 = GetMean(scalar.GetNodeArray("r"))
        
        e = MLab.max(Numeric.absolute(r1 - r2))
        self.assert_(e < 1e-8)

        if vector.GetVariable("ms") is not None:
            self.assert_(abs(GetMean(vector.GetVariable("ms"))
                             - GetMean(scalar.GetVariable("ms"))) < 1e-8)
        else:
            self.assert_(abs(GetMean(vector.GetVariable("m0"))
                             - GetMean(scalar.GetVariable("m0"))) < 1e-8)

        self.assert_(abs(GetMean(vector.GetVariable("vs"))
                         - GetMean(scalar.GetVariable("vs"))) < 1e-8)
        
    def test1(self):
        """Vector/scalar comparison with static model."""

        net1 = self.buildstatnet1()
        net2 = self.buildstatnet2()
        
        for i in range(50):
            net1.UpdateAll()
            net2.UpdateAll()
            
            #print "%d : %f : %f" % (i, net1.Cost(), net2.Cost())

        self.comparenets(net1, net2)

    def builddynnet1(self):
        net = PyNet.PyNet(self.data.shape[0])
        nf = PyNet.PyNodeFactory(net)

        c0 = nf.GetConstant("c0", 0.0)
        c_5 = nf.GetConstant("c_5", -5.0)

        m0 = nf.GetGaussian("ms", c0, c_5)
        vs = nf.GetGaussian("vs", c0, c_5)

        proxy = nf.GetProxy("proxy", "s")
        delay = nf.GetDelayV("delay", m0, proxy)

        s = nf.GetGaussRectV("s", delay, vs)
        r = nf.GetRectificationV("r", s)

        vx = nf.GetConstant("vx", -2*Numeric.log(0.01))
        #vx = nf.GetGaussian("vx", c0, c_5)
        x = nf.GetGaussianV("x", r, vx)

        x.Clamp(Helpers.Array2DV(self.data))

        net.ConnectProxies()

        return net

    def builddynnet2(self):
        net = PyNet.PyNet(1)
        nf = PyNet.PyNodeFactory(net)

        c0 = nf.GetConstant("c0", 0.0)
        c_5 = nf.GetConstant("c_5", -5.0)

        m0 = nf.GetGaussian("ms", c0, c_5)
        vs = nf.GetGaussian("vs", c0, c_5)

        vx = nf.GetConstant("vx", -2*Numeric.log(0.01))

        prev = m0

        for i in range(self.data.shape[0]):
            s = nf.GetGaussRect(Label("s", i), prev, vs)
            prev = s
            r = nf.GetRectification(Label("r", i), s)
            x = nf.GetGaussian(Label("x", i), r, vx)
            x.Clamp(self.data[i])

        return net

    def test2(self):
        """Vector/scalar comparison with dynamic model."""

        net1 = self.builddynnet1()
        net2 = self.builddynnet2()
        
        for i in range(100):
            net1.UpdateAll()
            net2.UpdateAll()
            
            #print "%d : %f : %f" % (i, net1.Cost(), net2.Cost())

        self.comparenets(net1, net2)

    def test3(self):
        """Load/save vector."""
        from PickleHelpers import load, save
        
        ref = self.buildstatnet2()
        net = self.buildstatnet1()
        for i in range(100):
            ref.UpdateAll()
            net.UpdateAll()
        self.comparenets(net, ref)

        save(net, "/tmp/bblocks_tmpnet")
        alt = load("/tmp/bblocks_tmpnet")
        self.comparenets(alt, ref)

        ref = self.builddynnet2()
        net = self.builddynnet1()
        for i in range(100):
            ref.UpdateAll()
            net.UpdateAll()
        self.comparenets(net, ref)

        save(net, "/tmp/bblocks_tmpnet")
        alt = load("/tmp/bblocks_tmpnet")
        self.comparenets(alt, ref)
        
    def test4(self):
        """Load/save scalar."""
        from PickleHelpers import load, save

        ref = self.buildstatnet1()
        net = self.buildstatnet2()
        for i in range(100):
            ref.UpdateAll()
            net.UpdateAll()
        self.comparenets(ref, net)

        save(net, "/tmp/bblocks_tmpnet")
        alt = load("/tmp/bblocks_tmpnet")
        self.comparenets(ref, alt)

        ref = self.builddynnet1()
        net = self.builddynnet2()
        for i in range(100):
            ref.UpdateAll()
            net.UpdateAll()
        self.comparenets(ref, net)

        save(net, "/tmp/bblocks_tmpnet")
        alt = load("/tmp/bblocks_tmpnet")
        self.comparenets(ref, alt)

if __name__ == '__main__':
    cases = [s for s in dir() if s.endswith("TestCase")]
    suites = []

    if len(sys.argv) > 1:
        for i in range(1, len(sys.argv)):
            if sys.argv[i] in cases:
                suites.append(unittest.makeSuite(eval(sys.argv[i]), "test"))
            else:
                print sys.argv[i] + "?"
    else:
        for case in cases:
            suites.append(unittest.makeSuite(eval(case), "test"))

    bigsuite = unittest.TestSuite(suites)
    runner = unittest.TextTestRunner()

    runner.run(bigsuite)

##    import pdb
##    pdb.run("runner.run(bigsuite)")
##