qubx_about_plug.qubx

1 """Base classes for QUB Express plugins, and sample implementations. 2 3 This unit is not required to make a plugin (see qubx.pyenv), 4 but is intended to simplify and explain the process. 5 6 To make a plugin: 7 * make a folder in qub-express's Plugins folder (e.g. "my_plugin") 8 * make a subfolder with a unique name for your plugin's package (e.g. "my_plugin/my_uniquestr_package"); 9 the my_plugin folder will be part of sys.path; if there is any other my_uniquestr_package in sys.path, it 10 may be loaded instead. 11 * copy this file into my_uniquestr_package; by putting it there, you're guaranteed to import this copy, 12 and not one from some other plugin. 13 * make a python module, e.g. "my_uniquestr_package/my_module.py", which starts with "from . import qubx_plugin", 14 and defines a subclass of one of the plugin base classes. 15 * make a python module, "my_uniquestr_package/__init__.py" along the lines of: 16 __all__ = ['qubx_plugin', 'my_module'] 17 or if you want, leave it blank, but make sure it exists. 18 * make these text files in the my_plugin folder: 19 - qubx_plugin_about.txt - information for the Admin:Plugins panel 20 - qubx_plugin_imports.txt - newline-separated list of modules, in the order they should be re-loaded 21 - qubx_plugin_init.py - python script executed on load; 22 typically creates and adds one or more "Face" to the tool/about panels. 23 - qubx_plugin_fini.py - python script executed on unload; typically removes the "Face"(s) 24 - qubx_plugin_name.txt - name of plugin, as shown in Admin:Plugins panel 25 26 For examples, see qubx_sample_plugin, qubx_macroscopic, and qubx_single_molecule. 27 """ 28 29 import gobject 30 import gtk 31 import pango 32 import traceback 33 34 import qubx.notebook 35 import qubx.optimize 36 import qubx.optimizeGTK 37 import qubx.pyenv 38 import qubx.trialset 39 from qubx.accept import * 40 from qubx.faces import Face 41 from qubx.GTK import Requestable, pack_item, pack_space, pack_label, pack_scrolled, pack_check, pack_button 42 from qubx.model import QubModel, ModelToQ, ModelToP0, ModelToClazz, QtoPe, K012ToModel, ModelToRateMatricesP 43 from qubx.task import Robot, Tasks 44 from qubx.util_types import Anon, Product, Reffer, WeakEvent 45 46 from itertools import izip, count 47 from numpy import array, matrix, zeros, any, sum 48 import scipy.linalg.matfuncs 49 expm = scipy.linalg.matfuncs.expm 50

51 -class AboutPlugin(Face):

52 - def __init__(self, label, global_name=""):

53 super(AboutPlugin, self).__init__(label, global_name=global_name)

54

55 -class ToolPlugin(Face):

56 """Base class for top-pane plugin interface with optional worker thread (robot). 57 58 @ivar robot: L{qubx.task.Robot} managing a background work thread 59 """ 60 61 # these names are visible in the Explore dialog (Admin:Scripts): 62 __explore_featured = ['robot', 'on_robot_begin', 'on_robot_end', 'on_exception', 'on_interrupt'] 63

64 - def __init__(self, label, global_name="", has_robot=True):

65 """ 66 @param label: name on this tab in the gui 67 @param has_robot: if True, starts a background thread called self.robot 68 """ 69 super(ToolPlugin, self).__init__(label, global_name=global_name) 70 self.__ref = Reffer() 71 self.set_size_request(200, 50) 72 73 if has_robot: 74 self.robot = Robot(label, self.__ref(self.__on_robot_begin), self.__ref(self.__on_robot_end)) 75 self.robot.OnException += self.__ref(self.on_exception) 76 self.robot.OnInterrupt += self.__ref(self.on_interrupt) 77 QubX = qubx.pyenv.env.globals['QubX'] 78 QubX.OnQuit += self.__ref(self.__qubx_finalize) 79 else: 80 self.robot = None

81 - def __qubx_finalize(self):

82 self.robot.stop()

83 - def __on_robot_begin(self):

84 Tasks.add_task(self.robot) 85 self.on_robot_begin()

86 - def __on_robot_end(self):

87 Tasks.remove_task(self.robot) 88 self.on_robot_end()

89 - def on_robot_begin(self):

90 """Override this method to respond when the robot starts working.""" 91 pass

92 - def on_robot_end(self):

93 """Override this method to respond when the robot finishes working.""" 94 pass

95 - def on_exception(self, task, typ, val, tb):

96 """Override this method to do something with exceptions thrown on robot thread; default: print to stdout.""" 97 traceback.print_exception(typ, val, tb)

98 - def on_interrupt(self, task, cancel):

99 """Override this method to customize the effect of the Stop button (e.g. set a flag); default: KeyboardInterrupt raised in robot thread.""" 100 # if you don't call cancel, qubx will throw a KeyboardInterrupt into the robot thread. 101 pass

102

103 -class OptToolPlugin(ToolPlugin):

