Merge pull request #113 from choderalab/quickmodel

Further small updates to examples, quickmodel, and pymcmodels

AssayTools/parser.py

@@ -83,7 +83,7 @@ def get_data_using_inputs(inputs):
                 except:
                     pass
 
-        for i in range(0,7,2):
+        for i in range(0,len(inputs['ligand_order']*2),2):
             protein_row = ALPHABET[i]
             buffer_row = ALPHABET[i+1]
 

AssayTools/pymcmodels.py

@@ -779,8 +779,8 @@ def run_mcmc(pymc_model, nthin=20, nburn=0, niter=20000, map=True, db='ram', dbn
         mcmc.use_step_method(pymc.Metropolis, getattr(pymc_model, 'epsilon_ex'), proposal_sd=10000.0, proposal_distribution='Normal')
     if hasattr(pymc_model, 'epsilon_em'):
         mcmc.use_step_method(pymc.Metropolis, getattr(pymc_model, 'epsilon_em'), proposal_sd=10000.0, proposal_distribution='Normal')
-    # log_F_PL is highly correlated with DeltaG
-    mcmc.use_step_method(pymc.AdaptiveMetropolis, [pymc_model.log_F_PL, pymc_model.DeltaG], scales={ pymc_model.log_F_PL : 0.1, pymc_model.DeltaG : 0.1  })
+    # Uncomment below to use log_F_PL highly correlated with DeltaG to improve sampling somewhat.
+    #mcmc.use_step_method(pymc.AdaptiveMetropolis, [pymc_model.log_F_PL, pymc_model.DeltaG], scales={ pymc_model.log_F_PL : 0.1, pymc_model.DeltaG : 0.1  })
 
     mcmc.sample(iter=(nburn+niter), burn=nburn, thin=nthin, progress_bar=False, tune_throughout=True)
 

examples/direct-fluorescence-assay/README.md

@@ -23,11 +23,11 @@ To plot raw data using `xml2png`, navigate into the data folder and run:
 
 To get a deltaG estimate from direct fluorescence data using our Bayesian modeling framework, via the `quickmodel` script run:
 
-  `env PYTHONPATH="./" quickmodel --inputs 'inputs_p38_singlet'`
+  `quickmodel --inputs 'inputs_p38_singlet'`
 
 or
 
-  `env PYTHONPATH="./" quickmodel --inputs 'inputs_Gef_WIP' --type 'spectra'`
+  `quickmodel --inputs 'inputs_Gef_WIP'`
 
 The two inputs.py files here are required to run these, but they can be changed if you want to run other types of analyses (different sections of the xml file, different wavelenghts, etc.).
 

scripts/inputs_example.py

@@ -10,6 +10,8 @@
     'ligand_order'  :  ['Bosutinib','Bosutinib Isomer','Erlotinib','Gefitinib','Bosutinib','Bosutinib Isomer','Erlotinib','Gefitinib'],
     'Lstated'       :  np.array([20.0e-6,14.0e-6,9.82e-6,6.88e-6,4.82e-6,3.38e-6,2.37e-6,1.66e-6,1.16e-6,0.815e-6,0.571e-6,0.4e-6,0.28e-6,0.196e-6,0.138e-6,0.0964e-6,0.0676e-6,0.0474e-6,0.0320e-6,0.0240e-6,0.0160e-6,0.0120e-6,0.008e-6,0.0], np.float64), # ligand concentration, M
     'Pstated'       :  0.5e-6 * np.ones([24],np.float64), # protein concentration, M
+    'P_error'       :  0.15,
+    'L_error'       :  0.08,
     'assay_volume'  :  50e-6, # assay volume, L
     'well_area'     :  0.1369, # well area, cm^2 for 4ti-0203 [http://4ti.co.uk/files/3113/4217/2464/4ti-0201.pdf]
     }

scripts/quickmodel.py

@@ -46,8 +46,15 @@ def quick_model(inputs, nsamples=1000, nthin=20):
             #these need to be changed so they are TAKEN FROM INPUTS!!!
 
             # Uncertainties in protein and ligand concentrations.
-            dPstated = 0.35 * inputs['Pstated'] # protein concentration uncertainty
-            dLstated = 0.08 * inputs['Lstated'] # ligand concentraiton uncertainty (due to gravimetric preparation and HP D300 dispensing)
+            try:
+                 dPstated = inputs['P_error'] * inputs['Pstated'] # protein concentration uncertainty               
+            except:
+                 dPstated = 0.35 * inputs['Pstated'] # protein concentration uncertainty
+
+            try:
+                 dLstated = inputs['L_error'] * inputs['Lstated']
+            except:
+                 dLstated = 0.08 * inputs['Lstated'] # ligand concentraiton uncertainty (due to gravimetric preparation and HP D300 dispensing)
 
             from assaytools import pymcmodels
             pymc_model = pymcmodels.make_model(inputs['Pstated'], dPstated, inputs['Lstated'], dLstated,