dials.algorithms.symmetry

Methods for symmetry determination.

This module provides a base class for symmetry determination algorithms.

dials.algorithms.symmetry.median_unit_cell(experiments)[source]
dials.algorithms.symmetry.resolution_filter_from_array(intensities, min_i_mean_over_sigma_mean, min_cc_half)[source]

Run the resolution filter using miller array data format.

dials.algorithms.symmetry.resolution_filter_from_reflections_experiments(reflections, experiments, min_i_mean_over_sigma_mean, min_cc_half)[source]

Run the resolution filter using native dials data formats.

class dials.algorithms.symmetry.symmetry_base(intensities, normalisation='ml_aniso', lattice_symmetry_max_delta=2.0, d_min=<libtbx.AutoType object>, min_i_mean_over_sigma_mean=4, min_cc_half=0.6, relative_length_tolerance=None, absolute_angle_tolerance=None, best_monoclinic_beta=True)[source]

Bases: object

Base class for symmetry analysis.

__init__(intensities, normalisation='ml_aniso', lattice_symmetry_max_delta=2.0, d_min=<libtbx.AutoType object>, min_i_mean_over_sigma_mean=4, min_cc_half=0.6, relative_length_tolerance=None, absolute_angle_tolerance=None, best_monoclinic_beta=True)[source]

Initialise a symmetry_base object.

Parameters:
  • intensities (cctbx.miller.array) – The intensities on which to perform symmetry analysis.

  • normalisation (str) – The normalisation method to use. Possible choices are ‘kernel’, ‘quasi’, ‘ml_iso’ and ‘ml_aniso’. Set to None to switch off normalisation altogether.

  • lattice_symmetry_max_delta (float) – The maximum value of delta for determining the lattice symmetry using the algorithm of Le Page (1982).

  • d_min (float) – Optional resolution cutoff to be applied to the input intensities. If set to libtbx.Auto then d_min will be automatically determined according to the parameters min_i_mean_over_sigma_mean and min_cc_half.

  • min_i_mean_over_sigma_mean (float) – minimum value of \(|I|/|sigma(I)|\) for automatic determination of resolution cutoff.

  • min_cc_half (float) – minimum value of CC½ for automatic determination of resolution cutoff.

  • relative_length_tolerance (float) – Relative length tolerance in checking consistency of input unit cells against the median unit cell.

  • absolute_angle_tolerance (float) – Absolute angle tolerance in checking consistency of input unit cells against the median unit cell.

  • best_monoclinic_beta (bool) – If True, then for monoclinic centered cells, I2 will be preferred over C2 if it gives a less oblique cell (i.e. smaller beta angle).

static kernel_normalisation(intensities)[source]

Kernel normalisation of the input intensities.

Parameters:

intensities (cctbx.miller.array) – The intensities to be normalised.

Returns:

The normalised intensities.

Return type:

cctbx.miller.array

static ml_aniso_normalisation(intensities)[source]

Anisotropic maximum-likelihood normalisation of the input intensities.

Parameters:

intensities (cctbx.miller.array) – The intensities to be normalised.

Returns:

The normalised intensities.

Return type:

cctbx.miller.array

static ml_iso_normalisation(intensities)[source]

Isotropic maximum-likelihood normalisation of the input intensities.

Parameters:

intensities (cctbx.miller.array) – The intensities to be normalised.

Returns:

The normalised intensities.

Return type:

cctbx.miller.array

Algorithms for determination of Laue group symmetry.

class dials.algorithms.symmetry.laue_group.CorrelationCoefficientAccumulator(x=None, y=None)[source]

Bases: object

Class for incremental computation of correlation coefficients.

Uses the single-pass formula for Pearson correlation coefficient:

https://en.wikipedia.org/wiki/Pearson_correlation_coefficient#For_a_sample

__init__(x=None, y=None)[source]

Initialise a CorrelationCoefficientAccumulator object.

