lenstronomy.Sampling package

Submodules

lenstronomy.Sampling.likelihood module

class lenstronomy.Sampling.likelihood.LikelihoodModule(kwargs_data_joint, kwargs_model, param_class, image_likelihood=True, check_bounds=True, check_matched_source_position=False, astrometric_likelihood=False, image_position_likelihood=False, source_position_likelihood=False, image_position_uncertainty=0.004, check_positive_flux=False, source_position_tolerance=0.001, source_position_sigma=0.001, force_no_add_image=False, source_marg=False, linear_prior=None, restrict_image_number=False, max_num_images=None, bands_compute=None, time_delay_likelihood=False, image_likelihood_mask_list=None, flux_ratio_likelihood=False, kwargs_flux_compute={}, prior_lens=[], prior_source=[], prior_extinction=[], prior_lens_light=[], prior_ps=[], prior_special=[], prior_lens_kde=[], prior_source_kde=[], prior_lens_light_kde=[], prior_ps_kde=[], prior_special_kde=[], prior_extinction_kde=[], prior_lens_lognormal=[], prior_source_lognormal=[], prior_extinction_lognormal=[], prior_lens_light_lognormal=[], prior_ps_lognormal=[], prior_special_lognormal=[], custom_logL_addition=None, kwargs_pixelbased=None)[source]

Bases: object

this class contains the routines to run a MCMC process the key components are: - imSim_class: an instance of a class that simulates one (or more) images and returns the likelihood, such as ImageModel(), Multiband(), MultiExposure() - param_class: instance of a Param() class that can cast the sorted list of parameters that are sampled into the conventions of the imSim_class

Additional arguments are supported for adding a time-delay likelihood etc (see __init__ definition)

static check_bounds(args, lowerLimit, upperLimit, verbose=False)[source]

checks whether the parameter vector has left its bound, if so, adds a big number

effective_num_data_points(**kwargs)[source]

returns the effective number of data points considered in the X2 estimation to compute the reduced X2 value

likelihood(a)[source]
logL(args, verbose=False)[source]

routine to compute X2 given variable parameters for a MCMC/PSO chain

log_likelihood(kwargs_return, verbose=False)[source]
negativelogL(a)[source]

for minimizer function, the negative value of the logl value is requested

Parameters:a – array of parameters
Returns:-logL
num_data
Returns:number of independent data points in the combined fitting
param_limits

lenstronomy.Sampling.parameters module

class lenstronomy.Sampling.parameters.Param(kwargs_model, kwargs_fixed_lens=None, kwargs_fixed_source=None, kwargs_fixed_lens_light=None, kwargs_fixed_ps=None, kwargs_fixed_special=None, kwargs_fixed_extinction=None, kwargs_lower_lens=None, kwargs_lower_source=None, kwargs_lower_lens_light=None, kwargs_lower_ps=None, kwargs_lower_special=None, kwargs_lower_extinction=None, kwargs_upper_lens=None, kwargs_upper_source=None, kwargs_upper_lens_light=None, kwargs_upper_ps=None, kwargs_upper_special=None, kwargs_upper_extinction=None, kwargs_lens_init=None, linear_solver=True, joint_lens_with_lens=[], joint_lens_light_with_lens_light=[], joint_source_with_source=[], joint_lens_with_light=[], joint_source_with_point_source=[], joint_lens_light_with_point_source=[], joint_extinction_with_lens_light=[], joint_lens_with_source_light=[], mass_scaling_list=None, point_source_offset=False, num_point_source_list=None, image_plane_source_list=None, solver_type='NONE', Ddt_sampling=None, source_size=False, num_tau0=0, lens_redshift_sampling_indexes=None, source_redshift_sampling_indexes=None, source_grid_offset=False, num_shapelet_lens=0, log_sampling_lens=[])[source]

Bases: object

class that handles the parameter constraints. In particular when different model profiles share joint constraints.

Options between same model classes:

‘joint_lens_with_lens’:list [[i_lens, k_lens, [‘param_name1’, ‘param_name2’, …]], […], …], joint parameter between two lens models

‘joint_lens_light_with_lens_light’:list [[i_lens_light, k_lens_light, [‘param_name1’, ‘param_name2’, …]], […], …], joint parameter between two lens light models, the second adopts the value of the first

