Source Code for Module qubx_mac.rates

  1  """Provides a visual panel (L{qubx.faces.Face}) to optimize rate constants from ensemble measurement. 
  2   
  3  Uses math from U{qubx_macroscopic<http://www.qub.buffalo.edu/src/qub-express/api/qubx_macroscopic>} 
  4   
  5  /* Copyright 2008-2014 Research Foundation State University of New York   */ 
  6   
  7  /* This file is part of QUB Express.                                      */ 
  8   
  9  /* QUB Express is free software; you can redistribute it and/or modify    */ 
 10  /* it under the terms of the GNU General Public License as published by   */ 
 11  /* the Free Software Foundation, either version 3 of the License, or      */ 
 12  /* (at your option) any later version.                                    */ 
 13   
 14  /* QUB Express is distributed in the hope that it will be useful,         */ 
 15  /* but WITHOUT ANY WARRANTY; without even the implied warranty of         */ 
 16  /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          */ 
 17  /* GNU General Public License for more details.                           */ 
 18   
 19  /* You should have received a copy of the GNU General Public License,     */ 
 20  /* named LICENSE.txt, in the QUB Express program directory.  If not, see  */ 
 21  /* <http://www.gnu.org/licenses/>.                                        */ 
 22   
 23  """ 
 24   
 25  import cStringIO 
 26  import itertools 
 27  import gc 
 28  import operator 
 29  import os 
 30  import sys 
 31  import time 
 32  import traceback 
 33   
 34  from qubx.data_types import * 
 35  from qubx.util_types import * 
 36  from qubx.faces import * 
 37  from qubx.GTK import * 
 38  import qubx.fast.data 
 39  import qubx.fast.model 
 40  import qubx.notebook 
 41  import qubx.optimize 
 42  import qubx.optimizeGTK 
 43  import qubx.pyenv 
 44  import qubx.tree 
 45  from qubx.model import * 
 46  from qubx.task import * 
 47  from qubx.settings import Propertied, Property 
 48  from qubx.trialset import add_trial_rates, add_trial_amps, add_trial_nChannel 
 49  import max_mac_ll 
 50  from . import qubx_plugin 
 51   
 52  import scipy.optimize 
 53  from numpy import * 
 54   
 55  import ctypes 
 56  from ctypes import c_int, c_float, c_double, c_void_p, c_char_p, Structure, POINTER, CFUNCTYPE, byref, sizeof, py_object, pointer, cast 
 57  pybuf = ctypes.pythonapi.PyBuffer_FromReadWriteMemory 
 58  pybuf.restype = py_object 
 59   
 60   
 61  MIN_SECONDS_BETWEEN_FIT = 3 
 62  MAXITER = 300 
 63   
 64 -class DataPrep(Anon): 
 65      """Represents data ready for computation.""" 
 66 -    def __init__(self, segs, datas): 
 67          Anon.__init__(self, segs=segs, datas=datas) 
 68   
 69   
 70   
 71  #class MacData(Requestable): 
 72  #    """Tracks data source, prepares data and stimulus as needed for LL, opt.  Data is 
 73  #[L{max_mac_ll.Data}], each item decorated with additional fields: 
 74  #    - iVoltage: index of voltage signal 
 75  #    - signal_names: list of strings 
 76  #    - file: L{qubx.data_types.QubData} 
 77  #    - ff_ll_counts: [(firsts, lasts, count)] -- resampled interval bounds, per data chunk 
 78  # 
 79  #    @ivar resample: True if data should be resampled by delta next time 
 80  #    @ivar resample_std: the "delta" of resample-by-delta 
 81  #    @ivar use_group: a Group index in the Data table, to fit multi files, or None to use the Data Source 
 82  #""" 
 83   
 84  @Propertied(Property('use_Peq', True, "True to start each segment from equilibrium; False to use States:Pr"), 
 85              Property('use_group', False, "True to fit a group of data files", in_tree=False), 
 86              Property('use_group_index', 1, "Group of interest in Data table"), 
 87              Property('resample', False, "True to resample by delta y"), 
 88              Property('resample_delta', 0.5, 'an interval with std < delta is resampled as one data point'), 
 89              Property('Gleak', 0.0, 'Leakage conductance (pS)'), 
 90              Property('Vleak', 0.0, 'Leak voltage (mV)'), 
 91              Property('opt_Nchannel', False, "True to optimize the model's channel count"), 
 92              Property('opt_cond', False, "True to optimize the model's Amps/Conds"), 
 93              Property('opt_var', False, "True to optimize the model's Stds/CondStds"), 
 94              # not persistent: 
 95              Property('live', False, "True to redo fit curve when model changes", in_tree=False), 
 96              Property('use_GPU', False, "True to use OpenCL-capable graphics card, if avail", in_tree=False), 
 97              Property('multi_model', False, "True to fit the sum of a group of models", in_tree=False), 
 98              Property('multi_model_group', 1, "Group for multi_model fitting (column in Models table)") 
 99              ) 
