Changeset 165

Show
Ignore:
Timestamp:
04/29/08 17:21:15 (7 months ago)
Author:
edsuom
Message:

ProjectManager? tested OK!

Files:

Legend:

Unmodified
Added
Removed
Modified
Copied
Moved

  • projects/AsynCluster/trunk/svpmc/model.py

    r164 r165  
    1 # svpmc: Stochastic Volatility Inference via Population Monte Carlo 
     1# svPMC: Stochastic Volatility Inference via Population Monte Carlo 
    22# 
    33# Copyright (C) 2007-2008 by Edwin A. Suominen, http://www.eepatents.com 
    …… 
    7474 
    7575    def likelihood(paramContainer, wiggle): 
    76         L, v = M.likelihood(paramContainer, wiggle) 
    77         return Packer(L, v)() 
     76        L, v, h = M.likelihood(paramContainer, wiggle) 
     77        return Packer(L, v, h)() 
    7878 
    7979    ########################################################################### 
    8080    """ 
    8181    def __init__(self, paramNames, tsList): 
    82         self.mgr = projectManager 
    8382        self.modelObj = Model(paramNames=paramNames, tsList=tsList) 
    8483        self.queue = NullQueue(1) 
    …… 
    124123        undergoing some nested MCMC. 
    125124 
    126         Returns a deferred that fires with the log-likelihood when the 
    127         likelihood computation is done and the paramContainer has been updated. 
     125        Returns a deferred that fires (with C{None}) when the likelihood 
     126        computation is done and the paramContainer has been updated. 
    128127        """ 
    129128        def gotPackedLikelihood(result): 
    …… 
    132131                # infinitely low likelihood 
    133132                paramContainer.Lx = -s.inf 
    134             up = pack.Unpacker(result) 
    135             L = paramContainer.L = up() 
    136             paramContainer.z = up() 
    137             return L 
     133                paramContainer.z = [] 
     134                paramContainer.h = [] 
     135            else: 
     136                return gotLikelihood(list(pack.Unpacker(result))) 
    138137 
    139138        def gotLikelihood(result): 
    140             L = paramContainer.L = result[0] 
    141             paramContainer.z = result[1] 
    142             return L 
     139            for k, name in enumerate(('Lx', 'z', 'h')): 
     140                setattr(paramContainer, name, result[k]) 
    143141         
    144142        if (remote or self.remoteMode) and not local: 
    …… 
    338336        where each conditional draw is done with a Metropolis-Hastings step. 
    339337         
    340         Returns the likelihood of the innovations in I{x} given the simulated 
    341         log-volatilities. Updates the IID normal variates in place. 
     338        Updates the IID normal variates in place. Returns the likelihood of the 
     339        innovations in I{x}, given the simulated log-volatilities, along with 
     340        the last sample's multivariate log-volatility value. 
    342341 
    343342        The method will account for the effect of log-volatility proposals on 
    …… 
    379378            h0 = LV[2][:,0] 
    380379            L_total += LV[0] 
    381         return L_total 
     380        return L_total, h0 
    382381     
    383382 
    …… 
    408407        the simulated log-volatilities. 
    409408 
    410         The value of the log-likelihood at the end of the log-volatility 
    411         simulation is returned in a tuple, followed by the IID normal 
    412         variates underlying the log-volatilities drawn at that point. 
     409        The method returns a list containing the following values reached at 
     410        the end of the log-volatility simulation: 
     411         
     412            1. The value of the log-likelihood. 
     413 
     414            2. The IID normal variates underlying the log-volatilities. 
     415 
     416            3. The last sample's multivariate log-volatility value. 
    413417         
    414418        The returned likelihood does not consider the prior probability of the 
    …… 
    435439        for k in xrange(self.Mv-1): 
    436440            self.hybridGibbs(z, x, rv, sigma) 
    437         return self.hybridGibbs(z, x, rv, sigma), z 
    438      
     441        return list(self.hybridGibbs(z, x, rv, sigma)) + [z] 
     442     

  • projects/AsynCluster/trunk/svpmc/params.py

    r164 r165  
    391391        newVersion.Lj = Lj 
    392392        return newVersion 
    393  
    394  
    395  

  • projects/AsynCluster/trunk/svpmc/project.py

    r164 r165  
    6565        paramTitles, dimensions = self._setupParams() 
    6666        self._setupCDF( 
    67             os.path.join(specDir, ncFile), tsList, seriesTitles, paramTitles) 
     67            os.path.join(specDir, ncFileName), 
     68            tsList, seriesTitles, paramTitles, dimensions) 
    6869        self.mgr = model.ModelManager(self.paramNames, tsList) 
    6970     
    …… 
    219220            raise ValueError( 
    220221                "Each iteration record must have %d parameters" % self.m) 
     222        # Write the parameters to their variables 
    221223        for name in self.paramNames: 
    222224            var = self.cdf.variables[name] 
     225            # Iterate over a 1-D FlexArray of the population's parameter arrays 
    223226            for member, paramArray in enumerate(getattr(parameters, name)): 
    224227                var[self.iteration, member, :] = paramArray.astype('f') 
     228        # Write the last sample's multivariate log-volatility value for each 
     229        # population member 
     230        var = self.cdf.variables['log_volatility'] 
     231        for member, h in enumerate(parameters.h): 
     232            var[self.iteration, member, :] = h.astype('f') 
     233        # Write the log-volatility value for each population member 
     234        var = self.cdf.variables['log_likelihood'] 
     235        var[self.iteration, :] = parameters.Lx.astype('f') 
     236        # Done writing this iteration of the CDF record 
    225237        self.iteration += 1 
    226      
    227          
    228          
     238 
     239    def done(self): 
     240        """ 
     241        Call this when the project is done. 
     242        """ 
     243        self.cdf.close() 

  • projects/AsynCluster/trunk/svpmc/test/test_params.py

    r161 r165  
    4848    def test_get_set_2d(self): 
    4949        x = self._make_stringArray(3, 4) 
    50         self.failUnlessEqual(len(x), 12) 
     50        self.failUnlessEqual(len(x), 3) 
    5151        self.failUnlessEqual(x[0,0], '0:0') 
    5252        self.failUnlessEqual(x[2,3], '2:3') 
    …… 
    105105        def __init__(self, c): 
    106106            self.c = c 
     107            self.array = s.ones((2,3)) 
    107108        def __add__(self, other): 
    108109            self.c += other 
    …… 
    111112        def addThese(self, x, y): 
    112113            return self.c + x + y 
     114        def addArray(self, x): 
     115            return self.array + x 
    113116 
    114117    class Scaler(object): 
    …… 
    133136                self.failUnlessEqual(z[j, k].c, 100*j+k+y[j, k]) 
    134137 
    135     def test_getattr(self): 
     138    def test_getattr_scalar(self): 
    136139        y = self.x.c 
    137140        for j in xrange(4): 
    138141            for k in xrange(4): 
    139142                self.failUnlessEqual(y[j, k], 100*j+k) 
     143 
     144    def test_getattr_array(self): 
     145        y = self.x.array 
     146        for j in xrange(4): 
     147            for k in xrange(4): 
     148                self.failUnlessEqual(y[j, k].shape, (2,3)) 
    140149 
    141150    def test_setattr_scalar(self): 
    …… 
    165174                    "Error doing addThese method for j=%d, k=%d" % (j,k)) 
    166175 
     176    def test_call_arrayArgAndResult(self): 
     177        yExpected = s.ones((2,3)) + s.arange(6).reshape(2,3) 
     178        y = self.x.call('addArray', s.arange(6).reshape(2,3)) 
     179        for j in xrange(4): 
     180            for k in xrange(4): 
     181                self.failUnless( 
     182                    s.equal(y[j, k], yExpected).all(), 
     183                    "Error doing addArray method for j=%d, k=%d" % (j,k)) 
     184 
    167185    def test_call_flexArrayArg(self): 
    168186        arrayScale = params.FlexArray(4, 4) 
    …… 
    186204    def test_uniformPrior(self): 
    187205        p = params.Prior(dname='uniform', loc=1, scale=2) 
    188         X = p.rvs(1000) 
     206        X = s.array([p.rvs() for k in xrange(1000)]) 
    189207        self.failUnlessEqual(X.shape, (1000,)) 
    190208        self.failUnless(s.greater(X, 1).all()) 
    191209        self.failUnless(s.less(X, 3).all()) 
    192         self.failUnless(s.equal(p.pdf(X), 0.5).all()) 
     210        self.failUnless(s.equal(s.array([p.pdf(x) for x in X]), 0.5).all()) 
    193211 
    194212 

  • projects/AsynCluster/trunk/svpmc/test/test_project.py

    r164 r165  
    2222import os.path 
    2323import scipy as s 
    24 from scipy import stats, random, linalg, signal 
     24from Scientific.IO import NetCDF 
    2525 
    26 from twisted.internet import reactor, defer 
    27 from twisted.spread import pb 
    28 from asynqueue import ThreadQueue 
    29 from twisted_goodies.pybywire import pack, params 
    30  
    31 import project, tseries 
     26import project, tseries, params 
    3227import util 
    3328 
    …… 
    10297        # Make sure we can write to the open file 
    10398        var = cdf.variables['a'] 
    104         var[0,0,:] = s.rand(self.mgr.p).astype('f') 
    105         var[0,1,:] = s.rand(self.mgr.p).astype('f') 
    106         print var[0,:4,:] 
     99        x = s.rand(self.mgr.p, 2).astype('f') 
     100        var[0,0,:] = x[:,0] 
     101        var[0,1,:] = x[:,1] 
     102        self.failUnless(s.equal(s.transpose(var[0,:2,:]), x).all()) 
    107103        # So far, so good. Now check out the file 
    108104        cdf.close() 
    109         self.fail("TODO") 
     105        cdf = NetCDF.NetCDFFile(cdfPath, 'r') 
     106        self.failUnlessElementsEqual( 
     107            cdf.variables.keys(), 
     108            ['observations', 'log_volatility', 'log_likelihood', 
     109             'a', 'b', 'cr', 'd', 'e', 'fs', 'fr', 'g']) 
     110        var = cdf.variables['a'] 
     111        self.failUnless(s.equal(s.transpose(var[0,:2,:]), x).all()) 
     112        # All done 
     113        cdf.close() 
     114 
     115 
     116class Mock_ParameterContainer(util.Mock): 
     117    p, xcorrs = 3, 6 
     118    paramNames =  ['a', 'b', 'cr', 'd', 'e', 'fs', 'fr', 'g'] 
     119    paramShapes = [(p,), (p,p), (xcorrs,), (p,), (p,p), (p,), (xcorrs,), (p,)] 
     120 
     121    def __init__(self, scale): 
     122        for k, name in enumerate(self.paramNames): 
     123            array = k * scale * s.ones(self.paramShapes[k]) 
     124            setattr(self, name, array) 
     125        self.h = self.p * scale * s.ones(self.p) 
     126        self.Lx = 0.0 
     127         
     128 
     129class Test_ProjectManager_writeParams(util.TestCase): 
     130    def setUp(self): 
     131        projectDir = os.path.dirname(__file__) 
     132        self.cdfPath = os.path.join(projectDir, "svpmc.nc") 
     133        self.specFile = os.path.join(projectDir, "project-spec.txt") 
     134        self.mgr = project.ProjectManager(self.specFile, m=3) 
     135 
     136    def _makeParameters(self, *scales): 
     137        parameters = params.FlexArray(len(scales)) 
     138        for k, scale in enumerate(scales): 
     139            parameters[k] = Mock_ParameterContainer(scale) 
     140        return parameters 
     141         
     142    def test_writeParams(self): 
     143        # Write them 
     144        parameters = self._makeParameters(7, 11, 13) 
     145        self.mgr.writeParams(parameters) 
     146        self.mgr.cdf.close() 
     147        # Now read them 
     148        cdf = NetCDF.NetCDFFile(self.cdfPath, 'r') 
     149        for member, paramContainer in enumerate(parameters): 
     150            for name in self.mgr.paramNames: 
     151                array = cdf.variables[name][0, member, :] 
     152                arrayExpected = getattr(paramContainer, name) 
     153                self.failUnless(s.equal(array, arrayExpected).all()) 
     154            array = cdf.variables['log_volatility'][0, member, :] 
     155            self.failUnless(s.equal(array, paramContainer.h).all()) 
     156        L = cdf.variables['log_likelihood'][0,:] 
     157        self.failUnless(s.equal(L, parameters.Lx).all()) 
     158        # Make sure we are ready for the next iteration 
     159        self.failUnlessEqual(self.mgr.iteration, 1) 
     160 
     161             
     162