‘joint_source_with_source’:list [[i_source, k_source, [‘param_name1’, ‘param_name2’, …]], […], …], joint parameter between two source surface brightness models, the second adopts the value of the first

Options between different model classes:

‘joint_lens_with_light’: list [[i_light, k_lens, [‘param_name1’, ‘param_name2’, …]], […], …], joint parameter between lens model and lens light model

‘joint_source_with_point_source’: list [[i_point_source, k_source], […], …], joint position parameter between lens model and source light model

‘joint_lens_light_with_point_source’: list [[i_point_source, k_lens_light], […], …], joint position parameter between lens model and lens light model

‘joint_extinction_with_lens_light’: list [[i_lens_light, k_extinction, [‘param_name1’, ‘param_name2’, …]], […], …], joint parameters between the lens surface brightness and the optical depth models

‘joint_lens_with_source_light’: [[i_source, k_lens, [‘param_name1’, ‘param_name2’, …]], […], …], joint parameter between lens model and source light model. Samples light model parameter only.

hierarchy is as follows: 1. Point source parameters are inferred 2. Lens light joint parameters are set 3. Lens model joint constraints are set 4. Lens model solver is applied 5. Joint source and point source is applied

Alternatively to the format of the linking of parameters with IDENTICAL names as listed above as: [[i_1, k_2, [‘param_name1’, ‘param_name2’, …]], […], …] the following format of the arguments are supported to join parameters with DIFFERENT names: [[i_1, k_2, {‘param_old1’: ‘param_new1’, ‘ra_0’: ‘center_x’}], […], …] Log10 sampling of the lens parameters : ‘log_sampling_lens’: [[i_lens, [‘param_name1’, ‘param_name2’, …]], […], …], Sample the log10 of the lens model parameters.

args2kwargs(args, bijective=False)[source]
Parameters:
  • args – tuple of parameter values (float, strings, …)
  • bijective – boolean, if True (default) returns the parameters in the form as they are sampled (e.g. if image_plane_source_list is set =True it returns the position in the image plane coordinates), if False, returns the parameters in the form to render a model (e.g. image_plane_source_list positions are ray-traced back to the source plane).
Returns:

keyword arguments sorted in lenstronomy conventions

check_solver(kwargs_lens, kwargs_ps)[source]

test whether the image positions map back to the same source position :param kwargs_lens: lens model keyword argument list :param kwargs_ps: point source model keyword argument list :return: Euclidean distance between the ray-shooting of the image positions

image2source_plane(kwargs_source, kwargs_lens, image_plane=False)[source]

maps the image plane position definition of the source plane

Parameters:
  • kwargs_source – source light model keyword argument list
  • kwargs_lens – lens model keyword argument list
  • image_plane – boolean, if True, does not up map image plane parameters to source plane
Returns:

source light model keyword arguments with mapped position arguments from image to source plane

kwargs2args(kwargs_lens=None, kwargs_source=None, kwargs_lens_light=None, kwargs_ps=None, kwargs_special=None, kwargs_extinction=None)[source]

inverse of getParam function :param kwargs_lens: keyword arguments depending on model options :param kwargs_source: keyword arguments depending on model options :param kwargs_lens_light: lens light model keyword argument list :param kwargs_ps: point source model keyword argument list :param kwargs_special: special keyword arguments :param kwargs_extinction: extinction model keyword argument list :return: numpy array of parameters

num_param()[source]
Returns:number of parameters involved (int), list of parameter names
num_param_linear()[source]
Returns:number of linear basis set coefficients that are solved for
num_point_source_images
param_limits()[source]
Returns:lower and upper limits of the arguments being sampled
print_setting()[source]

prints the setting of the parameter class

Returns:
update_kwargs_model(kwargs_special)[source]

updates model keyword arguments with redshifts being sampled

Parameters:kwargs_special – keyword arguments from SpecialParam() class return of sampling arguments
Returns:kwargs_model, bool (True if kwargs_model has changed, else False)
update_lens_scaling(kwargs_special, kwargs_lens, inverse=False)[source]

multiplies the scaling parameters of the profiles