Parameters:
  • x (list) – Optional list of x values to initialise the accumulator.

  • y (list) – Optional list of y values to initialise the accumulator.

accumulate(x, y)[source]

Accumulate the x and y values provided.

Parameters:
  • x (list) – The list of x values to accumulate.

  • y (list) – The list of y values to accumulate.

coefficient()[source]

Calculate the correlation coefficient.

Returns:

The correlation coefficient.

Return type:

float

denominator()[source]

Calculate the denominator of the correlation coefficient formula.

\[\sqrt{n \sum{x^2} - \sum{x}^2} \sqrt{n \sum{y^2} - \sum{y}^2}\]
Returns:

The value of the denominator.

Return type:

float

n()[source]

Return the number of values contributing to the correlation coefficient.

Returns:

n (int)

numerator()[source]

Calculate the numerator of the correlation coefficient formula.

\[n \sum{x y} - \sum{x} \sum{y}\]
Returns:

The value of the numerator.

Return type:

float

class dials.algorithms.symmetry.laue_group.LaueGroupAnalysis(intensities, normalisation='ml_aniso', lattice_symmetry_max_delta=2.0, d_min=<libtbx.AutoType object>, min_i_mean_over_sigma_mean=4, min_cc_half=0.6, relative_length_tolerance=None, absolute_angle_tolerance=None, best_monoclinic_beta=True)[source]

Bases: symmetry_base

Determination of Laue group symmetry using algorithms similar to POINTLESS.

__init__(intensities, normalisation='ml_aniso', lattice_symmetry_max_delta=2.0, d_min=<libtbx.AutoType object>, min_i_mean_over_sigma_mean=4, min_cc_half=0.6, relative_length_tolerance=None, absolute_angle_tolerance=None, best_monoclinic_beta=True)[source]

Initialise a LaueGroupAnalysis object.

Parameters:
  • intensities (cctbx.miller.array) – The intensities on which to perform symmetry analysis.

  • normalisation (str) – The normalisation method to use. Possible choices are ‘kernel’, ‘quasi’, ‘ml_iso’ and ‘ml_aniso’. Set to None to switch off normalisation altogether.

  • lattice_symmetry_max_delta (float) – The maximum value of delta for determining the lattice symmetry using the algorithm of Le Page (1982).

  • d_min (float) – Optional resolution cutoff to be applied to the input intensities. If set to libtbx.Auto then d_min will be automatically determined according to the parameters min_i_mean_over_sigma_mean and min_cc_half.

  • min_i_mean_over_sigma_mean (float) – minimum value of \(|I|/|sigma(i)|\) for automatic determination of resolution cutoff.

  • min_cc_half (float) – minimum value of CC1/2 for automatic determination of resolution cutoff.

  • relative_length_tolerance (float) – Relative length tolerance in checking consistency of input unit cells against the median unit cell.

  • absolute_angle_tolerance (float) – Absolute angle tolerance in checking consistency of input unit cells against the median unit cell.

  • best_monoclinic_beta (bool) – If True, then for monoclinic centered cells, I2 will be preferred over C2 if it gives a less oblique cell (i.e. smaller beta angle).

as_dict()[source]

Return a dictionary representation of the results.

Returns:

dict

as_json(filename=None, indent=2)[source]

Return a json representation of the results.

Parameters:
  • filename (str) – Optional filename to export the json representation of the results.

  • indent (int) – The indent level for pretty-printing of the json. If None is the most compact representation.

Return type:

str

class dials.algorithms.symmetry.laue_group.ScoreCorrelationCoefficient(cc, sigma_cc, expected_cc, lower_bound=-1, upper_bound=1, k=2)[source]

Bases: object

__init__(cc, sigma_cc, expected_cc, lower_bound=-1, upper_bound=1, k=2)[source]
property p_cc_given_not_s

Probability of observing this CC if the sym op is NOT present, p(CC; !S).