104 """Base class for optimization plugins. 105 106 Override at minimum: robot_calculate_score, robot_get_pars, robot_set_pars. 107 Optionally, override robot_optimize to replace the built-in optimizers. 108 109 They are called robot_* because they run on the worker (robot) thread. When called, 110 the latest model has been copied into self.robot.model, a L{QubModel}. 111 112 Add any additional widgets to self.top_panel and self.right_panel, or rearrange layout as needed. 113 114 model_prep has .source (the QubModel in GUI), .model (private copy for robot computation), 115 and .set_rates2(K0, K1, K2) [faster than making new model_prep] 116 117 @ivar right_panel: L{gtk.VBox} along right edge, with optimize controls already at the bottom 118 @ivar top_panel: L{gtk.HBox} along top edge, with score display already at far right 119 @ivar txtStatus: L{gtk.TextView} with text output; convenience methods self.append_output(s), self.clear_output() 120 @ivar model_prep: L{ModelPrep} with source model and private copy for robot use 121 122 """ 123 124 __explore_featured = ['right_panel', 'top_panel', 'txtStatus', 'model_prep', 'grads', 'hessian', 125 'OnLL', 'OnIteration', 'OnOptimized', 'model_tracker', 'data', 'data_prep', 126 'btnGetScore', 'txtScore', 'lblScore', 'btnAddToTrials', 'txtStatus', 127 'btnOpt', 'optimizing', 'score', 'model', 'get_model_prep', 'get_data_prep', 'request_data_prep', 128 'append_output', 'clear_output', 'add_to_trials', 'make_trial', 'get_ndata', 'add_trial_fields', 129 'calculate_score', 'optimize', 'robot_setup_data', 'update_score', 'robot_calculate_score', 130 'robot_get_pars', 'robot_set_pars', 'robot_optimize'] 131

132 - def __init__(self, label, global_name="", opt_name=None, opt_options={}, have_grads=False):

133 """ 134 @param label: tab label, becomes self.face_name 135 @param opt_name: optimizer preferences selector, defaults to label 136 @param opt_options: {option_label : default} to override specific defaults of the built-in optimizers (see qubx.optimize) 137 @param have_grads: True if your robot_calculate_score can also calculate dScore/dPar[i] 138 """ 139 super(OptToolPlugin, self).__init__(label, global_name=global_name) 140 self.__ref = Reffer() 141 self.__score = 0.0 142 self.__optimizing = False 143 self.have_grads = have_grads 144 self.opt_name = opt_name or label 145 qubx.optimize.setup_opt(self.opt_name, opt_options) 146 147 self.grads = [] 148 self.hessian = None 149 self.OnLL = WeakEvent() 150 self.OnIteration = WeakEvent() 151 self.OnOptimized = WeakEvent() 152 153 self.model_tracker = ModelTracker() 154 self.model_tracker.OnSwitch += self.__ref(self.__onSwitchModel) 155 self.model_tracker.OnChange += self.__ref(self.__onSetModel) 156 self.model_tracker.OnSetRate += self.__ref(self.__onSetRate) 157 self.model_prep = Product(self.get_model_prep) 158 self.model_prep.OnInvalidate += self.__ref(self.__onInvalidateModel) 159 self.rates = Product(self.__get_rates) 160 self.rates.OnInvalidate += self.__ref(self.__onInvalidateModel) 161 self.__model = None 162 self.data = OptToolData(self.get_data_prep) 163 self.data.OnInvalidate += self.__ref(self.__onInvalidateData) 164 self.data_prep = None 165 166 self.top_panel = pack_item(gtk.HBox(), self) 167 v = pack_item(gtk.VBox(), self.top_panel, at_end=True) 168 vh = pack_item(gtk.HBox(), v) 169 self.btnGetScore = pack_button('Get Score', vh, at_end=True, on_click=self.__onClickGetScore) 170 self.txtScore = pack_item(gtk.Entry(), vh, at_end=True) 171 self.txtScore.set_width_chars(10) 172 self.lblScore = pack_label('Score:', vh, at_end=True) 173 self.btnAddToTrials = pack_button('Add to Trials', v, expand=True, on_click=self.__onClickAddToTrials) 174 pack_space(self.top_panel, x=20, at_end=True) 175 176 box = pack_item(gtk.HBox(), self, expand=True) 177 self.txtStatus = gtk.TextView() 178 self.txtStatus.set_editable(False) 179 self.txtStatus.modify_font(pango.FontDescription('monospace 10')) 180 pack_scrolled(self.txtStatus, box, expand=True) 181 182 self.right_panel = pack_item(gtk.VBox(), box) 183 h = pack_item(gtk.HBox(True), self.right_panel, at_end=True) 184 self.btnOpt = pack_button('Optimize', h, on_click=self.__onClickOptimize) 185 186 self.optOptions = qubx.optimizeGTK.OptOptionsFace(label, self.opt_name, 'QubX.About.%s'%label) 187 # set self.optOptions to None if you're not using qubx.optimize 188 189 try: 190 qubx.global_namespace.QubX.Modeling.Utils.add_optimizer(label, self.get_model_prep, self.get_data_prep, 191 self.calculate_score, self.optimize, 192 self.OnLL, self.OnIteration, self.OnOptimized) 193 except: 194 traceback.print_exc() 195 196 self.__onSwitchModel(qubx.global_namespace.QubX.Models.file)

197 optimizing = property(lambda self: self.__optimizing)

198 - def set_score(self, x):

199 self.txtScore.set_text("%.5g" % x) 200 self.__score = x

201 score = property(lambda self: self.__score, lambda self, x: self.set_score(x), 202 "quantity to maximize, e.g. log-likelihoo; auto-displayed at top-right")

203 - def onShow(self, showing):

204 # if using built-in optimizers, show their options in the About panel 205 if not self.optOptions: 206 return 207 QubX = qubx.pyenv.env.globals['QubX'] 208 if showing: 209 QubX.About.append_face(self.optOptions) 210 QubX.About.show_face(self.optOptions.face_name) 211 else: 212 QubX.About.remove_face(QubX.About.index(self.optOptions.face_name))

213 - def __onSwitchModel(self, model):

