Note: This page was generated from a Jupyter notebook which can be found at docs/tutorial_notebooks/5_parameter_space.ipynb in your DYNAMITE directory. —-

5. Parameter Space

To run a DYNAMITE model, one must specify a number of parameters for the gravitational potential. The aim of this notebook is to demonstrate how to specify these parameters and to highlight features that we have implemented in order to help you explore parameter space.

We’ll start as before by reading the same configuration file as previously,

import dynamite as dyn

# read the config file
fname = 'NGC6278_config.yaml'
c = dyn.config_reader.Configuration(fname, reset_logging=True)

[INFO] 21:49:20 - dynamite.config_reader.Configuration - Config file NGC6278_config.yaml read.
[INFO] 21:49:20 - dynamite.config_reader.Configuration - io_settings...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - Output directory tree: NGC6278_output/.
[INFO] 21:49:20 - dynamite.config_reader.Configuration - system_attributes...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - model_components...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - system_parameters...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - orblib_settings...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - weight_solver_settings...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - Will attempt to recover partially run models.
[INFO] 21:49:20 - dynamite.config_reader.Configuration - parameter_space_settings...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - multiprocessing_settings...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - ... using 4 CPUs for orbit integration.
[INFO] 21:49:20 - dynamite.config_reader.Configuration - ... using 4 CPUs for weight solving.
[INFO] 21:49:20 - dynamite.config_reader.Configuration - legacy_settings...
[INFO] 21:49:20 - dynamite.config_reader.Configuration - System assembled
[INFO] 21:49:20 - dynamite.config_reader.Configuration - Configuration validated
[INFO] 21:49:20 - dynamite.config_reader.Configuration - Instantiated parameter space
[INFO] 21:49:20 - dynamite.model.AllModels - Previous models have been found: Reading NGC6278_output/all_models.ecsv into AllModels.table
[INFO] 21:49:20 - dynamite.config_reader.Configuration - Instantiated AllModels object
[INFO] 21:49:20 - dynamite.model.AllModels - No all_models table update required.

When the configuration object is created, internally, a parameter space object is created. This parspace object is a list, and every entry of this list is a parameter in our model, Lets extract this and have a look

# extract the parameter space
parspace = c.parspace
print('type of parspace is', type(parspace))
print('length of parspace is', len(parspace))
print('the parameter names are:')
for par in parspace:
    print('   -', par.name)

type of parspace is <class 'dynamite.parameter_space.ParameterSpace'>
length of parspace is 8
the parameter names are:
   - m-bh
   - a-bh
   - c-dh
   - f-dh
   - q-stars
   - p-stars
   - u-stars
   - ml

Several properties are specified for each parameter in the configuration file. Let’s look at the value,

print('Parameter / value in config file:')
for par in c.parspace:
    print(f'   {par.name} = {par.raw_value}')

Parameter / value in config file:
   m-bh = 5.0
   a-bh = 0.001
   c-dh = 8.0
   f-dh = 1.0
   q-stars = 0.54
   p-stars = 0.99
   u-stars = 0.9999
   ml = 5.0

These are the starting values from which we would like to run a model.

One complication in specifying these values is that, for some parameters, we would like to take logarithmically spaced steps through parameter space, i.e. the ones which are specificed as

parameters -> XXX -> logarithmic : True

Logarithmic spacing can be useful for mass parameters. For other parameters (e.g. length scales) linearly spaced steps may be more appropriate. For other types of parameters (e.g. angles) a different spacing altogether may be preferable.

To handle these possibilities, we introduce the concept of raw parameter values, distinct from the values themselves. All values associated with parameters in the configuration file are given in raw units. When we step through parameter space, we take linear steps in raw values. The conversion from raw values to the parameter values is handled by the Parameter class and the parameter values are accessible via the

Parameter.par_value

property. So to convert the above list from raw values, we can do the following,

print('Parameter / value in linear units:')
for par in c.parspace:
    print(f'   {par.name} = {par.par_value}')

Parameter / value in linear units:
   m-bh = 100000.0
   a-bh = 0.001
   c-dh = 8.0
   f-dh = 10.0
   q-stars = 0.54
   p-stars = 0.99
   u-stars = 0.9999
   ml = 5.0