property p_cc_given_s

Probability of observing this CC if the sym op is present, p(CC; S).

Modelled by a Cauchy distribution centred on cc_true and width gamma = sigma_cc

property p_s_given_cc

The likelihood of this symmetry element being present.

p(S; CC) = p(CC; S) / (p(CC; S) + p(CC; !S))

class dials.algorithms.symmetry.laue_group.ScoreSubGroup(subgroup, sym_op_scores)[source]

Bases: object

Score the probability of a given subgroup being the true subgroup.

  1. Calculates the combined correlation coefficients for symmetry operations present/absent from the subgroup.

  2. Calculates overall Zcc scores for symmetry elements present/absent from the subgroup.

  3. Calculates the overall likelihood for this subgroup.

See appendix A2 of Evans, P. R. (2011). Acta Cryst. D67, 282-292.

__init__(subgroup, sym_op_scores)[source]

Initialise a ScoreSubGroup object.

Parameters:
as_dict()[source]

Return a dictionary representation of the subgroup scoring.

The dictionary will contain the following keys:
  • patterson_group: The current subgroup

  • likelihood: The likelihood of the subgroup being correct

  • confidence: The confidence of the subgroup being correct

  • z_cc_for: The combined Z-scores for all symmetry elements present in the subgroup

  • z_cc_against: The combined Z-scores for all symmetry elements present in the lattice group but not in the subgroup

  • z_cc_net: The net Z-score, i.e. z_cc_for - z_cc_against

  • cc_for: The overall correlation coefficient for all symmetry elements present in the subgroup

  • cc_against: The overall correlation coefficient for all symmetry elements present in the lattice group but not in the subgroup

  • max_angular_difference: The maximum angular difference between the symmetrised unit cell and the P1 unit cell.

  • cb_op: The change of basis operation from the input unit cell to the ‘best’ unit cell.

Return type:

dict

class dials.algorithms.symmetry.laue_group.ScoreSymmetryElement(intensities, sym_op, cc_true, cc_sig_fac)[source]

Bases: object

Analyse intensities for presence of a given symmetry operation.

  1. Calculate the correlation coefficient, CC, for the given sym op.

  2. Calculate the probability of observing this CC if the sym op is present, p(CC; S), modelled by a Cauchy distribution centred on cc_true and width gamma = sigma_cc.

  3. Calculate the probability of observing this CC if the sym op is NOT present, p(CC; !S).

  4. Calculate the likelihood of symmetry element being present, p(S; CC) = p(CC; S) / (p(CC; S) + p(CC; !S))

See appendix A1 of Evans, P. R. (2011). Acta Cryst. D67, 282-292.

__init__(intensities, sym_op, cc_true, cc_sig_fac)[source]

Initialise a ScoreSymmetryElement object.

Parameters:
  • intensities (cctbx.miller.array) – The intensities on which to perform symmetry analysis.

  • sym_op (cctbx.sgtbx.rt_mx) – The symmetry operation for analysis.

  • cc_true (float) – the expected value of CC if the symmetry element is present, E(CC; S)

  • cc_sig_fac (float) – Estimation of sigma(CC) as a function of sample size.

as_dict()[source]

Return a dictionary representation of the symmetry element scoring.

The dictionary will contain the following keys:
  • likelihood: The likelihood of the symmetry element being present

  • z_cc: The Z-score for the correlation coefficient

  • cc: The correlation coefficient for the symmetry element

  • n_ref: The number of reflections contributing to the correlation coefficient

  • operator: The xyz representation of the symmetry element

Return type:

dict

dials.algorithms.symmetry.laue_group.trunccauchy_pdf(x, a, b, loc=0, scale=1)[source]

Calculate a truncated Cauchy probability density function.

Parameters:
  • x (float) – The point at which to calculate the PDF.

  • a (float) – The lower bound of the truncated distribution.

  • b (float) – The upper bound of the truncated distribution.

  • loc (float) – The location parameter for the Cauchy distribution.

  • scale (float) – The scale parameter for the Cauchy distribution.