100 -class RatesFace(qubx_plugin.OptToolPlugin): 
101      """Calculates log-likelihood and optimizes model parameters, a la Milescu 2005. 
102   
103  @ivar OnIteration: L{WeakEvent}(iteration, LL, grads, K0, K1, amp, std, amp0, std0, Nchannel) 
104  @ivar OnIteration2: L{WeakEvent}(iteration, LL, grads, K0, K1, K2, amp, std, amp0, std0, Nchannel) 
105  @ivar OnOptimized: L{WeakEvent}(iterations, LL, grads, Hessian) 
106  @ivar OnLiveLL: L{WeakEvent}(LL) 
107  """ 
108       
109      __explore_featured = ['accel', 'rates', 'iterations', 'LL', 'last_LL', 'grads', 'OnLiveLL', 'fitbot', 'get_ndata', 
110                            'add_trial_fields', 'wait_for_fits'] 
111       
112 -    def __init__(self, name='MacRates', global_name='QubX.Modeling.MacRates'): 
113          super(RatesFace, self).__init__(name, global_name, 
114                                          "qubx_mac_rates_optimize", 
115                                          {"max step" : 0.1, 
116                                           "num grad dx" : 1e-4, 
117                                           "conv ll" : 1e-3, 
118                                           "conv grad" : 1e-3}) 
119          self.__ref = Reffer() 
120          self.__invalid_score = True 
121          self.__invalid_live = True 
122          self.accel = 1 # (0: single; 1: openmp; 2: opencl) 
123          self.QubX = qubx.pyenv.env.globals['QubX'] 
124          self.QubX.OnQuit += self.__ref(self.mac_finalize) 
125          self.propertied_connect_settings('MacRates') 
126          self.data.use_group = self.use_group_index if self.use_group else None 
127          self.rates = Product(self.__get_rates) # overriding OptToolPlugin.rates 
128          self.rates.OnInvalidate += self.__ref(self.__onInvalidateRates) 
129          self.OnIteration = WeakEvent() 
130          self.OnIteration2 = WeakEvent() 
131          self.OnOptimized = WeakEvent() 
132   
133          self.robot.working = False 
134          self.robot.workspaces = None 
135          self.robot.optNchannel = self.robot.optCond = self.robot.optVar = False 
136          self.iterations = 0 
137          self.LL = self.last_LL = 0.0 
138          self.grads = [] 
139          self.__Ndata = 0 
140          self.OnLiveLL = WeakEvent() 
141   
142          # When calculation is fast enough, showing fit curves can be the slowest part of an iteration. 
143          # Ideally we'd skip displaying an iteration's curves if the optimizer has already finished the next iteration. 
144          # Showing curves on a separate thread makes it easier: 
145          self.fitbot = Robot('Apply Fit') # not in task list 
146          self.fitbot.OnException += self.__ref(self.on_exception) 
147          self.fitbot.serial = self.fitbot.last_serial = 0 
148          # We also skip iterations when too few seconds have elapsed: 
149          self.__last_fits_time = datetime.datetime.now() 
150          # rewind the clock so fit curves can show right away (there must be a simpler way): 
151          if self.__last_fits_time.hour: 
152              self.__last_fits_time = self.__last_fits_time.replace(hour=self.__last_fits_time.hour-1) 
153          # else: it will not show any fits for the first few seconds oh well 
154   
155          # these figures will be published through Admin:Notebook 
156          qubx.notebook.Notebook.register_auto('MacRates_eval.Messages', 'Messages, on MacRates: Get LL', False) 
157          qubx.notebook.Notebook.register_auto('MacRates_eval.Model.Picture', 'Model picture, on MacRates: Get LL', False) 
158          qubx.notebook.Notebook.register_auto('MacRates_eval.Model.PictureNO', 'Model picture, no overlays, on MacRates: Get LL', False) 
159          qubx.notebook.Notebook.register_auto('MacRates_eval.Data.Picture', 'Data picture, on MacRates: Get LL', False) 
160          qubx.notebook.Notebook.register_auto('MacRates_eval.Model.Classes', 'Classes table, on MacRates: Get LL', False) 
161          qubx.notebook.Notebook.register_auto('MacRates_eval.Model.Rates', 'Rates table, on MacRates: Get LL', False) 
162           
163          qubx.notebook.Notebook.register_auto('MacRates_opt.Iterations', 'Iteration log, on MacRates: Optimize', True) 
164          qubx.notebook.Notebook.register_auto('MacRates_opt.Messages', 'Messages, on MacRates: Optimize', True) 
165          qubx.notebook.Notebook.register_auto('MacRates_opt.Model.Picture', 'Model picture, on MacRates: Optimize', True) 
166          qubx.notebook.Notebook.register_auto('MacRates_opt.Model.PictureNO', 'Model picture, no overlays, on MacRates: Optimize', False) 
167          qubx.notebook.Notebook.register_auto('MacRates_opt.Data.Picture', 'Data picture, on MacRates: Optimize', True) 
168          qubx.notebook.Notebook.register_auto('MacRates_opt.Model.Classes', 'Classes table, on MacRates: Optimize', True) 
169          qubx.notebook.Notebook.register_auto('MacRates_opt.Model.Rates', 'Rates table, on MacRates: Optimize', False) 
170          qubx.notebook.Notebook.register_auto('MacRates_opt.Model.Constraints', 'Kinetic Constraints table, on MacRates: Optimize', True) 
171   
172          col = pack_item(gtk.VBox(), self.top_panel) 
173          self.chkPeq = pack_check('start at equilibrium p', col) 
174          self.chkPeq.set_tooltip_text('otherwise use States:Pr') 
175          self.propertied_connect_check('use_Peq', self.chkPeq) 
176          h = pack_item(gtk.HBox(), col) 
177          self.chkResample = pack_check('resample by delta:', h) 
178          self.chkResample.set_tooltip_text('skip slow changes:\ncompress data with std.dev < delta into a single point') 
179          self.propertied_connect_check('resample', self.chkResample) 
180          self.txtResample = pack_item(NumEntry(self.resample_delta, acceptFloatGreaterThan(0.0), '%.3g', width_chars=6), h) 
181          self.propertied_connect_NumEntry('resample_delta', self.txtResample) 
182           
183          pack_space(self.top_panel, x=20) 
184          col = pack_item(gtk.VBox(), self.top_panel) 
185          h = pack_item(gtk.HBox(), col) 
186          pack_label('Gleak:', h) 
187          self.txtGleak = pack_item(NumEntry(self.Gleak, acceptFloat, '%.3g', width_chars=6), h, at_end=True) 
188          self.propertied_connect_NumEntry('Gleak', self.txtGleak) 
189          h = pack_item(gtk.HBox(), col) 
190          pack_label('Vleak:', h) 
191          self.txtVleak = pack_item(NumEntry(self.Vleak, acceptFloat, '%.3g', width_chars=6), h, at_end=True) 
192          self.propertied_connect_NumEntry('Vleak', self.txtVleak) 
193   
194          pack_space(self.top_panel, expand=True) 
195          col = pack_item(gtk.VBox(True), self.top_panel) 
196          lbl_d = pack_label('Multi-data:', col) 
197          lbl_m = pack_label('Multi-model:', col) 
198          col = pack_item(gtk.VBox(True), self.top_panel) 
199          h = pack_item(gtk.HBox(), col) 
200          self.chkUseGroup = pack_check('Group', h) 
201          self.propertied_connect_check('use_group', self.chkUseGroup) 
202          self.txtUseGroup = pack_item(qubx.GTK.NumEntry(self.use_group_index, acceptIntGreaterThan(-1), width_chars=3), h) 
203          self.txtUseGroup.set_sensitive(self.use_group) 
204          self.propertied_connect_NumEntry('use_group_index', self.txtUseGroup) 
205          for widget in (lbl_d, self.chkUseGroup, self.txtUseGroup): 
206              widget.set_tooltip_text('To select data files, set their Group number in the Data table.\nThey will be fit to common parameters') 
207          h = pack_item(gtk.HBox(), col) 
208          self.chkModelsGroup = pack_check('Group', h) 
209          self.propertied_connect_check('multi_model', self.chkModelsGroup) 
210          self.txtModelsGroup = pack_item(qubx.GTK.NumEntry(self.multi_model_group, acceptIntGreaterThan(-1), width_chars=3), h) 
211          self.txtModelsGroup.set_sensitive(self.multi_model) 
212          self.propertied_connect_NumEntry('multi_model_group', self.txtModelsGroup) 
213          for widget in (lbl_m, self.chkModelsGroup, self.txtModelsGroup): 
214              widget.set_tooltip_text("To select models, set their Group number in the Models table.\nData will be fit to the sum of all models' currents.") 
215           
216          pack_space(self.top_panel, expand=True) 
217          col = pack_item(gtk.VBox(True), self.top_panel) 
218          self.chkLive = pack_check('live   ', col) 
219          self.chkLive.set_tooltip_text('update LL and fit curve when the model changes') 
220          self.propertied_connect_check('live', self.chkLive) 
221          self.chkGPU = pack_check('use GPU   ', col) 
222          self.chkGPU.set_tooltip_text('use graphics card if available; faster but less accurate') 
223          self.propertied_connect_check('use_GPU', self.chkGPU) 
224           
225          self.lblScore.set_text('LL:') 
226          self.btnGetScore.set_label('Get LL') 
227   
228          pack_label('Optimize:', self.right_panel) 
229          self.chkOptNchannel = pack_check('channel count', self.right_panel) 
230          self.propertied_connect_check('opt_Nchannel', self.chkOptNchannel) 
231          self.chkOptCond = pack_check('conductance', self.right_panel) 
232          self.propertied_connect_check('opt_cond', self.chkOptCond) 
233          self.chkOptVar = pack_check('variance', self.right_panel) 
234          self.propertied_connect_check('opt_var', self.chkOptVar) 
235          self.chkOptRates = pack_check('rates', self.right_panel, active=True) 
236          self.chkOptRates.set_tooltip_text('To hold a rate constant, right-click it, "Add Constraint...," "Fix..."') 
237          self.chkOptRates.set_sensitive(False) 
238           
239          self.outStatus = TextViewAppender(self.txtStatus) 
240          self.outStatus.write('Based on "Maximum liklihood estimation of ion channel kinetics from macroscopic currents," Milescu LS, Akk G, Sachs F, Biophys J. 2005 Apr; 88(4):2494-515.\n') 
241          self.outStatus.OnWrite += self.__ref(self.__onWriteStatus) 
242          self.__nb_iter_buf = None 
243   
244          self.modelgroup = ModelGroup(self.multi_model_group) 
245          self.modelgroup.OnChange += self.__ref(self.__onChangeModelGroup) 
246          self.modelgroup.OnChangeRates += self.__ref(self.__onChangeModelGroupRates) 
247 -    def mac_finalize(self): 
248          """Stops the worker threads; call this in qubx_plugin_fini.""" 
249          self.robot.stop() 
250          self.fitbot.stop() 
251          self.data = self.data_prep = None 
252 -    def on_exception(self, robot, typ, val, trace): 
253          traceback.print_exception(typ, val, trace) 
254          traceback.print_exception(typ, val, trace, file=self.outStatus) 
255 -    def __onWriteStatus(self, msg): 
256          if not (self.__nb_iter_buf is None): 
257              self.__nb_iter_buf.write(msg) 
258 -    def propertied_set(self, value, name): 
259          if (name, value) == ('use_group_index', None): 
260              self.use_group = False 
261              return 
262          if (name, value) == ('multi_model_group', None): 
263              self.multi_model = False 
264              return 
265          super(RatesFace, self).propertied_set(value, name) 
266          if name in ('use_Peq', 'Gleak', 'Vleak'): 
267              self.invalidate_controls() 
268          if name in ('resample', 'resample_delta'): 
269              self.data.invalidate() 
270          if name in ('use_group', 'use_group_index'): 
271              self.data.use_group = self.use_group_index if self.use_group else None 
272              self.txtUseGroup.set_sensitive(self.use_group) 
273          if (name == 'use_group_index') and not self.use_group: 
274              self.use_group = True 
275          if name == 'multi_model': 
276              self.txtModelsGroup.set_sensitive(self.multi_model) 
277              self.model_prep.invalidate() 
278              self.rates.invalidate() 
279          if (name == 'multi_model_group'): 
280              if not self.multi_model: 
281                  self.multi_model = True 
282              self.modelgroup.group = value 
283          if name == 'live': 
284              if value: 
285                  self.invalidate_controls() 
286          if name == 'use_GPU': 
287              self.invalidate_controls() 
288 -    def __onInvalidateRates(self, prod): 
289          self.invalidate_controls() 
290 -    def invalidate_controls(self, gray=True): 
291          """Changes LL text color gray (outdated) or black (up-to-date), and starts computation if showing in live mode.""" 
292          if gray != self.__invalid_score: 
293              self.__invalid_score = gray 
294              color = gdk.color_parse(gray and "#999999" or "#000000") 
295              self.txtScore.modify_text(gtk.STATE_NORMAL, color) 
296          if gray and self.chkLive.get_active(): 
297              self.__queue_robot_live() 
298 -    def __queue_robot_live(self): 
299          self.__invalid_live = True 
300          if not (self.robot.working or self.optimizing): 
301              self.__capture_robot_settings() 
302              self.robot.do(self.robot_live) 
303 -    def __sensitize(self, sens): 
304          try: 
305              self.chkPeq.set_sensitive(sens) 
306          except: 
307              pass 
308 -    def on_robot_begin(self): 
309          super(RatesFace, self).on_robot_begin() 
310          gobject.idle_add(self.__sensitize, False) 
311 -    def on_robot_end(self): 
312          super(RatesFace, self).on_robot_end() 
313          gobject.idle_add(self.__sensitize, True) 
314 -    def fitbot_send_datapic(self, auto_id): 
315          gobject.idle_add(lambda: qubx.notebook.Notebook.send(self.QubX.Data.view.hires.nbPicture, auto_id=auto_id)) 
316 -    def get_ndata(self, data_prep=None): 
317          if data_prep and data_prep.workspaces: 
318              return reduce(operator.add, (ws.data.Ndata for ws in data_prep.workspaces), 0) 
319          else: 
320              return self.__Ndata 
321 -    def __onChangeModelGroup(self): 
322          if self.multi_model: 
323              self.model_prep.invalidate() 
324              self.rates.invalidate() 
325 -    def __onChangeModelGroupRates(self): 
326          self.rates.invalidate() 
327 -    def __get_rates(self): 
328          model_prep = self.model_prep.val 
329          signal_names = iunion((prep.source or prep.model).variables for prep in model_prep.preps) 
330          return [ModelToRateMatricesP(prep.source or prep.model, signal_names)[:3] for prep in model_prep.preps] 
331 -    def add_trial_fields(self, row, mtree, prefix='', model_prep=None, data_prep=None): 
332          add_trial_rates(row, mtree, prefix) 
333          if (not prefix) or self.robot.optNchannel: 
334              add_trial_nChannel(row, mtree, prefix) 
335          if (not prefix) or self.robot.optCond or self.robot.optVar: 
336              add_trial_amps(row, mtree, prefix) 
337 -    def __onGotLL(self, LL): 
338          self.LL = LL 
339          if not self.optimizing: 
340              qubx.notebook.Notebook.send(qubx.notebook.NbText("""MacRates: Get LL: %.6g""" % LL), auto_id='MacRates_eval.Messages') 
341              if self.QubX.Models.file in self.robot.model_prep.sources: ### TODO: send all sources 
342                  qubx.notebook.Notebook.send(self.QubX.Models.view.nbPicture, auto_id='MacRates_eval.Model.Picture') 
343                  qubx.notebook.Notebook.send(self.QubX.Models.view.nbPictureNO, auto_id='MacRates_eval.Model.PictureNO') 
344                  qubx.notebook.Notebook.send(self.QubX.Tables.find_view('Classes').nbTable, auto_id='MacRates_eval.Model.Classes') 
345                  qubx.notebook.Notebook.send(self.QubX.Tables.find_view('Rates').nbTable, auto_id='MacRates_eval.Model.Rates') 
346              # fit robot send datapic with updated fit curves 
347              self.fitbot.do(self.fitbot_send_datapic, 'MacRates_eval.Data.Picture') 
348              if self.__invalid_live: 
349                  self.__queue_robot_live() 
350              else: 
351                  self.invalidate_controls(False) 
352              self.set_score(LL) 
353              ##self.last_LL = LL 
354 -    def calculate_score(self, model_prep=None, data_prep=None, wait=True, receiver=lambda score: None, on_trial=None, **kw): 
355          self.__capture_robot_settings() 
356          return super(RatesFace, self).calculate_score(model_prep, data_prep, wait, receiver, on_trial, **kw) 
357 -    def optimize(self, model_prep=None, data_prep=None, wait=True, receiver=None, ext_return=False, on_trial=None): 
358          self.__nb_iter_buf = cStringIO.StringIO() 
359          self.__capture_robot_settings() 
360          self.QubX.Models.view.optimizing = True 
361          return super(RatesFace, self).optimize(model_prep, data_prep, wait, receiver, ext_return, on_trial) 
362 -    def __onOptDone(self, LL, iterations, grads, Hessian, workspaces, model_prep, data_prep): 
363          self.QubX.Models.view.optimizing = False 
364          self.LL = LL 
365          self.iterations = iterations 
366          self.grads = grads 
367          iteration_messages = self.__nb_iter_buf.getvalue() if self.__nb_iter_buf else "" 
368          self.__nb_iter_buf = cStringIO.StringIO() # to capture final output for notebook 
369          for source in model_prep.sources: 
370              source.undoStack.seal_undo('set rates') 
371          Ndata = self.__Ndata 
372          self.append_output('===================================================\n') 
373          self.outStatus.write("""MacRates: finished optimizing %d points after %d iterations. 
374  Data:           %s 
375  Model:          %s 
376  P0:             %s 
377  Resample:       %s 
378  Gleak:          %.4g 
379  Vleak:          %.4g 
380  GPU:            %s 
381  Optimize:       %s 
382  Gradient:       %s 
383  """ % (Ndata, 
384                        iterations, 
385                        '\n           '.join(os.path.split(ws.data.file.path)[1] for ws in workspaces), 
386                        model_prep.sources and ', '.join([os.path.split(source.path)[1] for source in model_prep.sources]) or "from script", 
387                        self.robot.Peq and "equilibrium" or "fixed (States table)", 
388                        self.data.resample and ("by delta = %.4g" % self.data.resample_std) or "no", 
389                        self.robot.Gleak, 
390                        self.robot.Vleak, 
391                        self.chkGPU.get_active() and "yes" or "no", 
392                        ((self.robot.optNchannel and "channel count, " or "") + 
393                         (self.robot.optCond and "conductance, " or "") + 
394                         (self.robot.optVar and "variance, " or "") + 
395                         "rates"), 
396                        ', '.join('%8.5g' % x for x in grads))) 
397          if workspaces: 
398              self.__set_model_confidence(workspaces[0], model_prep) 
399          self.outStatus.write('\nLL: ') 
400          self.outStatus.write_bold('%.5g' % LL) 
401          if self.last_LL: 
402              self.outStatus.write('        last LL: %.5g' % self.last_LL) 
403          self.append_output('\n===================================================\n') 
404          self.invalidate_controls(False) 
405          qubx.notebook.Notebook.send(qubx.notebook.NbText(iteration_messages), auto_id='MacRates_opt.Iterations') 
406          qubx.notebook.Notebook.send(qubx.notebook.NbText(self.__nb_iter_buf.getvalue()), auto_id='MacRates_opt.Messages') 
407          self.__nb_iter_buf = None 
408          if self.QubX.Models.file in model_prep.sources: ### TODO: send all source models 
409              qubx.notebook.Notebook.send(self.QubX.Models.view.nbPicture, auto_id='MacRates_opt.Model.Picture') 
410              qubx.notebook.Notebook.send(self.QubX.Models.view.nbPictureNO, auto_id='MacRates_opt.Model.PictureNO') 
411              qubx.notebook.Notebook.send(self.QubX.Tables.find_view('Classes').nbTable, auto_id='MacRates_opt.Model.Classes') 
412              qubx.notebook.Notebook.send(self.QubX.Tables.find_view('Rates').nbTable, auto_id='MacRates_opt.Model.Rates') 
413              qubx.notebook.Notebook.send(self.QubX.Tables.find_view('Kinetic Constraints').nbTable, auto_id='MacRates_opt.Model.Constraints') 
414          # fit robot send datapic with updated fit curves 
415          self.fitbot.do(self.fitbot_send_datapic, 'MacRates_opt.Data.Picture') 
416          self.last_opt_model = workspaces[0].param.modelinfo[0].fastmodel if workspaces else None 
417          self.last_LL = LL 
418          self.OnOptimized(iterations, LL, grads, Hessian) 
419 -    def __set_model_confidence(self, ws, model_prep): 
420          for mi in ws.param.modelinfo: 
421              self.outStatus.write("\n", scroll=False) # model name? 
422              if self.robot.optNchannel: 
423                  self.outStatus.write("Channel count:  ", scroll=False) 
424                  self.outStatus.write_bold('%i' % mi.fastmodel.Nchannel, scroll=False) 
425                  self.outStatus.write(" +/- %.1f\n" % mi.fastmodel.Nchannel_err, scroll=False) 
426              if self.robot.optCond or self.robot.optVar: 
427                  amp_lbl, std_lbl = ('Amp', 'Variance') if (not mi.fastmodel.IeqFv) else ('Cond', 'CondVar') 
428                  if mi.fastmodel.IeqFv and mi.fastmodel.iVoltage: 
429                      if self.robot.optCond: 
430                          self.outStatus.write('Baseline Amp:  ', scroll=False) 
431                          self.outStatus.write_bold('%.5g' % mi.fastmodel.amp0, scroll=False) 
432                          self.outStatus.write(' +/- %.5g\n' % mi.fastmodel.amp0_err, scroll=False) 
433                      if self.robot.optVar: 
434                          self.outStatus.write('Baseline Variance: ', scroll=False) 
435                          self.outStatus.write_bold('%.5g' % mi.fastmodel.var0, scroll=False) 
436                          self.outStatus.write(' +/- %.5g\n' % mi.fastmodel.var0_err, scroll=False) 
437                  self.outStatus.write('%7s' % 'State', scroll=False) 
438                  if self.robot.optCond: 
439                      self.outStatus.write_bold('%14s' % amp_lbl, scroll=False) 
440                      self.outStatus.write('%14s' % (amp_lbl+' err'), scroll=False) 
441                  if self.robot.optVar: 
442                      self.outStatus.write_bold('%14s' % std_lbl, scroll=False) 
443                      self.outStatus.write('%14s' % (std_lbl+' err'), scroll=False) 
444                  self.outStatus.write('\n', scroll=False) 
445                  for c in xrange(mi.fastmodel.Nclass): 
446                      self.outStatus.write('%7d' % c) 
447                      if self.robot.optCond: 
448                          self.outStatus.write_bold('%14.5g' % mi.fastmodel.amp[c], scroll=False) 
449                          self.outStatus.write('%14.5g' % mi.fastmodel.amp_err[c], scroll=False) 
450                      if self.robot.optVar: 
451                          self.outStatus.write_bold('%14.5g' % mi.fastmodel.var[c], scroll=False) 
452                          self.outStatus.write('%14.5g' % mi.fastmodel.var[c], scroll=False) 
453                      self.outStatus.write('\n', scroll=False) 
454              self.outStatus.write('%10s%7s' % ("Rate From", "To"), scroll=False) 
455              self.outStatus.write_bold('%14s' % "k0", scroll=False) 
456              self.outStatus.write('%11s' % "+/-", scroll=False) 
457              self.outStatus.write_bold('%14s' % 'k1', scroll=False) 
458              self.outStatus.write('%11s' % "+/-", scroll=False) 
459              self.outStatus.write_bold('%14s' % 'k2', scroll=False) 
460              self.outStatus.write('%11s' % "+/-", scroll=False) 
461          for mi, source in itertools.izip(ws.param.modelinfo, model_prep.sources): 
462              self.outStatus.write('\n', scroll=False) 
463              for rate in source.rates: 
464                  a,b = rate.From, rate.To 
465                  source.rates[rate.Index, 'k0 std'] = mi.fastmodel.K0_err[a,b] 
466                  source.rates[rate.Index, 'k1 std'] = mi.fastmodel.K1_err[a,b] 
467                  source.rates[rate.Index, 'k2 std'] = mi.fastmodel.K2_err[a,b] 
468                  self.outStatus.write('%10d%7d' % (a, b), scroll=False) 
469                  self.outStatus.write_bold('%14.5g' % rate.k0, scroll=False) 
470                  self.outStatus.write('%11.3g' % mi.fastmodel.K0_err[a,b], scroll=False) 
471                  self.outStatus.write_bold('%14.5g' % rate.k1, scroll=False) 
472                  self.outStatus.write('%11.3g' % mi.fastmodel.K1_err[a,b], scroll=False) 
473                  self.outStatus.write_bold('%14.5g' % rate.k2, scroll=False) 
474                  self.outStatus.write('%11.3g' % mi.fastmodel.K2_err[a,b], scroll=False) 
475                  self.outStatus.write('\n') 
476 -    def __onIteration(self, iteration, LL, grads, optNchannel, workspaces, modelinfo, K0s, K1s, K2s, amps, stds, amp0s, std0s, Nchannels): 
477          self.LL = LL 
478          self.iterations = iteration 
479          self.grads = grads 
480          for mi, K0, K1, K2, amp, std, amp0, std0, Nchannel in itertools.izip(modelinfo, K0s, K1s, K2s, amps, stds, amp0s, std0s, Nchannels): 
481              model = mi.prep.source 
482              if not model: 
483                  continue 
484              K012ToModel(K0, K1, K2, model) 
485              if model.IeqFv: 
486                  for c,a in enumerate(amp): 
487                      model.classes.set(c, 'Cond', a) 
488                      model.classes.set(c, 'CondStd', std[c]) 
489                  model.classes.set(0, 'Amp', amp0) 
490                  model.classes.set(0, 'Std', std0) 
491              else: 
492                  for c,a in enumerate(amp): 
493                      model.classes.set(c, 'Amp', a) 
494                      model.classes.set(c, 'Std', std[c]) 
495              model.set_channelCount(Nchannel, seal=False) # seal: don't want undo-per-iteration 
496              if (self.__Ndata > 1) and optNchannel and (len(modelinfo) == 1): 
497                  for i, ws in enumerate(workspaces): 
498                      datafile = qubx.global_namespace.QubX.Data.views[ws.param.data.file_index].file 
499                      try: 
500                          cc_ix = datafile.constants.index('Channel Count') 
501                          datafile.constants[cc_ix, 'Value'] = ws.param.modelinfo[0].fastmodel.Nchannel 
502                      except KeyError: 
503                          datafile.constants.append({'Name' : 'Channel Count', 'Value' : ws.param.modelinfo[0].fastmodel.Nchannel}) 
504          self.invalidate_controls(False) 
505          self.outStatus.write('Iter: %i\tLL: %8.5g\tgrads: %s\n' % (iteration, LL, ', '.join("%8.5g"%g for g in grads))) 
506          self.OnIteration(iteration, LL, grads, K0, K1, amp, std, amp0, std0, Nchannel) 
507          self.OnIteration2(iteration, LL, grads, K0, K1, K2, amp, std, amp0, std0, Nchannel) 
508 -    def on_interrupt(self, task, cancel): 
509          cancel() # the flag is good enough, while the KeyboardInterrupt can bewilder it 
510          for ws in self.robot.workspaces: 
511              for mi in ws.param.modelinfo: 
512                  mi.fastmodel.obj[0].stop_flag = 1 
513 -    def get_model_prep(self, model=None): 
514          if (model is None) and self.multi_model: 
515              # make it look like a single-model prep for compat, but add .preps 
516              mp = Anon(preps=[], models=[], sources=[]) 
517              for i, mv in enumerate(qubx.global_namespace.QubX.Models.views): 
518                  if qubx.global_namespace.QubX.Models.table[i, 'Group'] == self.multi_model_group: 
519                      mp1 = qubx_plugin.ModelPrep(mv.file) 
520                      mp.preps.append(mp1) 
521                      mp.models.append(mp1.model) 
522                      mp.sources.append(mp1.source) 
523              if mp.preps: 
524                  mp.model = mp.models[0] 
525                  mp.source = mp.sources[0] 
526          else: 
527              mp = super(RatesFace, self).get_model_prep(model) 
528              # either way, list all models and sources for loop processing 
529              mp.models = [mp.model] 
530              mp.sources = [mp.source] if mp.source else [] 
531              mp.preps = [mp] # please please please garbage collect despite this cycle 
532          mp.set_rates2 = set_modelprep_rates(mp) 
533          return mp 
534 -    def request_data_prep(self, receiver, segs, model_prep, last_data_prep=None, resample=None, resample_std=None): 
535          QubX = qubx.global_namespace.QubX 
536          view_of_file = dict([(view.file, view) for view in QubX.Data.views]) 
537          view_segm = [] 
538          vs_by_file = {} 
539          for seg in segs: 
540              if not (seg.file in vs_by_file): 
541                  vs = (view_of_file[seg.file], []) 
542                  view_segm.append(vs) 
543                  vs_by_file[seg.file] = vs 
544              vs_by_file[seg.file][1].append(seg) 
545          if not view_segm: 
546              if receiver: 
547                  gobject.idle_add(receiver, DataPrep(segs, [])) 
548              return 
549   
550          resamp = QubX.Modeling.MacRates.resample if (resample is None) else resample 
551          resamp_std = QubX.Modeling.MacRates.resample_delta if (resample_std is None) else resample_std 
552   
553          datas = [] 
554          data_prep = DataPrep(segs, datas) 
555          self.last_LL = 0.0 
556          # segment is divided in chunks, alternating included/excluded (data / no data) 
557          # each segment, we request idl of each signal, each chunk 
558          # this could involve re-idealization, indefinite time, so result comes via callback 
559          # we provide generator instance as callback, repeatedly, to process segments in a loop 
560           
561          # gobject thread:                               stimulus robot: 
562          #   rcv_stims.next() 
563          #      repeat: 
564          #         stimulus.call_with_stimulus_idl 
565          #         yield... 
566          #   [get(self, ...) returns] 
567          #                                                  idealize stim if needed 
568          #                                                  idle_add(rcv_stims.send, ...) 
569          #         ...process chunks 
570          #      yield None 
571           
572          rcv_stims = None # will have instance of receive_stims() 
573          def rcv_hack(*args): # can send only one argument into a generator 
574              rcv_stims.send(args) 
575          def receive_stims(): 
576              stim_names = iunion((prep.source or prep.model).variables for prep in model_prep.preps) 
577              for dataview, segsrc in view_segm: 
578                  # re-use L{max_mac_ll.Data} object: (single-file only) 
579                  if last_data_prep and len(last_data_prep.datas) and (len(view_segm) == 1): 
580                      data = last_data_prep.datas[0] 
581                      data.reset(len(stim_names)) 
582                  else: 
583                      data = max_mac_ll.Data(len(stim_names), file_index=QubX.Data.views.index(dataview)) # , data_cap_samples) 
584               
585                  sampling_sec = dataview.file.sampling 
586                  if 'Voltage' in stim_names: 
587                      data.iVoltage = stim_names.index('Voltage') + 1 
588                  else: 
589                      data.iVoltage = 0 
590                  data.signal_names = ['Current'] + stim_names 
591                  data.ff_ll_counts = [] 
592   
593                  iseg = -1 
594                  for seg in segsrc: 
595                      data.file = seg.file 
596                      data.iSignal = seg.signal 
597                      QubX.Modeling.Stimulus.call_with_stimulus_idl(rcv_hack, [(chunk.f, chunk.l) for chunk in seg.chunks], datafile=dataview.file, names=stim_names) 
598                      stim_names, stimidl = (yield) 
599                      chunkidl = [ Anon(dwellCounts=[], classeses=[], durationses=[], ampses=[]) for chunk in seg.chunks] 
600                      missing = False 
601                      for segs_flcda in stimidl: # for each signal 
602                          for ichunk, flcda in enumerate(segs_flcda): # each chunk 
603                              ff, ll, cc, dd, aa  = flcda 
604                              if not len(cc): 
605                                  missing = True 
606                              chunkidl[ichunk].dwellCounts.append(len(cc)) 
607                              chunkidl[ichunk].classeses.append(cc) 
608                              chunkidl[ichunk].durationses.append(dd) 
609                              chunkidl[ichunk].ampses.append(aa) 
610                      if missing: 
611                          continue 
612                   
613                      iseg += 1 
614                      data.append_segment() 
615                      for ichunk, chunk in enumerate(seg.chunks): 
616                          if chunk.included: 
617                              Nd = chunk.n 
618                              chunkoff = chunk.f - seg.f 
619                              xx = arange(chunkoff, chunkoff+Nd, dtype='float32') 
620                              xx *= sampling_sec 
621                              yy = zeros(shape=(Nd,), dtype='float32') 
622                              ii = 0 
623                              for ch in generate_chunk_samples([chunk]): 
624                                  yy[ii:ii+ch.n] = ch.samples 
625                                  ii += ch.n 
626                              if resamp: 
627                                  means, stds, ff, ll, closest = qubx.fast.data.adaptive_resample(yy, resamp_std) 
628                                  xx = array(xx[closest], dtype='float32', copy=True) 
629                                  yy = array(yy[closest], dtype='float32', copy=True) 
630                                  ff += chunk.f 
631                                  ll += chunk.f 
632                                  Nd = len(means) 
633                              else: 
634                                  ff = arange(seg.f+chunkoff, seg.f+chunkoff+Nd, dtype='int32') 
635                                  ll = ff 
636                              data.ff_ll_counts.append( (ff, ll, Nd) ) 
637                          else: 
638                              Nd = 0 
639                              xx = yy = None 
640                          data.append_chunk(iseg, (chunk.l - chunk.f + 1) * sampling_sec, Nd, xx, yy, None, 
641                                            chunkidl[ichunk].dwellCounts, chunkidl[ichunk].classeses, 
642                                            chunkidl[ichunk].durationses, chunkidl[ichunk].ampses) 
643                  datas.append(data) 
644              gobject.idle_add(receiver, data_prep) 
645              yield None 
646   
647          rcv_stims = receive_stims() 
648          rcv_stims.next() 
649 -    def __onGPUFail(self): 
650          self.chkGPU.set_active(False) 
651          self.chkGPU.set_sensitive(False) 
652          self.chkGPU.set_tooltip_text('no OpenCL-compatible graphics processor is available') 
653 -    def wait_for_fits(self): 
654          """Waits until all pending work is done on the thread which updates fit curves, then returns.""" 
655          while self.fitbot.last_serial != self.fitbot.serial: 
656              qubx.pyenv.env.process_events() 
657 -    def __capture_robot_settings(self): 
658          """Should happen at the same time the modelprep is captured, e.g. at button press""" 
659          self.robot.Peq = self.use_Peq 
660          self.robot.GPU = self.chkGPU.get_active() 
661          self.robot.Gleak = self.Gleak 
662          self.robot.Vleak = self.Vleak 
663          self.robot.optNchannel = self.opt_Nchannel 
664          self.robot.optCond = self.opt_cond 
665          self.robot.optVar = self.opt_var 
666          self.data.resample = self.resample 
667          self.data.resample_std = self.resample_delta 
668          self.robot.K012 = self.rates.val 
669          self.robot.rates_source = qubx.global_namespace.QubX.Models.file 
670          if self.robot.GPU: 
671              self.robot.accel = self.accel and (self.accel + int(self.robot.GPU)) or 0 
672          else: 
673              self.robot.accel = self.accel 
674 -    def robot_setup(self): 
675          model_prep = self.robot.model_prep 
676          self.robot.models = model_prep.models 
677          self.robot.modelinfo = [Anon(model=prep.model, prep=prep) for prep in model_prep.preps if prep.model.states.size] 
678          self.data.model_prep = model_prep 
679           
680          for mi in self.robot.modelinfo: 
681              mi.clazz = ModelToClazz(mi.model) 
682              mi.Nclass = int(max(mi.clazz) + 1) if mi.clazz.shape[0] else 0 
683   
684          datas = self.robot.data = self.data_prep.datas # if Nclass else [] 
685          if self.robot.workspaces: 
686              for ws in self.robot.workspaces: 
687                  try: 
688                      ws.fastmodel.dispose() # work around callback-related memory leak 
689                  except: 
690                      pass 
691          workspaces = self.robot.workspaces = [] 
692          self.__Ndata = reduce(operator.add, (data.Ndata for data in datas), 0) 
693   
694          for data in datas: 
695              stimclasses = data.get_stimclasses() 
696              Nstim = data.Nstim 
697              Nsig = data.Nsig 
698              sc_frac = data.get_stimclass_frac() 
699              param = max_mac_ll.Param(len(self.robot.modelinfo)) 
700              param.modelinfo = [] 
701              param.data = data 
702              param.useVar = True ## 
703              param.optNchannel = self.robot.optNchannel 
704              param.accel = self.robot.accel 
705   
706              for im, mi in enumerate(self.robot.modelinfo): 
707                  fastmodel = ModelToFast(mi.model, data.signal_names, self.robot_on_status) 
708                  fastmodel.usePeq = self.robot.Peq 
709                  fastmodel.iVoltage = data.iVoltage 
710                  fastmodel.Gleak = self.robot.Gleak 
711                  fastmodel.Vleak = self.robot.Vleak 
712                  K0, K1, K2, L, V, P = ModelToRateMatricesP(mi.model, data.signal_names) 
713                  fastmodel.K0[:], fastmodel.K1[:], fastmodel.K2[:] = K0, K1, K2 = self.robot.K012[im] 
714   
715                  if len(datas) > 1: 
716                      try: 
717                          fastmodel.Nchannel = int(round(qubx.global_namespace.QubX.Data.views[data.file_index].file.constants.get_by_name('Channel Count', 'Value'))) 
718                      except: 
719                          traceback.print_exc() 
720                          pass # use value from model 
721           
722                  ratesm = qubx.fast.model.StimRates(fastmodel, Nstim, Nsig, stimclasses, sc_frac, 0.0, custom_eigen=None) # skip ratesm.AA generation 
723                  A, B = mi.model.constraints_kin.get_matrices_p(K0, K1, K2, L, V, P, auto=False) 
724                  if A != None: 
725                      ratesm.cns.Ain[:A.shape[0],:A.shape[1]] = A 
726                      ratesm.cns.Bin[:B.shape[0]] = B.flatten() 
727                      ratesm.cns.Ncns = B.shape[0] 
728                  else: 
729                      ratesm.cns.Ncns = 0 
730                       
731                  ampm = qubx.fast.model.StimAmps(fastmodel, Nstim, Nsig, stimclasses, sc_frac) 
732                  if not self.robot.optCond: 
733                      for c in xrange(fastmodel.Nclass): 
734                          ampm.fix_amp(c) 
735                      if fastmodel.IeqFv: 
736                          ampm.fix_amp0() 
737                  if not self.robot.optVar: 
738                      for c in xrange(fastmodel.Nclass): 
739                          ampm.fix_var(c) 
740                      if fastmodel.IeqFv: 
741                          ampm.fix_var0() 
742                  param.add_model(fastmodel, ampm, ratesm) 
743                  param.modelinfo.append(mi) 
744                  mi.fastmodel = fastmodel 
745                  mi.Nchannel = fastmodel.Nchannel 
746   
747              if param.models: 
748                  param.setup_constraints() 
749                  workspaces.append(Anon(param=param, data=data, modelinfo=param.modelinfo, L=L, V=V)) 
750 -    def robot_live(self, model_prep=None): 
751          if self.__invalid_live and not (self.robot.working or self.optimizing): 
752              if model_prep is None: 
753                  with self.robot.main_hold: 
754                      mp = self.model_prep.val 
755              else: 
756                  mp = model_prep 
757              self.robot_setup_data(mp, None) 
758              LL = self.robot_calculate_score() 
759              gobject.idle_add(self.OnLiveLL, LL) 
760 -    def robot_calculate_score(self, grads=None): 
761          LL = 0.0 
762          try: 
763              self.robot.working = True 
764              self.__invalid_live = False 
765              if not self.optimizing: 
766                  self.robot_setup() 
767              LL = 0.0 
768              if ('RejectMacModels' in qubx.global_namespace.__dict__) and qubx.global_namespace.RejectMacModels([ 
769                      mi.fastmodel for mi in self.robot.modelinfo 
770              ]): 
771                  LL = -1e20 # super-unlikely 
772              else: 
773                  for ws in self.robot.workspaces: 
774                      try: 
775                          LL += max_mac_ll.mac_ll(ws.param) 
776                      except: 
777                          gobject.idle_add(self.outStatus.write, traceback.format_exc() + '\n') 
778                          LL += -1.0 
779                      if ws.param.gpu_fail: 
780                          gobject.idle_add(self.__onGPUFail) 
781                          self.robot.accel -= 1 
782                          for ws2 in self.robot.workspaces: 
783                              ws2.param.accel = self.robot.accel 
784              self.robot_queue_fits(0, True) 
785          except: 
786              self.robot.send_exception() 
787          gobject.idle_add(self.__onGotLL, LL) 
788          self.robot.working = False 
789          return LL 
790 -    def robot_get_pars(self): 
791          """Returns the list of initial parameter values, on worker thread; override this method. 
792  Derive parameter values from self.robot.model and/or from your own widgets.""" 
793          self.robot_setup() 
794          if self.robot.workspaces: 
795              pars = [float(x) for x in self.robot.workspaces[0].param.fPars[:self.robot.workspaces[0].param.NfPar]] # get as python array 
796          else: 
797              pars = [] 
798          # fastmodel has room for one channel count, used for first file; append a param for the rest 
799          if (len(self.robot.workspaces) > 1) and self.robot.optNchannel and (len(self.robot.models) == 1): 
800              pars.extend([1.0] * (len(self.robot.workspaces)-1)) 
801              self.__scale_Nchannel = [log(ws.param.Nchannel) for ws in self.robot.workspaces[1:]] 
802          return array(pars) 
803 -    def robot_set_pars(self, pars): 
804          if (pars != None) and self.robot.workspaces: 
805              NfPar = self.robot.workspaces[0].param.NfPar 
806              self.robot.workspaces[0].param.fPars[:NfPar] = pars[:NfPar] 
807              self.robot.workspaces[0].param.do_fPars_to_model() 
808              if (len(self.robot.workspaces) > 1) and self.robot.optNchannel and (len(self.robot.models) == 1): 
809                  ccpar = pars[-(len(self.robot.workspaces)-1):] 
810                  for i, ws in enumerate(self.robot.workspaces[1:]): 
811                      float_Nchannel = max(1.0, exp(self.__scale_Nchannel[i] * ccpar[i])) 
812                      ws.param.float_Nchannel = float_Nchannel 
813                      ws.param.modelinfo[0].fastmodel.Nchannel = int(round(float_Nchannel)) 
814 -    def robot_optimize(self): 
815          try: 
816              self.robot.working = True 
817              self.robot.LL = 0.0 
818              self.__last_fits_iter = -1 
819              iterations = 0 
820              pars = grads = [] 
821              Hessian = None 
822              self.robot.LL, pars, iterations, grads, Hessian = super(RatesFace, self).robot_optimize() 
823              self.robot.workspaces[0].param.calc_std(Hessian) 
824              if self.accel and (self.__last_fits_iter < iterations): 
825                  self.robot_queue_fits(iterations, True) 
826          except: 
827              self.robot.send_exception() 
828          gobject.idle_add(self.__onOptDone, self.robot.LL, iterations, grads, Hessian, self.robot.workspaces, self.robot.model_prep, self.data_prep) 
829          self.robot.working = False 
830          return self.robot.LL, pars, iterations, grads, Hessian 
831 -    def robot_on_iteration(self, pars, grads, iterations, score): 
832          if self.robot.workspaces[0].param.modelinfo[0].fastmodel.obj[0].stop_flag: 
833              raise KeyboardInterrupt() 
834          else: 
835              self.robot_set_pars(pars) 
836              self.robot_queue_fits(iterations) 
837              self.robot.LL = score 
838              gobject.idle_add(self.__onIteration, iterations, - score, grads, self.robot.optNchannel, 
839                               self.robot.workspaces, self.robot.workspaces[0].param.modelinfo, 
840                               [array(mi.fastmodel.K0, copy=True) for mi in self.robot.workspaces[0].param.modelinfo], 
841                               [array(mi.fastmodel.K1, copy=True) for mi in self.robot.workspaces[0].param.modelinfo], 
842                               [array(mi.fastmodel.K2, copy=True) for mi in self.robot.workspaces[0].param.modelinfo], 
843                               [array(mi.fastmodel.amp, copy=True) for mi in self.robot.workspaces[0].param.modelinfo], 
844                               [numpy.sqrt(array(mi.fastmodel.var, copy=True)) for mi in self.robot.workspaces[0].param.modelinfo], 
845                               [mi.fastmodel.amp0 for mi in self.robot.workspaces[0].param.modelinfo], 
846                               [sqrt(mi.fastmodel.var0) for mi in self.robot.workspaces[0].param.modelinfo], 
847                               [mi.fastmodel.Nchannel for mi in self.robot.workspaces[0].param.modelinfo]) 
848 -    def robot_queue_fits(self, iteration=0, force=False): 
849          if (not force) and self.optimizing and self.accel: 
850              if (datetime.datetime.now() - self.__last_fits_time).seconds < MIN_SECONDS_BETWEEN_FIT: 
851                  return 
852          self.__last_fits_time = datetime.datetime.now() 
853          self.__last_fits_iter = iteration 
854          self.fitbot.serial += 1 # shared serial so it updates all files or none 
855          for data in self.robot.data: 
856              #print 'Ndata:',self.robot.data.Ndata 
857              chunks = [] 
858              samples = data.get_buffer() 
859              Ndata = data.Ndata 
860              start = 0 
861              for firsts, lasts, count in data.ff_ll_counts: 
862                  means = numpy.array(samples[start:start+count, 0], copy=True) 
863                  stds = numpy.sqrt(array(samples[start:start+count, 1])) # what if some are neg/nan? 
864                  #means = array([sample.fit for sample in samples[start:start+count]], dtype='float32') 
865                  #stds = array([sqrt(sample.fit_var) for sample in samples[start:start+count]], dtype='float32') 
866                  chunks.append((firsts, lasts, means, stds)) 
867                  start += count 
868                  means[numpy.isinf(means)] = 0.0 
869                  means[numpy.isnan(means)] = 0.0 
870                  stds[numpy.isinf(stds)] = 0.0 
871                  stds[numpy.isnan(stds)] = 0.0 
872              self.fitbot.do(self.fitbot_set_fits, data.file, data.iSignal, chunks, self.fitbot.serial) 
873              #self.robot_set_fits(self.robot.data.file, self.robot.data.iSignal, chunks, self.fitbot.serial) 
874 -    def robot_on_status(self, s): 
875          gobject.idle_add(self.outStatus.write, s+"\n") 
876 -    def fitbot_set_fits(self, datafile, iSignal, chunks, serial): 
877          # separate thread to wait for gui to finish drawing, discard stale iterations' fit curves, apply latest 
878          if serial == self.fitbot.serial: 
879              gobject.idle_add(self.do_set_fits, datafile, iSignal, chunks) 
880          self.fitbot.last_serial = serial 
881 -    def do_set_fits(self, datafile, iSignal, chunks): 
882          fits = datafile.fits[iSignal].idl 
883          fits.clear() 
884          for ff, ll, mm, ss in chunks: 
885              fits.set_fit(ff, ll, mm, ss) 
886          datafile.OnChangeFits(datafile, iSignal) 
887   
888   
889 -def set_modelprep_rates(mp): 
890      # override of ModelPrep.set_rates2 for multi_model 
891      def set_rates2(*K012s): 
892          for prep, K012 in itertools.izip(mp.preps, K012s): 
893              prep.K0, prep.K1, prep.K2 = K012 
894              K012ToModel(prep.K0, prep.K1, prep.K2, prep.model) 
895      return set_rates2 
896   
897 -class ModelGroup(object): 
898 -    def __init__(self, group=0): 
899          self.__group = 0 
900          self.__ref = Reffer() 
901          self.OnChange = WeakEvent() 
902          self.OnChangeRates = WeakEvent() 
903          self.Models = qubx.global_namespace.QubX.Models 
904          self.Models.table.OnSet += self.__ref(self.__onSetTable) 
905          self.Models.table.OnAfterInsert += self.__ref(self.__onAfterInsertTable) 
906          self.Models.table.OnRemoving += self.__ref(self.__onRemovingTable) 
907          for view in self.Models.views: 
908              self.__subscribe(view.file) 
909 -    def set_group(self, x): 
910          if self.__group == x: 
911              return 
912          self.__group = x 
913          self.OnChange() 
914      group = property(lambda self: self.__group, lambda self, x: self.set_group(x)) 
915 -    def __subscribe(self, model, adding=True): 
916          events = (model.OnSetChannelCount, 
917                    model.OnSetIeqFv, 
918                    model.OnSetVRev, 
919                    model.states.OnSet, 
920                    model.states.OnInsert, 
921                    model.states.OnRemoved, 
922                    model.classes.OnSet, 
923                    model.rates.OnInsert, 
924                    model.rates.OnRemoved, 
925                    model.OnChangeBalanceLoops, 
926                    model.constraints_kin.OnEdit) 
927          if adding: 
928              for event in events: 
929                  event += self.__ref(self.__onChangeModel) 
930              model.rates.OnSet += self.__ref(self.__onSetRate) 
931          else: 
932              for event in events: 
933                  event -= self.__ref(self.__onChangeModel) 
934              model.rates.OnSet -= self.__ref(self.__onSetRate) 
935          self.OnChange() 
936 -    def __onChangeModel(self, *args): 
937          self.OnChange() 
938 -    def __onSetRate(self, r, field, val, prev, undoing): 
939          self.OnChangeRates() 
940 -    def __onSetTable(self, index, field_name, val, prev, undoing): 
941          if field_name == 'Group': 
942              self.OnChange() 
943 -    def __onAfterInsertTable(self, index, undoing): 
944          self.__subscribe(self.Models.views[index].file) 
945          if self.Models.table[index, 'Group'] == self.__group: 
946              self.OnChange() 
947 -    def __onRemovingTable(self, index, undoing): 
948          self.__subscribe(self.Models.views[index].file, False) 
949          if self.Models.table[index, 'Group'] == self.__group: 
950              self.OnChange() 
951   
952   
953   
954  ## [x] start at equilibrium P           Gleak: [_________0.0_]  [x] live  LL: [___-33952.5_] ## 
955  ## [x] resample by delta [___1.0_]      Vleak: [_________0.0_]            [[    Get LL    ]] ## 
956  ## |==========================================================| Optimize:                    ## 
957  ## |MacRates started...                                       | [x] channel count            ##  
958  ## | 1:      -33129                                           | [x] amp                      ## 
959  ## | 2:      -32981                                           | [x] var                      ## 
960  ## | 3:      -32360                                           | [x] rates (grayed out; tooltip fix via constraints) 
961  ## | 4:      -32110                                           |                              ## 
962  ## |==========================================================| [ DFP       |v] [ Optimize ] ## 
963   
964  # repair or lose pascal max*ll; harmonize with QuB  # was: fix reduceconstraints in pascal (what else has changed?) 
965  # G/Vleak sensitive only for IeqFv models? 
966   
967  # other outputs (tabular results?) 
968  # averagecount, variance of current 
969  # fast/slow numgrad? 
970  # useVar always? 
971