214 if self.__optimizing: # put it off until optimization finishes 215 self.robot.do(gobject.idle_add, self.__onSwitchModel, model) 216 else: 217 self.__model = model 218 self.__onSetModel()

219 model = property(lambda self: self.__model)

220 - def __onInvalidateModel(self, product):

221 pass

222 #if not self.__optimizing: 223 # self.robot.do(gobject.idle_add, self.invalidate_controls)

224 - def __onSetModel(self, *args):

225 self.model_prep.invalidate() 226 self.rates.invalidate()

227 - def __onSetRate(self, r, field, val, prev, undoing):

228 self.rates.invalidate()

229 - def get_model_prep(self, model=None):

230 # override to customize the ModelPrep 231 m = model or self.__model 232 if not m: 233 return None 234 return ModelPrep(m)

235 - def __onInvalidateData(self, product):

236 self.model_prep.invalidate() 237 self.rates.invalidate()

238 - def get_data_prep(self, segs, model_prep, wait=True, receiver=lambda dp: None):

239 # leave this one alone; override the next one to customize DataPrep 240 if wait: 241 return qubx.pyenv.call_async_wait(self.request_data_prep, segs, model_prep, last_data_prep=self.data_prep) 242 else: 243 self.request_data_prep(receiver, segs, model_prep, last_data_prep=self.data_prep)

244 - def request_data_prep(self, receiver, segs=None, model_prep=None, last_data_prep=None):

245 # override to customize the DataPrep 246 data_prep = Anon(segs=segs, model_prep=model_prep) 247 receiver(data_prep)

248 - def __get_rates(self):

249 if self.__model == None: 250 return None, None 251 signal_names = qubx.global_namespace.QubX.Modeling.Stimulus.active_names 252 K0, K1, K2, L, V, P = ModelToRateMatricesP(self.__model, signal_names) 253 return K0, K1, K2

254 - def on_robot_begin(self):

255 def do_on_idle(): 256 self.btnGetScore.set_sensitive(False) 257 self.btnOpt.set_sensitive(False)

258 gobject.idle_add(do_on_idle)

259 - def on_robot_end(self):

260 def do_on_idle(): 261 self.btnGetScore.set_sensitive(True) 262 self.btnOpt.set_sensitive(True)

263 gobject.idle_add(do_on_idle)

264 - def __onClickGetScore(self, btn):

265 qubx.pyenv.env.scriptable_if_matching('%s.calculate_score()' % self.global_name, 266 [(self.global_name, self)]) 267 self.calculate_score(wait=False, receiver=self.set_score)

268 - def __onClickAddToTrials(self, btn):

269 self.add_to_trials(*self.make_trial())

270 - def __onClickOptimize(self, btn):

271 qubx.pyenv.env.scriptable_if_matching('%s.optimize(on_trial=%s.add_to_trials)' % (self.global_name, self.global_name), 272 [(self.global_name, self)]) 273 self.__mtree_pre_opt = self.__model.as_tree() 274 self.optimize(wait=False, on_trial=self.add_to_trials)

275 - def append_output(self, s):

276 """Appends a string to the output panel.""" 277 buf = self.txtStatus.get_buffer() 278 buf.insert(buf.get_end_iter(), s) 279 self.txtStatus.scroll_to_mark(buf.get_insert(), 0)

280 - def clear_output(self):

281 """Clears the text output panel.""" 282 self.txtStatus.set_text("")

283 - def on_exception(self, task, typ, val, tb):

284 self.append_output('\n'.join(traceback.format_exception(typ, val, tb)))

285 - def add_to_trials(self, mtree, row):

286 qubx.pyenv.env.OnScriptable('%s.add_to_trials(*%s.make_trial())' % (self.global_name, self.global_name)) 287 Trials = qubx.global_namespace.QubX.Trials 288 ts = Trials.get_trialset(self.face_name, in_qsf=True) 289 ts.append_trial(mtree, row) 290 Trials.trialset = ts

291 - def make_trial(self, model_prep=None, data_prep=None, fields={}, initial_tree=None):

292 try: 293 QubX = qubx.global_namespace.QubX 294 if model_prep: 295 model = model_prep.source or model_prep.model 296 else: 297 model = self.QubX.Models.file 298 if data_prep and data_prep.segs: 299 data = data_prep.segs[0].file 300 else: 301 data = self.QubX.Data.file 302 mtree = model.as_tree() 303 nParam = fields['nParam'] if ('nParam' in fields) else 0 304 if 'nData' in fields: 305 nData = fields['nData'] 306 else: 307 nData = self.get_ndata(data_prep) 308 row = qubx.trialset.make_trial_row(model.path, data.path, self.score, nParam, nData=nData) 309 row['nData'] = nData 310 if initial_tree: 311 self.add_trial_fields(row, initial_tree, 'Initial ', model_prep, data_prep) 312 self.add_trial_fields(row, mtree, "", model_prep, data_prep) 313 for k,v in fields.iteritems(): 314 row[k] = v 315 return mtree, row 316 except: 317 traceback.print_exc()

318 - def get_ndata(self, data_prep=None):

319 """Override to tell how many datapoints/events used, for Trials.""" 320 return 0

321 - def add_trial_fields(self, row, mtree, prefix='', model_prep=None, data_prep=None):

322 """Override to pick info from the model.""" 323 pass # e.g. qubx.trialset.add_trial_rates(row, mtree, prefix)

324 - def calculate_score(self, model_prep=None, data_prep=None, wait=True, receiver=lambda score: None, on_trial=None, **kw):