Returns:

The value of the probability density function.

Return type:

float

Functions to help with reindexing against a reference dataset.

dials.algorithms.symmetry.reindex_to_reference.determine_reindex_operator_against_reference(test_miller_set, reference_miller_set)[source]

Reindex a miller set to match a reference miller set.

This function takes two miller arrays, a reference and a test array. The space group is checked to see if any reindexing may be required to give consistent indexing between both datasets. If possible twin operators exist, the different indexing options are tested against the reference set, using the correlation between datasets as the test.

Parameters:
Returns:

The change of basis operator which should be applied to the test dataset to give consistent indexing with the reference.

Return type:

cctbx.sgtbx.change_of_basis_op

dials.algorithms.symmetry.cosym

Methods for symmetry determination from partial datasets.

This module implements the methods of Gildea, R. J. & Winter, G. (2018). Acta Cryst. D74, 405-410 for determination of Patterson group symmetry from sparse multi-crystal data sets in the presence of an indexing ambiguity.

class dials.algorithms.symmetry.cosym.CosymAnalysis(intensities, params, seed_dataset: int | None = None)[source]

Bases: symmetry_base, Subject

Perform cosym analysis.

Perform cosym analysis on the input intensities using the methods of Gildea, R. J. & Winter, G. (2018). Acta Cryst. D74, 405-410 for determination of Patterson group symmetry from sparse multi-crystal data sets in the presence of an indexing ambiguity.

__init__(intensities, params, seed_dataset: int | None = None)[source]

Initialise a CosymAnalysis object.

Parameters:
  • intensities (cctbx.miller.array) – The intensities on which to perform cosym analysis.

  • params (libtbx.phil.scope_extract) – Parameters for the analysis.

  • seed_dataset (int) – The index into the intensities list to use when choosing a seed dataset for the reindexing analysis (the x,y,z reindex mode will be used for this dataset). If None, a high density cluster point is chosen.

as_dict()[source]

Return a dictionary representation of the results.

Returns:

dict

as_json(filename=None, indent=2)[source]

Return a json representation of the results.

Parameters:
  • filename (str) – Optional filename to export the json representation of the results.

  • indent (int) – The indent level for pretty-printing of the json. If None is the most compact representation.

Return type:

str

run()[source]
class dials.algorithms.symmetry.cosym.ScoreSubGroup(subgroup, sym_op_scores)[source]

Bases: object

Score the probability of a given subgroup being the true subgroup.

  1. Calculates overall Zcc scores for symmetry elements present/absent from the subgroup.

  2. Calculates the overall likelihood for this subgroup.

See appendix A2 of Evans, P. R. (2011). Acta Cryst. D67, 282-292.

__init__(subgroup, sym_op_scores)[source]

Initialise a ScoreSubGroup object.

Parameters:
as_dict()[source]

Return a dictionary representation of the subgroup scoring.

The dictionary will contain the following keys:
  • patterson_group: The current subgroup

  • likelihood: The likelihood of the subgroup being correct

  • confidence: The confidence of the subgroup being correct

  • z_cc_for: The combined Z-scores for all symmetry elements present in the subgroup

  • z_cc_against: The combined Z-scores for all symmetry elements present in the lattice group but not in the subgroup

  • z_cc_net: The net Z-score, i.e. z_cc_for - z_cc_against

  • max_angular_difference: The maximum angular difference between the symmetrised unit cell and the P1 unit cell.

  • cb_op: The change of basis operation from the input unit cell to the ‘best’ unit cell.

Return type:

dict

property stars
class dials.algorithms.symmetry.cosym.ScoreSymmetryElement(cc, sigma_cc, cc_true)[source]

Bases: object

