Package prompt :: Module monte_carlo
[hide private]
[frames] | no frames]

Source Code for Module prompt.monte_carlo

  1  ############################################################################### 
  2  #                                                                             # 
  3  # Copyright (C) 2004-2005 Edward d'Auvergne                                   # 
  4  #                                                                             # 
  5  # This file is part of the program relax.                                     # 
  6  #                                                                             # 
  7  # relax is free software; you can redistribute it and/or modify               # 
  8  # it under the terms of the GNU General Public License as published by        # 
  9  # the Free Software Foundation; either version 2 of the License, or           # 
 10  # (at your option) any later version.                                         # 
 11  #                                                                             # 
 12  # relax is distributed in the hope that it will be useful,                    # 
 13  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
 14  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
 15  # GNU General Public License for more details.                                # 
 16  #                                                                             # 
 17  # You should have received a copy of the GNU General Public License           # 
 18  # along with relax; if not, write to the Free Software                        # 
 19  # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA   # 
 20  #                                                                             # 
 21  ############################################################################### 
 22   
 23  import sys 
 24   
 25  import help 
 26   
 27   
28 -class Monte_carlo:
29 - def __init__(self, relax):
30 # Help. 31 self.__relax_help__ = \ 32 """Class containing the functions for Monte Carlo and related simulations.""" 33 34 # Add the generic help string. 35 self.__relax_help__ = self.__relax_help__ + "\n" + help.relax_class_help 36 37 # Place relax in the class namespace. 38 self.__relax__ = relax
39 40
41 - def create_data(self, run=None, method='back_calc'):
42 """Function for creating simulation data. 43 44 Keyword Arguments 45 ~~~~~~~~~~~~~~~~~ 46 47 run: The name of the run. 48 49 method: The simulation method. 50 51 52 Description 53 ~~~~~~~~~~~ 54 55 The method argument can either be set to 'back_calc' or 'direct', the choice of which 56 determines the simulation type. If the values or parameters of a run are calculated rather 57 than minimised, this option will have no effect, hence, 'back_calc' and 'direct' are 58 identical. 59 60 For error analysis, the method argument should be set to 'back_calc' which will result in 61 proper Monte Carlo simulations. The data used for each simulation is back calculated from 62 the minimised model parameters and is randomised using Gaussian noise where the standard 63 deviation is from the original error set. When the method is set to 'back_calc', this 64 function should only be called after the model or run is fully minimised. 65 66 The simulation type can be changed by setting the method argument to 'direct'. This will 67 result in simulations which cannot be used in error analysis and which are no longer Monte 68 Carlo simulations. However, these simulations are required for certain model selection 69 techniques (see the documentation for the model selection function for details), and can be 70 used for other purposes. Rather than the data being back calculated from the fitted model 71 parameters, the data is generated by taking the original data and randomising using Gaussian 72 noise with the standard deviations set to the original error set. 73 """ 74 75 # Function intro text. 76 if self.__relax__.interpreter.intro: 77 text = sys.ps3 + "monte_carlo.create_data(" 78 text = text + "run=" + `run` 79 text = text + ", method=" + `method` + ")" 80 print text 81 82 # The run argument. 83 if type(run) != str: 84 raise RelaxStrError, ('run', run) 85 86 # The method argument. 87 if type(method) != str: 88 raise RelaxStrError, ('method', method) 89 90 # Execute the functional code. 91 self.__relax__.generic.monte_carlo.create_data(run=run, method=method)
92 93
94 - def error_analysis(self, run=None, prune=0.0):
95 """Function for calculating parameter errors from the Monte Carlo simulations. 96 97 Keyword Arguments 98 ~~~~~~~~~~~~~~~~~ 99 100 run: The name of the run. 101 102 prune: Legacy argument corresponding to 'trim' in Art Palmer's Modelfree program. 103 104 105 Description 106 ~~~~~~~~~~~ 107 108 Parameter errors are calculated as the standard deviation of the distribution of parameter 109 values. This function should never be used if parameter values are obtained by minimisation 110 and the simulation data are generated using the method 'direct'. The reason is because only 111 true Monte Carlo simulations can give the true parameter errors. 112 113 The prune argument is legacy code which corresponds to the 'trim' option in Art Palmer's 114 Modelfree program. To remove failed simulations, the eliminate function should be used 115 prior to this function. Eliminating the simulations specifically identifies and removes the 116 failed simulations whereas the prune argument will only, in a few cases, positively identify 117 failed simulations but only if severe parameter limits have been imposed. Most failed 118 models will pass through the pruning process and hence cause a catastrophic increase in the 119 parameter errors. If the argument must be used, the following must be taken into account. 120 If the values or parameters of a run are calculated rather than minimised, the prune 121 argument must be set to zero. The value of this argument is proportional to the number of 122 simulations removed prior to error calculation. If prune is set to 0.0, all simulations are 123 used for calculating errors, whereas a value of 1.0 excludes all data. In almost all cases 124 prune must be set to zero, any value greater than zero will result in an underestimation of 125 the error values. If a value is supplied, the lower and upper tails of the distribution of 126 chi-squared values will be excluded from the error calculation. 127 """ 128 129 # Function intro text. 130 if self.__relax__.interpreter.intro: 131 text = sys.ps3 + "monte_carlo.error_analysis(" 132 text = text + "run=" + `run` 133 text = text + ", prune=" + `prune` + ")" 134 print text 135 136 # The run argument. 137 if type(run) != str: 138 raise RelaxStrError, ('run', run) 139 140 # The prune argument. 141 if type(prune) != int and type(prune) != float: 142 raise RelaxNumError, ('prune', prune) 143 144 # Execute the functional code. 145 self.__relax__.generic.monte_carlo.error_analysis(run=run, prune=prune)
146 147
148 - def initial_values(self, run=None):
149 """Function for setting the initial simulation parameter values. 150 151 Keyword Arguments 152 ~~~~~~~~~~~~~~~~~ 153 154 run: The name of the run. 155 156 157 Description 158 ~~~~~~~~~~~ 159 160 This function only effects runs where minimisation occurs and can therefore be skipped if 161 the values or parameters of a run are calculated rather than minimised. However, if 162 accidentally run in this case, the results will be unaffected. It should only be called 163 after the model or run is fully minimised. Once called, the functions 'grid_search' and 164 'minimise' will only effect the simulations and not the model parameters. 165 166 The initial values of the parameters for each simulation is set to the minimised parameters 167 of the model. A grid search can be undertaken for each simulation instead, although this 168 is computationally expensive and unnecessary. The minimisation function should be executed 169 for a second time after running this function. 170 """ 171 172 # Function intro text. 173 if self.__relax__.interpreter.intro: 174 text = sys.ps3 + "monte_carlo.initial_values(" 175 text = text + "run=" + `run` + ")" 176 print text 177 178 # The run argument. 179 if type(run) != str: 180 raise RelaxStrError, ('run', run) 181 182 # Execute the functional code. 183 self.__relax__.generic.monte_carlo.initial_values(run=run)
184 185
186 - def off(self, run=None):
187 """Function for turning simulations off. 188 189 Keyword Arguments 190 ~~~~~~~~~~~~~~~~~ 191 192 run: The name of the run. 193 """ 194 195 # Function intro text. 196 if self.__relax__.interpreter.intro: 197 text = sys.ps3 + "monte_carlo.off(" 198 text = text + "run=" + `run` + ")" 199 print text 200 201 # The run argument. 202 if type(run) != str: 203 raise RelaxStrError, ('run', run) 204 205 # Execute the functional code. 206 self.__relax__.generic.monte_carlo.off(run=run)
207 208
209 - def on(self, run=None):
210 """Function for turning simulations on. 211 212 Keyword Arguments 213 ~~~~~~~~~~~~~~~~~ 214 215 run: The name of the run. 216 """ 217 218 # Function intro text. 219 if self.__relax__.interpreter.intro: 220 text = sys.ps3 + "monte_carlo.on(" 221 text = text + "run=" + `run` + ")" 222 print text 223 224 # The run argument. 225 if type(run) != str: 226 raise RelaxStrError, ('run', run) 227 228 # Execute the functional code. 229 self.__relax__.generic.monte_carlo.on(run=run)
230 231
232 - def setup(self, run=None, number=500):
233 """Function for setting up Monte Carlo simulations. 234 235 Keyword Arguments 236 ~~~~~~~~~~~~~~~~~ 237 238 run: The name of the run. 239 240 number: The number of Monte Carlo simulations. 241 242 243 Description 244 ~~~~~~~~~~~ 245 246 This function must be called prior to any of the other Monte Carlo functions. The effect is 247 that the number of simulations for the given run will be set and that simulations will be 248 turned on. 249 """ 250 251 # Function intro text. 252 if self.__relax__.interpreter.intro: 253 text = sys.ps3 + "monte_carlo.setup(" 254 text = text + "run=" + `run` 255 text = text + ", number=" + `number` + ")" 256 print text 257 258 # The run argument. 259 if type(run) != str: 260 raise RelaxStrError, ('run', run) 261 262 # The number of simulations. 263 if type(number) != int: 264 raise RelaxIntError, ('number', number) 265 266 # Execute the functional code. 267 self.__relax__.generic.monte_carlo.setup(run=run, number=number)
268 269 270 271 # Modify all docstrings. 272 ######################## 273 274 __description__ = """ 275 Monte Carlo Simulation Overview 276 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 277 278 For proper error analysis using Monte Carlo simulations, a sequence of function calls is 279 required for running the various simulation components. The steps necessary for 280 implementing Monte Carlo simulations are: 281 282 1. The measured data set together with the corresponding error set should be loaded into 283 relax. 284 285 2. Either minimisation is used to optimise the parameters of the chosen model, or a 286 calculation is run. 287 288 3. To initialise and turn on Monte Carlo simulations, the number of simulations, n, needs 289 to be set. 290 291 4. The simulation data needs to be created either by back calculation from the fully 292 minimised model parameters from step 2 or by direct calculation when values are calculated 293 rather than minimised. The error set is used to randomise each simulation data set by 294 assuming Gaussian errors. This creates a synthetic data set for each Monte Carlo 295 simulation. 296 297 5. Prior to minimisation of the parameters of each simulation, initial parameter estimates 298 are required. These are taken as the optimised model parameters. An alternative is to use 299 a grid search for each simulation to generate initial estimates, however this is extremely 300 computationally expensive. For the case where values are calculated rather than minimised, 301 this step should be skipped (although the results will be unaffected if this is accidentally 302 run). 303 304 6. Each simulation requires minimisation or calculation. The same techniques as used in 305 step 2, excluding the grid search when minimising, should be used for the simulations. 306 307 7. Failed simulations are removed using the techniques of model elimination. 308 309 8. The model parameter errors are calculated from the distribution of simulation 310 parameters. 311 312 313 Monte Carlo simulations can be turned on or off using functions within this class. Once the 314 function for setting up simulations has been called, simulations will be turned on. The 315 effect of having simulations turned on is that the functions used for minimisation (grid 316 search, minimise, etc) or calculation will only affect the simulation parameters and not the 317 model parameters. By subsequently turning simulations off using the appropriate function, 318 the functions used in minimisation will affect the model parameters and not the simulation 319 parameters. 320 321 322 An example, for model-free analysis, which includes only the functions required for 323 implementing the above steps is: 324 325 relax> grid_search('m1', inc=11) # Step 2. 326 relax> minimise('newton', run='m1') # Step 2. 327 relax> monte_carlo.setup('m1', number=500) # Step 3. 328 relax> monte_carlo.create_data('m1', method='back_calc') # Step 4. 329 relax> monte_carlo.initial_values('m1') # Step 5. 330 relax> minimise('newton', run='m1') # Step 6. 331 relax> eliminate('m1') # Step 7. 332 relax> monte_carlo.error_analysis('m1') # Step 8. 333 334 An example for reduced spectral density mapping is: 335 336 relax> calc('600MHz') # Step 2. 337 relax> monte_carlo.setup('600MHz', number=500) # Step 3. 338 relax> monte_carlo.create_data('600MHz', method='back_calc') # Step 4. 339 relax> calc('600MHz') # Step 6. 340 relax> monte_carlo.error_analysis('600MHz') # Step 8. 341 """ 342 343 create_data.__doc__ = create_data.__doc__ + "\n\n" + __description__ 344 error_analysis.__doc__ = error_analysis.__doc__ + "\n\n" + __description__ 345 initial_values.__doc__ = initial_values.__doc__ + "\n\n" + __description__ 346 off.__doc__ = off.__doc__ + "\n\n" + __description__ 347 on.__doc__ = on.__doc__ + "\n\n" + __description__ 348 setup.__doc__ = setup.__doc__ + "\n\n" + __description__
349