325 """Starts the robot calculating the score. 326 @param receiver: func(score) called when score is finished calculating 327 """ 328 if wait: 329 qubx.pyenv.call_async_wait(self.__calculate_score, model_prep, data_prep, on_trial=on_trial, **kw) 330 receiver(self.__score) 331 return self.__score 332 else: 333 self.__calculate_score(receiver, model_prep, data_prep, on_trial=on_trial, **kw)

334 - def __calculate_score(self, receiver=lambda score: None, model_prep=None, data_prep=None, on_trial=None, **kw):

335 QubX = qubx.pyenv.env.globals['QubX'] 336 self.__set_keys(**kw) 337 mp = model_prep 338 if mp is None: 339 mp = self.model_prep.val 340 mp.set_rates2(*self.rates.val) 341 self.robot.do(self.__robot_calculate_score, receiver, mp, data_prep, on_trial)

342 - def optimize(self, model_prep=None, data_prep=None, wait=True, receiver=None, ext_return=False, on_trial=None, **kw):

343 """Starts the robot optimizing. 344 @param receiver: func() called when optimization is finished 345 """ 346 if wait: 347 qubx.pyenv.call_async_wait(self.__optimize, model_prep, data_prep, on_trial=on_trial, **kw) 348 if ext_return: 349 if receiver: 350 receiver(self.score, self.iterations, self.grads, self.hessian) 351 return self.score, self.iterations, self.grads, self.hessian 352 else: 353 if receiver: 354 receiver(self.score) 355 return self.score 356 else: 357 self.__optimize(receiver, model_prep, data_prep, ext_return, on_trial=on_trial, **kw)

358 - def __optimize(self, receiver=None, model_prep=None, data_prep=None, ext_return=False, on_trial=None, **kw):

359 QubX = qubx.pyenv.env.globals['QubX'] 360 self.__set_keys(**kw) 361 mp = model_prep 362 if mp is None: 363 mp = self.model_prep.val 364 if mp.source: 365 mp.source.undoStack.enable(False) 366 mp.set_rates2(*self.rates.val) 367 self.__mtree_pre_opt = mp.model.as_tree() 368 QubX.Trials.get_trialset("%s Iterations"%self.face_name, in_qsf=True).clear() 369 self.__optimizing = True 370 self.robot.do(self.__robot_optimize, receiver, mp, data_prep, ext_return, on_trial)

371 - def __set_keys(self, **kw):

372 for k, v in kw.iteritems(): 373 try: 374 setattr(self, k, v) 375 except: 376 print 'MacRates optimizer has no option named %s' % k

377 - def robot_setup_data(self, model_prep, data_prep):

378 self.robot.model_prep = model_prep 379 if data_prep: 380 self.data_prep = data_prep 381 else: 382 self.data.model_prep = model_prep 383 self.data_prep = self.robot.gui_call_recv(self.data.request, self.robot)[0]

384 - def __robot_calculate_score(self, on_finish=lambda: None, model_prep=None, data_prep=None, on_trial=None):

385 """Calculates score on worker thread, then passes it back to GUI thread.""" 386 self.robot_setup_data(model_prep, data_prep) 387 try: 388 score = self.robot_calculate_score() 389 except: 390 score = 0.0 391 self.robot.send_exception() 392 if on_trial: 393 mtree, row = self.make_trial(model_prep, data_prep) 394 row['LL'] = score 395 gobject.idle_add(on_trial, mtree, row) 396 gobject.idle_add(self.update_score, score, lambda: on_finish(score)) 397 gobject.idle_add(self.OnLL, score)

398 - def __robot_optimize(self, on_finish=None, model_prep=None, data_prep=None, ext_return=False, on_trial=None):

399 """Optimizes on worker thread, then notifies GUI thread.""" 400 self.robot_setup_data(model_prep, data_prep) 401 mtree_initial = model_prep.model.as_tree() 402 try: 403 score, pars, iterations, grads, Hessian = self.robot_optimize() 404 self.iterations = iterations 405 self.grads = grads 406 self.hessian = Hessian 407 except: 408 self.robot.send_exception() 409 score = 0.0 410 self.score = score 411 self.__optimizing = False 412 if on_trial: 413 row = {'nParam' : len(self.grads), 414 'Iterations' : self.iterations} 415 for i, g in enumerate(self.grads): 416 row['Gradient %i' % i] = g 417 gobject.idle_add(on_trial, *self.make_trial(model_prep, data_prep, row, mtree_initial)) 418 if ext_return: 419 gobject.idle_add(self.update_score, score, lambda: on_finish and on_finish((score, iterations, grads, Hessian))) 420 else: 421 gobject.idle_add(self.update_score, score, lambda: on_finish and on_finish(score)) 422 gobject.idle_add(self.__onOptimized, iterations, score, grads, Hessian)

423 - def __onOptimized(self, iterations, score, grads, Hessian):

424 model_prep = self.robot.model_prep 425 if model_prep.source: 426 model_prep.source.undoStack.enable() 427 model_prep.source.undoStack.push_undo(bind(model_prep.source.from_tree, self.__mtree_pre_opt), 428 bind(model_prep.source.from_tree, model_prep.source.as_tree())) 429 model_prep.source.undoStack.seal_undo('set rates') 430 self.OnOptimized(iterations, score, grads, Hessian)

431 - def update_score(self, score, call_me=lambda: None):

432 """Updates score and txtScore then calls call_me(). Intended as GUI callback from robot thread; e.g. 433 >>> gobject.idle_add(self.update_score, score, on_finish) 434 """ 435 self.score = score 436 call_me()

437 - def robot_calculate_score(self, grads=None):