Analyse intensities for presence of a given symmetry operation.

  1. Calculate the probability of observing this CC if the sym op is present, p(CC; S), modelled by a Cauchy distribution centred on cc_true and width gamma = sigma_cc.

  2. Calculate the probability of observing this CC if the sym op is NOT present, p(CC; !S).

  3. Calculate the likelihood of symmetry element being present, p(S; CC) = p(CC; S) / (p(CC; S) + p(CC; !S))

See appendix A1 of Evans, P. R. (2011). Acta Cryst. D67, 282-292.

__init__(cc, sigma_cc, cc_true)[source]

Initialise a ScoreSymmetryElement object.

Parameters:
  • cc (float) – the correlation coefficient for this symmetry element

  • sigma_cc (float) – the estimated error in the correlation coefficient

  • cc_true (float) – the expected value of CC if the symmetry element is present, E(CC; S)

as_dict()[source]

Return a dictionary representation of the symmetry element scoring.

The dictionary will contain the following keys:
  • likelihood: The likelihood of the symmetry element being present

  • z_cc: The Z-score for the correlation coefficient

  • cc: The correlation coefficient for the symmetry element

  • operator: The xyz representation of the symmetry element

Return type:

dict

property stars
class dials.algorithms.symmetry.cosym.SymmetryAnalysis(coords, sym_ops, subgroups, cb_op_inp_min)[source]

Bases: object

__init__(coords, sym_ops, subgroups, cb_op_inp_min)[source]
as_dict()[source]

Return a dictionary representation of the results.

Returns:

dict

static subgroups_table(d)[source]
static summary_table(d)[source]
static sym_ops_table(d)[source]
dials.algorithms.symmetry.cosym.change_of_basis_op_to_best_cell(experiments, max_delta, relative_length_tolerance, absolute_angle_tolerance, best_subgroup)[source]

Compute change of basis op to map experiments from P1 cell to the best cell that matches the best subgroup

dials.algorithms.symmetry.cosym.extract_reference_intensities(params: iotbx.phil.scope_extract, wavelength: float) miller.array[source]

Target function for cosym analysis.

class dials.algorithms.symmetry.cosym.target.Target(intensities, lattice_ids, weights=None, min_pairs=3, lattice_group=None, dimensions=None)[source]

Bases: object

Target function for cosym analysis.

dim

The number of dimensions used in the analysis.

Type:

int

__init__(intensities, lattice_ids, weights=None, min_pairs=3, lattice_group=None, dimensions=None)[source]

Initialise a Target object.

Parameters:
  • intensities (cctbx.miller.array) – The intensities on which to perform cosym analysis.

  • lattice_ids (np.ndarray) – An array of equal size to intensities which maps each reflection to a given lattice (dataset).

  • weights (str) – Optionally include weights in the target function. Allowed values are None, “count” and “standard_error”. The default is to use no weights. If “count” is set, then weights are equal to the number of pairs of reflections used in calculating each value of the rij matrix. If “standard_error” is used, then weights are defined as \(w_{ij} = 1/s\), where \(s = \sqrt{(1-r_{ij}^2)/(n-2)}\). See also http://www.sjsu.edu/faculty/gerstman/StatPrimer/correlation.pdf.

  • min_pairs (int) – Only calculate the correlation coefficient between two datasets if they have more than min_pairs of common reflections.

  • lattice_group (cctbx.sgtbx.space_group) – Optionally set the lattice group to be used in the analysis.

  • dimensions (int) – Optionally override the number of dimensions to be used in the analysis. If not set, then the number of dimensions used is equal to the greater of 2 or the number of symmetry operations in the lattice group.

compute_functional(x: ndarray) float[source]

Compute the target function at coordinates x.

Parameters:

x (np.ndarray) – a flattened list of the N-dimensional vectors, i.e. coordinates in the first dimension are stored first, followed by the coordinates in the second dimension, etc.

Returns:

The value of the target function at coordinates x.

Return type:

f (float)