Parameters:
  • kwargs_special – keyword arguments of the ‘special’ arguments
  • kwargs_lens – lens model keyword argument list
  • inverse – bool, if True, performs the inverse lens scaling for bijective transforms
Returns:

updated lens model keyword argument list

lenstronomy.Sampling.sampler module

class lenstronomy.Sampling.sampler.Sampler(likelihoodModule)[source]

Bases: object

class which executes the different sampling methods Available are: MCMC with emcee and comsoHammer and a Particle Swarm Optimizer. This are examples and depending on your problem, you might find other/better solutions. Feel free to sample with your convenient sampler!

mcmc_emcee(n_walkers, n_run, n_burn, mean_start, sigma_start, mpi=False, progress=False, threadCount=1, initpos=None, backup_filename=None, start_from_backup=False)[source]

Run MCMC with emcee. For details, please have a look at the documentation of the emcee packager.

Parameters:
  • n_walkers (integer) – number of walkers in the emcee process
  • n_run (integer) – number of sampling (after burn-in) of the emcee
  • n_burn (integer) – number of burn-in iterations (those will not be saved in the output sample)
  • mean_start (numpy array of length the number of parameters) – mean of the parameter position of the initialising sample
  • sigma_start (numpy array of length the number of parameters) – spread of the parameter values (uncorrelated in each dimension) of the initialising sample
  • mpi (bool) – if True, initializes an MPIPool to allow for MPI execution of the sampler
  • progress (bool) – if True, prints the progress bar
  • threadCount (integer) – number of threats in multi-processing (not applicable for MPI)
  • initpos (numpy array of size num param x num walkser) – initial walker position to start sampling (optional)
  • backup_filename (string) – name of the HDF5 file where sampling state is saved (through emcee backend engine)
  • start_from_backup (bool) – if True, start from the state saved in backup_filename. Otherwise, create a new backup file with name backup_filename (any already existing file is overwritten!).
Returns:

samples, ln likelihood value of samples

Return type:

numpy 2d array, numpy 1d array

pso(n_particles, n_iterations, lower_start=None, upper_start=None, threadCount=1, init_pos=None, mpi=False, print_key='PSO')[source]

Return the best fit for the lens model on catalogue basis with particle swarm optimizer.

Parameters:
  • n_particles – number of particles in the sampling process
  • n_iterations – number of iterations of the swarm
  • lower_start – numpy array, lower end parameter of the values of the starting particles
  • upper_start – numpy array, upper end parameter of the values of the starting particles
  • threadCount – number of threads in the computation (only applied if mpi=False)
  • init_pos – numpy array, position of the initial best guess model
  • mpi – bool, if True, makes instance of MPIPool to allow for MPI execution
  • print_key – string, prints the process name in the progress bar (optional)
Returns:

kwargs_result (of best fit), [lnlikelihood of samples, positions of samples, velocity of sampels)

simplex(init_pos, n_iterations, method, print_key='SIMPLEX')[source]
Parameters:
  • init_pos – starting point for the optimization
  • n_iterations – maximum number of iterations
  • method – the optimization method, default is ‘Nelder-Mead’
Returns:

the best fit for the lens model using the optimization routine specified by method

lenstronomy.Sampling.special_param module

class lenstronomy.Sampling.special_param.SpecialParam(Ddt_sampling=False, mass_scaling=False, num_scale_factor=1, kwargs_fixed=None, kwargs_lower=None, kwargs_upper=None, point_source_offset=False, source_size=False, num_images=0, num_tau0=0, num_z_sampling=0, source_grid_offset=False)[source]

Bases: object

class that handles special parameters that are not directly part of a specific model component. These includes cosmology relevant parameters, astrometric errors and overall scaling parameters.

get_params(args, i)[source]
Parameters:
  • args – argument list
  • i – integer, list index to start the read out for this class
Returns:

keyword arguments related to args, index after reading out arguments of this class

num_param()[source]
Returns:integer, number of free parameters sampled (and managed) by this class, parameter names (list of strings)
set_params(kwargs_special)[source]
Parameters:kwargs_special – keyword arguments with parameter settings
Returns:argument list of the sampled parameters extracted from kwargs_special

Module contents