438 """Calculates and returns the score, on worker thread; override this method. 439 self.model_prep.model has just been updated. 440 @param grads: either None (no gradient requested) or an array to fill with dScore/dPar[i] 441 """ 442 if grads != None: 443 grads[0] = 0.0 444 grads[1] = 0.0 445 return 1.0

446 - def robot_get_pars(self):

447 """Returns the list of initial parameter values, on worker thread; override this method. 448 Derive parameter values from self.robot.model_prep.model and/or from your own widgets.""" 449 return [1.0, -1.0]

450 - def robot_set_pars(self, pars):

451 """Updates self.robot.model_prep.model and/or other state from new param values, on worker thread; override this method.""" 452 pass

453 - def robot_optimize(self):

454 """Finds improved param values, starting from robot_get_pars(); override if you don't want the built-in optimizers. 455 If you override, also do self.optOptions = None. 456 self.robot.model_prep.model has just been updated. 457 """ 458 pars = self.robot_get_pars() 459 self.robot.iter = 0 460 self.robot.maxiter = qubx.optimize.get_options(self.opt_name)['max iterations'] 461 if self.optOptions: 462 optimize = qubx.optimize.get_opts(self.opt_name, 463 on_start=lambda rec, opts: gobject.idle_add(self.optOptions.optimizers.onStartOpt, rec, opts), 464 on_end=lambda rec, opts: gobject.idle_add(self.optOptions.optimizers.onEndOpt, rec, opts)) 465 else: 466 optimize = qubx.optimize.get_opts(self.opt_name) 467 score, pars, iterations, grads, Hessian = optimize(array(pars), self.__robot_func, self.__robot_iter, self.have_grads and self.__robot_grad_func or None) 468 # important: use idle_add to pass GUI-changing stuff to the main thread: 469 gobject.idle_add(self.append_output, "Optimization finished after %d iterations.\nScore = %.9g\n" % (iterations, -score)) 470 # (on linux, it can be any thread, as long as the main thread is paused [with robot.main_hold:], 471 # but windows will crash unless all GUI API calls come from one thread.) 472 473 return -score, pars, iterations, grads, Hessian

474 # callbacks for built-in optimizers:

475 - def __robot_func(self, pars):

476 self.robot_set_pars(pars) 477 return - self.robot_calculate_score()

478 - def __robot_grad_func(self, pars, grads):

479 self.robot_set_pars(pars) 480 return - self.robot_calculate_score(grads)

481 - def __robot_iter(self, pars, grads, iterations, score):

482 self.robot.iter = iterations 483 self.robot.progress = self.robot.iter * 100.0 / self.robot.maxiter 484 gobject.idle_add(self.update_score, -score) 485 gobject.idle_add(self.OnIteration, self.robot.iter, score, grads) 486 self.robot_on_iteration(pars, grads, iterations, score) 487 mtree, row = self.make_trial(model_prep=Anon(model=self.robot.model_prep.model), data_prep=self.data_prep, fields={"Iteration": iterations}) 488 gobject.idle_add(self.__add_iter_trial, mtree, row)

489 - def robot_on_iteration(self, pars, grads, iterations, score):

490 """Responds to a completed iteration; override this function to send model params etc back to gui.""" 491 pass

492 - def __add_iter_trial(self, mtree, row):

493 Trials = qubx.global_namespace.QubX.Trials 494 Trials.get_trialset("%s Iterations"%self.face_name, in_qsf=True).append_trial(mtree, row)

495 496

497 -class ModelPrep(Anon):

498 """Represents a model ready for computation with QubX.Modeling.Utils."""

499 - def __init__(self, model=None, modeltree=None): # one or the other, not both; modeltree is more efficient

500 Anon.__init__(self, source=model, modeltree=modeltree or model.as_tree(), model=QubModel()) 501 self.model.from_tree(self.modeltree) 502 self.K0, self.K1, self.K2, L, V, P = ModelToRateMatricesP(self.model)

503 - def set_rates2(self, K0, K1, K2):

504 self.K0, self.K1, self.K2 = K0, K1, K2 505 K012ToModel(K0, K1, K2, self.model)

506 507

508 -class OptToolData(Requestable):

509 """Caches last-used data until changed."""

510 - def __init__(self, get_data_prep):

511 Requestable.__init__(self) 512 self.get_data_prep = get_data_prep 513 self.__ref = Reffer() 514 self.__use_group = None 515 self.__data = None 516 QubX = qubx.global_namespace.QubX 517 QubX.DataSource.OnChangeIdealization += self.__ref(self.__onChangeIdl) 518 QubX.Data.table.OnInsert += self.__ref(self.__onChangeFileList) 519 QubX.Data.table.OnRemoved += self.__ref(self.__onChangeFileList) 520 QubX.Data.table.OnSet += self.__ref(self.__onSetFileList) 521 QubX.Data.OnSwitch += self.__ref(self.__onSwitchDataFile) 522 self.__datafile = None 523 self.__onSwitchDataFile(QubX.Data, QubX.Data.file) 524 QubX.Modeling.Stimulus.OnInsert += self.__ref(self.__onChangeStimulus) 525 QubX.Modeling.Stimulus.OnChange += self.__ref(self.__onChangeStimulus) 526 QubX.Modeling.Stimulus.OnRemoving += self.__ref(self.__onChangeStimulus) 527 self.invalidate()

528 - def set_use_group(self, x):

529 if self.__use_group != x: 530 self.__use_group = x 531 self.invalidate()

532 use_group = property(lambda self: self.__use_group, lambda self, x: self.set_use_group(x))

533 - def __onChangeIdl(self, *args):

534 self.invalidate()

535 - def __onChangeFileList(self, *args):