Notice how only those parameters which have been specified with logarithmic : True have been modified.

Another property that we specifie for each parameter is whether or not it is fixed, a boolean value,

for par in parspace:
    if par.fixed:
        fix_string = ' is fixed.'
    if not par.fixed:
        fix_string = ' is NOT fixed.'
    print(f'{par.name}{fix_string}')

m-bh is fixed.
a-bh is fixed.
c-dh is fixed.
f-dh is NOT fixed.
q-stars is fixed.
p-stars is fixed.
u-stars is fixed.
ml is NOT fixed.

The only parameters which are not fixed for this example are the dark matter fraction f-dh and the mass-to-light ratio ml. For these free parameters, additional properties about how search through parameter space are stored in the par_generator_settings attribute,

for par in parspace:
    if not par.fixed:
        tmp = par.par_generator_settings
        lo, hi, step = tmp['lo'], tmp['hi'], tmp['step']
        print(f'{par.name} takes step-size {step} and bounds ({lo,hi})')

f-dh takes step-size 0.5 and bounds ((-1.5, 1.5))
ml takes step-size 4.0 and bounds ((1.0, 9.0))

How do we search over these free parameters? Running models (especially calcuating the orbit library) is expensive, so we will want to search through parameter space in the most efficient way possible.

In general, an algorithm to search through parameter space will take as input 1. the output of all models which have been run so far (e.g. \(\chi^2\) values) 2. setting for the free parameters (e.g. step-size and bounds) The algorithm will then output a new list of parameters for which we want to run models.

In DYNAMITE, we implement this generic idea in the class dyn.parameter_space.ParameterGenerator. In the configuration file, you specify which parameter generator you would like to use, at the location

parameter_space_settings -> generator_type

The current choice is

c.settings.parameter_space_settings['generator_type']

'LegacyGridSearch'

This parameter generator requires an additional setting which is set at

parameter_space_settings -> generator_settings -> threshold_del_chi2_abs

or

parameter_space_settings -> generator_settings -> threshold_del_chi2_as_frac_of_sqrt2nobs

(the options are mutually exclusive, set one or the other). Internally, the setting is converted to the appropriate threshold_del_chi2 and is accessed in the following way,

threshold_del_chi2_as_frac_of_sqrt2nobs = \
    c.settings.parameter_space_settings['generator_settings']['threshold_del_chi2_as_frac_of_sqrt2nobs']
threshold_del_chi2 = c.settings.parameter_space_settings['generator_settings']['threshold_del_chi2']
print(f'threshold_del_chi2_as_frac_of_sqrt2nobs = {threshold_del_chi2_as_frac_of_sqrt2nobs}')
print(f'threshold_del_chi2 = {threshold_del_chi2}')

threshold_del_chi2_as_frac_of_sqrt2nobs = 0.1
threshold_del_chi2 = 3.4871191548325395

The algorithm implemented to generate parameters in LegacyGridSearch is the following,

iteration = 0
if iteration == 0
    all parameters take `value` specified in the config
else:
    1. find the model with the lowest chi-squared
    2. find all models with chi-squared within threshold_del_chi2 of the lowest value
    3. for all models satisfying that criteria:
        - for all free parameters:
            - generate a new parameter set +/-1 step-size from the current value
    4. Remove any models with parameters outside specified bounds
    5. iteration = iteration + 1
stop if no new models are added, or any other stopping criteria are met

For those of you who have used the previous version of the trixial Schwarzschild modelling code (aka schwpy), this is the same algorithm which was implemented there.

The last line of the algorithm mentions stopping criteria. Settings which control the stopping criteria are also speicified in the configuration file, under

parameter_space_settings -> stopping_criteria

The current settings which are the following,

stopping_crierita = c.settings.parameter_space_settings['stopping_criteria']
for key in stopping_crierita:
    print(f'{key} = {stopping_crierita[key]}')

min_delta_chi2_abs = 0.5
n_max_mods = 15
n_max_iter = 4

These have the following meaning,

• if no new model impoves the chi-squared by at least min_delta_chi2, then stop

• if we have already run n_max_mods models, then stop

• if we have already run n_max_iter iterations, then stop

:)