Source Code for Module gmisclib.mcmcS2

 1   
 2  """Markov-Chain Monte-Carlo algorithms. 
 3   
 4  Here, we do MCMC when -log(p) is a sum of squares.""" 
 5   
 6  import mcmcSQ 
 7  import die 
 8   
 9  import mcmc 
10  import g_implements 
11  import findleak 
12   
13  _s = None 
14   
15 -def go_step(x, showvec): 
16          global _s 
17   
18          findleak.vm('go_step start') 
19          assert g_implements.impl(x.current(), mcmc.position_base), "Bad pseudoposition: %s" % g_implements.why(x.current(), mcmc.position_base) 
20          if _s is None: 
21                  print "INITIALIZING mcmc" 
22                  _s = mcmc.bootstepper(None, x.current(), x.archive.mvn.cov, None) 
23                  findleak.vm('go_step initialize end') 
24          print "SQ STEP" 
25          accepted = x.step(showvec=showvec, atmisc={'stepper': 'SQ'}) 
26          findleak.vm('go_step accepted start') 
27          if accepted: 
28                  cur = x.current() 
29                  print 'Cur from SQ:', cur.prms() 
30                  _s._set_current( cur ) 
31                  findleak.vm('go_step accepted end') 
32          findleak.vm('go_step Bootstep start') 
33          if(len(_s.archive) > _s.np/_s.F): 
34                  print "BOOT STEP" 
35                  accepted = _s.step() 
36                  print 'Accepted:', accepted 
37                  cur = x.current() 
38                  print 'Cur from BOOT:', cur.prms() 
39                  SQpos = _s.current() 
40                  assert isinstance(SQpos, mcmcSQ.position) 
41                  x.archive.add( SQpos ) 
42                  if accepted: 
43                          findleak.vm('go_step BOOTstep-set_current start') 
44                          x._set_current( SQpos ) 
45                          findleak.vm('go_step BOOTstep-set_current end') 
46                  x.print_at( misc = {'stepper': 'Boot', 'BootScale': _s.a[0].vs()} ) 
47                  findleak.vm('go_step BOOTstep end') 
48          #findleak.print_top_N(20) 
49          findleak.vm('go_step end') 
50   
51  mcmcSQ.go_step = go_step 
52   
53   
54   
55   
56   
57  if __name__ == '__main__': 
58          import sys 
59          if len(sys.argv) <= 1: 
60                  print __doc__ 
61                  die.exit(0) 
62          mcmcSQ.go(sys.argv) 
63