536 if self.__use_group != None: 537 self.invalidate()

538 - def __onSetFileList(self, i, field_name, val, prev, undoing):

539 if (self.__use_group != None) and (field_name == 'Group'): 540 self.invalidate()

541 - def __onChangeStimulus(self, *args):

542 self.invalidate()

543 - def __onSwitchDataFile(self, Data, file):

544 if file == self.__datafile: return 545 if self.__datafile: 546 self.__datafile.constants.OnSet -= self.__ref(self.__onSetConstant) 547 self.__datafile = file 548 if self.__datafile: 549 self.__datafile.constants.OnSet += self.__ref(self.__onSetConstant)

550 - def __onSetConstant(self, r, field, val, prev, undoing):

551 if self.__datafile.constants[r,'Name'] == 'Tdead (ms)': 552 self.invalidate()

553 - def get(self, receiver):

554 QubX = qubx.global_namespace.QubX 555 if self.__use_group != None: 556 view_segm = QubX.DataSource.get_segmentation_files(group=self.__use_group) 557 else: 558 view_segm = [(QubX.Data.view, QubX.DataSource.get_segmentation())] 559 segs = [] 560 for view, segm in view_segm: 561 segs.extend(segm) 562 self.get_data_prep(segs, self.model_prep, wait=False, receiver=receiver)

563 564

565 -class AboutPluginExample(AboutPlugin):

566 __explore_featured = ['lblModelStatus', 'modelTracker']

567 - def __init__(self, name='Example', global_name='QubX.About.Example'):

568 super(AboutPluginExample, self).__init__(name, global_name=global_name) 569 570 # qubx uses WeakEvent for synchronous notification of listeners. 571 # WeakEvent uses weak references, so it's not a big deal if a listener never unsubscribes, 572 # but this means the references must be held somewhere else; 573 # it requires special handling for lambda expressions, as well as bound methods, 574 # which are ordinarily created and destroyed as needed. Thus the following weak 575 # reference becomes invalid by the end of the line: 576 # QubX.OnSomeEvent += self.__onSomeEvent 577 # We keep a permanent reference to event handlers in a Reffer object. 578 # self.__ref(x) == x 579 # and the act of calling self.__ref(x) adds it to the Reffer's internal dict (as a strong reference). 580 # Thus the idiom: 581 # QubX.OnSomeEvent += self.__ref(self.__onSomeEvent) 582 self.__ref = Reffer() 583 584 # add components to self (it is a subclass of gtk.VBox): 585 self.lblModelStatus = pack_label('', self, expand=True) 586 self.lblDataStatus = pack_label('', self, expand=True) 587 # they will show updated info about data and model, as they change: 588 self.modelTracker = ModelTracker() # defined below 589 self.modelTracker.OnSwitch += self.__ref(self.__onSwitchModel) 590 self.modelTracker.OnChange += self.__ref(self.__onChangeModel) 591 self.modelTracker.OnSetRate += self.__ref(self.__onSetModelRate) 592 self.__onChangeModel(self.modelTracker.model) # initialize 593 # QubX is the main window. DataSource represents the active portion of data 594 # (frontmost file, screen or whole file, which signal) 595 QubX = qubx.pyenv.env.globals['QubX'] 596 self.DataSource = QubX.DataSource 597 QubX.DataSource.OnChangeSamples += self.__ref(self.__onChangeSamples) 598 QubX.DataSource.OnChangeIdealization += self.__ref(self.__onChangeIdealization) 599 self.__onChangeData() # initialize

600 - def __onSwitchModel(self, models, model):

601 self.__onChangeModel(model)

602 - def __onSetModelRate(self, r, field, val, prev, undoing):

603 self.__onChangeModel(self.__model)

604 - def __onChangeModel(self, model):

605 self.lblModelStatus.set_text("%d states" % model.states.size)

606 - def __onChangeSamples(self):

607 self.__onChangeData()

608 - def __onChangeIdealization(self):

609 self.__onChangeData()

610 - def __onChangeData(self):

611 segs = self.DataSource.get_segmentation() 612 nseg = len(segs) 613 npoint = 0 614 nideal = 0 615 for seg in segs: 616 for chunk in seg.chunks: 617 if chunk.included: 618 npoint += chunk.n 619 ff, ll, cc = chunk.get_idealization() 620 nideal += len(ff) 621 self.lblDataStatus.set_text("%d points in %d segments\n(%d events)" % (npoint, nseg, nideal))

622 623

624 -class ToolPluginExample(ToolPlugin):

625 - def __init__(self, name='Test', global_name='QubX.Modeling.Test'):

626 super(ToolPluginExample, self).__init__(name, global_name=global_name) 627 # unused: self.__ref = Reffer() 628 629 self.lblMsg = pack_label('This plugin finds the max and min value in the Data Source.', self, expand=True) 630 self.panButtons = pack_item(gtk.HBox(), self) 631 self.btnRun = pack_button('Run', self.panButtons, at_end=True, on_click=self.__onClickRun)

632 - def on_robot_begin(self):

633 gobject.idle_add(self.btnRun.set_sensitive, False)

634 - def on_robot_end(self):

635 gobject.idle_add(self.btnRun.set_sensitive, True)

636 - def __onClickRun(self, btn):

637 # calls self.__robot_compute() on the robot (worker) thread: 638 self.robot.do(self.__robot_compute)

639 - def __robot_compute(self):