compute_gradients(x: ndarray) ndarray[source]

Compute the gradients of the target function at coordinates x.

Parameters:

x (np.ndarray) – a flattened list of the N-dimensional vectors, i.e. coordinates in the first dimension are stored first, followed by the coordinates in the second dimension, etc.

Returns:

f: The value of the target function at coordinates x. grad: The gradients of the target function with respect to the parameters.

Return type:

Tuple[float, np.ndarray]

compute_gradients_fd(x: ndarray, eps=1e-06) ndarray[source]

Compute the gradients at coordinates x using finite differences.

Parameters:
  • x (np.ndarray) – a flattened list of the N-dimensional vectors, i.e. coordinates in the first dimension are stored first, followed by the coordinates in the second dimension, etc.

  • eps (float) – The value of epsilon to use in finite difference calculations.

Returns:

The gradients of the target function with respect to the parameters.

Return type:

grad (np.ndarray)

curvatures(x: ndarray) ndarray[source]

Compute the curvature of the target function at coordinates x.

Parameters:

x (np.ndarray) – a flattened list of the N-dimensional vectors, i.e. coordinates in the first dimension are stored first, followed by the coordinates in the second dimension, etc.

Returns:

The curvature of the target function with respect to the parameters.

Return type:

curvs (np.ndarray)

curvatures_fd(x: ndarray, eps=1e-06) ndarray[source]

Compute the curvatures at coordinates x using finite differences.

Parameters:
  • x (np.ndarray) – a flattened list of the N-dimensional vectors, i.e. coordinates in the first dimension are stored first, followed by the coordinates in the second dimension, etc.

  • eps (float) – The value of epsilon to use in finite difference calculations.

Returns:

The curvature of the target function with respect to the parameters.

Return type:

curvs (np.ndarray)

set_dimensions(dimensions)[source]

Set the number of dimensions for analysis.

Parameters:

dimensions (int) – The number of dimensions to be used.

LBFGS refinement engine for cosym analysis.

class dials.algorithms.symmetry.cosym.engine.lbfgs_with_curvs(target, coords, use_curvatures=True, termination_params=None)[source]

Bases: object

Minimise a target function using the LBFGS minimiser.

Implementation of an LBFGS minimiser using curvature information, according to the interface defined by scitbx.lbfgs.

__init__(target, coords, use_curvatures=True, termination_params=None)[source]

Initialise an lbfgs_with_curvs object.

Parameters:
  • target (dials.algorithms.target.Target) – The target function to minimise.

  • coords (np.ndarray) – The starting coordinates for minimisation.

  • use_curvatures (bool) – Whether or not to use curvature information in the minimisation. Defaults to True.

  • termination_params (scitbx.lbfgs.termination_parameters) – Override the default termination parameters for the minimisation.

callback_after_step(minimizer)[source]

Log progress after each successful step of the minimisation.

compute_functional_and_gradients()[source]

Compute the functional and gradients.

Returns:

A tuple of the functional and gradients.

Return type:

tuple

compute_functional_gradients_and_curvatures()[source]

Compute the functional, gradients and curvatures.

Returns:

A tuple of the functional, gradients and curvatures.

Return type:

tuple

compute_functional_gradients_diag()[source]

Compute the functional, gradients and diagonal.

Returns:

A tuple of the functional, gradients and diagonal, where the diagonal is the reciprocal of the curvatures.

Return type:

tuple

dials.algorithms.symmetry.cosym.engine.minimize_scipy(target, coords, method='L-BFGS-B', max_iterations=None, max_calls=None)[source]

Thin wrapper around scipy.optimize.minimize.

Parameters:
  • target (dials.algorithms.target.Target) – The target function to minimise.

  • coords (np.array) – The starting coordinates for minimisation.

dials.algorithms.symmetry.cosym.engine.minimize_scitbx_lbfgs(target, coords, use_curvatures=True, max_iterations=100, max_calls=None)[source]