640 QubX = qubx.pyenv.env.globals['QubX'] 641 with self.robot.main_hold: # pause the main thread to read/write volatile info: 642 units = QubX.Data.file.signals.get(QubX.DataSource.signal, 'Units') 643 lo = 1e20 644 hi = -lo 645 segments = QubX.DataSource.get_segmentation() 646 for s, seg in enumerate(segments): 647 # A segment is divided in chunks, alternating included/excluded. 648 # The chunks don't have sampled data yet. 649 # gen_samples yields sampled chunks, and also subdivides any chunks bigger than maxlen. 650 # It also takes care of pausing the main thread to read data. 651 for chunk in QubX.DataSource.gen_samples(seg.chunks, self.robot.main_hold, maxlen=(1<<20)): 652 if chunk.included: 653 lo = min(lo, min(chunk.samples)) 654 hi = max(hi, max(chunk.samples)) 655 self.robot.progress = s * 100.0 / len(segments) 656 # for short chunks, gen_samples is equivalent to: 657 # for chunk in seg.chunks: 658 # if chunk.included: 659 # with self.robot.main_hold: 660 # datachunk = chunk.get_samples() 661 # lo = min(lo, min(datachunk.samples)) 662 # hi = max(hi, max(datachunk.samples)) 663 664 # gobject.idle_add(a, b, c, ...) calls a(b, c, ...) on the main thread, sometime later. 665 # It's good for displaying intermediate and final results, messages, etc. which must 666 # be done on the main thread for Win32 compatibility. 667 gobject.idle_add(self.__show_results, 'Lo: %.6g %s\nHi: %.6g %s' % (lo, units, hi, units))

668 - def __show_results(self, results):

669 self.lblMsg.set_text(results)

670 671

672 -class OptToolPluginExample(OptToolPlugin):

673 """Sample L{OptToolPlugin} which optimizes rate constants given idealized single-molecule data. 674 The algorithm is a simplification of MIL (Qin et al). 675 """ 676 __explore_featured = ['chkPeq']

677 - def __init__(self, label='Example', global_name='QubX.Modeling.Example'):

678 super(OptToolPluginExample, self).__init__(label, global_name=global_name, have_grads=False) # it will auto-numerical-differentiate when required 679 self.__ref = Reffer() 680 681 # add a checkbox: equilibrium or fixed entry probability 682 self.chkPeq = pack_check('start at equilibrium (else use States:P0)', self.top_panel, active=True) 683 684 self.append_output('This plugin demonstrates log-likelihood optimization of rate constants from idealized data, in just 85 lines.\nThe algorithm is a simplification of MIL (Qin et al, 1996).\n')

685

686 - def robot_get_pars(self):

687 # this plugin doesn't handle stimulus-dependent rates, just k0 688 model = self.robot.model_prep.model 689 for i in xrange(model.rates.size): 690 if (model.rates.get(i, 'k1') or model.rates.get(i, 'k2') 691 or model.rates.get(i, 'Ligand') 692 or model.rates.get(i, 'Voltage') or model.rates.get(i, 'Pressure')): 693 gobject.idle_add(self.append_output, "Warning: Ignoring Ligand- and/or Voltage-dependence.\n") 694 # the params are log(k0), to constrain k0 > 0 695 return [log(model.rates.get(i, 'k0')) for i in xrange(model.rates.size)]

696 - def robot_set_pars(self, pars):

697 model = self.robot.model_prep.model 698 k0 = [exp(pars[i]) for i in xrange(model.rates.size)] 699 # update rates in the worker-thread's model: 700 for i in xrange(model.rates.size): 701 model.rates.set(i, 'k0', k0[i]) 702 # and in the visible one, on the main thread: 703 gobject.idle_add(self.__update_gui_rates, k0)

704 - def __update_gui_rates(self, k0):

705 if self.robot.model_prep.source: 706 model = self.robot.model_prep.source 707 for i in xrange(model.rates.size): 708 model.rates.set(i, 'k0', k0[i])

709

710 - def request_data_prep(self, receiver, segs=None, model_prep=None, last_data_prep=None):

711 # this one runs in the gobject thread, as needed (cached by OptToolPlugin) 712 idl = [] 713 nevent = 0 714 dt = segs[0].sampling if segs else 1e-4 715 for seg in segs: 716 # first, last, class, duration -- we need only class and duration arrays 717 # i.e. a segment is a sequence of dwells (events): 718 # - in a certain class (color; collection of indistinguishable states) 719 # - for some number of seconds (by assumption, durations are even multiples of dt) 720 ff, ll, cc, dd = seg.get_idealization(mark_excluded=True, get_durations=True) 721 # "excluded" and un-idealized portions will have class < 0 722 dd *= 1e-3 # convert from ms to sec 723 idl.append( (cc, dd) ) 724 if any(cc >= 0): # don't count segments with un-idealized or all-excluded data 725 nevent += len(cc) 726 receiver(Anon(segs=segs, model_prep=model_prep, idl=idl, nevent=nevent, dt=dt))

727

728 - def robot_calculate_score(self, grads=None):

729 QubX = qubx.pyenv.env.globals['QubX'] 730 idl = self.data_prep.idl 731 nevent = self.data_prep.nevent 732 dt = self.data_prep.dt 733 if nevent == 0: 734 gobject.idle_add(self.append_output, "No idealized data available (check the data source, or idealize).\n") 735 return 0.0 736 # read widgets while the main thread is paused: 737 with self.robot.main_hold: 738 usePeq = self.chkPeq.get_active() 739 # get rate matrix (Q) and state colors (class) 740 model = self.robot.model_prep.model 741 Q = ModelToQ(model) 742 N = Q.shape[0] 743 class_of_state = ModelToClazz(model) 744 # get equilibrium or fixed entry prob 745 if usePeq: 746 P0 = QtoPe(Q) 747 else: 748 P0 = ModelToP0(model) 749 750 # math subroutines for transition probability and log likelihood accumulation: 751 752 # Prob(state b at time t+dt |given| state a at time t) = A[a,b] = e^(Q*dt)[a,b] 753 A_base = expm(Q*dt) 754 # except we zero out columns which are ruled out by data (wrong class [color]), in this work matrix: 755 A_work = matrix(zeros(shape=(N,N), dtype='float64')) 756 # cache exponentiated rate matrices A(dt*nsample), since there will be a limited number of multiples of dt 757 A_cache = {} 758 def get_A(dur, dest_class): 759 """Returns A'^(dur/dt), where A' = e^(Q*dt), except with A'[r,c]==0.0 if c!=dest_class; 760 caches all answers in A_cache and returns answers from cache when possible.""" 761 if not ((dur, dest_class) in A_cache): 762 # start with one-sample transition prob. matrix: 763 A_work[:] = A_base[:] 764 # if data specifies destination class, then zero out other columns 765 # (once per sample, it is observed in dest_class, so prob(going anywhere else in one sample) is zero) 766 if dest_class >= 0: 767 A_work[:,class_of_state != c] = 0.0 768 # exponentiate to repeat by number of samples 769 A_cache[(dur, dest_class)] = A_work ** int(round(dur/dt)) 770 return A_cache[(dur, dest_class)]

771 772 # accumulating likelihood quickly results in underflow; following Qin et al (1996), 773 # we accumulate -log(scale) at each dwell, and rescale Pt 774 def Scale(Pt): 775 """Returns scale=1.0/sum(Pt); multiplies each component of Pt by scale.""" 776 mag = sum(Pt) 777 if mag == 0.0: 778 raise Exception('Impossible!') 779 scale = 1.0 / mag 780 Pt *= scale 781 return scale

782 783 # main loop: 784 LL = 0.0 785 for classes, durations in idl: # each segment 786 for i, c, d in izip(count(), classes, durations): # each dwell: index, class, duration 787 # transition probability into class c: 788 if i == 0: # first dwell of segment: Pt[i] = P0[i], if state i is in class c 789 Pt = matrix(zeros(shape=(1,N), dtype='float64')) 790 for i in xrange(N): 791 Pt[0,i] = (class_of_state[i] == c) and P0[i] or 0.0 792 else: # apply transition matrix for one sample duration (rough approximation) 793 transition_dur = min(dt, durations[i-1]) 794 Pt *= get_A(transition_dur, c) 795 # transition probability for dwell duration (minus transition of dt): 796 dwell_dur = d - dt 797 if dwell_dur > 0.0: 798 Pt *= get_A(dwell_dur, c) 799 # rescale Pt toward 1.0, and accumulate log likelihood: 800 LL -= log(Scale(Pt)) 801 return LL 802 803 804

805 -class ModelTracker(object):

806 """Follows the frontmost L{qubx.model.QubModel} and L{qubx.modelGTK.ModelView}, 807 with events when anything changes. 808 809 @ivar OnSwitch: L{WeakEvent}(QubX.Models, model) when a new file is brought to front 810 @ivar OnChange: L{WeakEvent}(model) when something in the model changes (except a rate constant) 811 @ivar OnSetRate: L{WeakEvent}(rate_index, field_name, value, prev, is_undoing) 812 """

813 - def __init__(self):

814 self.__ref = Reffer() 815 QubX = qubx.pyenv.env.globals['QubX'] 816 self.Models = QubX.Models 817 self.Models.OnSwitch += self.__ref(self.__onSwitchModel) 818 self.__model = self.__view = None 819 self.OnSwitch = WeakEvent() 820 self.OnChange = WeakEvent() 821 self.OnSetRate = WeakEvent() 822 self.__onSwitchModel(self.Models, self.Models.file)

823 - def __onSwitchModel(self, models, model):

824 if model == self.__model: return 825 if self.__model: 826 for event in (self.__model.OnSetChannelCount, 827 self.__model.OnSetIeqFv, 828 self.__model.OnSetVRev, 829 self.__model.states.OnSet, 830 self.__model.states.OnInsert, 831 self.__model.states.OnRemoved, 832 self.__model.classes.OnSet, 833 self.__model.rates.OnInsert, 834 self.__model.rates.OnRemoved, 835 self.__model.OnChangeBalanceLoops, 836 self.__model.constraints_kin.OnEdit): 837 event -= self.__ref(self.__onChangeModel) 838 self.__model.rates.OnSet -= self.__ref(self.__onSetRate) 839 self.__model = model 840 self.__view = models.view 841 if self.__model: 842 for event in (self.__model.OnSetChannelCount, 843 self.__model.OnSetIeqFv, 844 self.__model.OnSetVRev, 845 self.__model.states.OnSet, 846 self.__model.states.OnInsert, 847 self.__model.states.OnRemoved, 848 self.__model.classes.OnSet, 849 self.__model.rates.OnInsert, 850 self.__model.rates.OnRemoved, 851 self.__model.OnChangeBalanceLoops, 852 self.__model.constraints_kin.OnEdit): 853 event += self.__ref(self.__onChangeModel) 854 self.__model.rates.OnSet += self.__ref(self.__onSetRate) 855 self.OnSwitch(model)

856 model = property(lambda self: self.__model) 857 view = property(lambda self: self.__view)

858 - def __onChangeModel(self, *args):

859 self.OnChange(self.__model)

860 - def __onSetRate(self, r, field, val, prev, undoing):

861 self.OnSetRate(r, field, val, prev, undoing)

862

Source Code for Module qubx_about_plug.qubx_plugin