pyextal.jaxMultipole.fit¶

Weighted least-squares fit of a multipole model to structure factors.

The objective is the usual

\[\chi^2 = \sum_i \left(\frac{y_i^{\rm obs} - y_i^{\rm calc}(\mathbf x)}{\sigma_i}\right)^2 + \sum_r \left(\frac{g_r(\mathbf p) - t_r}{\sigma_r}\right)^2\]

over the observations (real parts, and imaginary parts where an observation carries one) and any soft restraints. It is minimised with scipy.optimize.least_squares (trust-region reflective, box bounds), with an exact Jacobian from jax.jacfwd – in float64, inside jax.enable_x64, so the rest of the process (the float32 Bloch engine) is unaffected.

What is reported¶

esd – from (J^T J)^-1, i.e. taking the sigmas at face value, and esd_scaled – the same times the goodness of fit, the conventional crystallographic choice when GoF > 1;
the correlation matrix of the free parameters and the pairs above 0.9;
R = sum |dy| / sum |y|, wR = sqrt(sum w dy^2 / sum w y^2), GoF = sqrt(chi^2 / (n - m));
the singular values and condition number of the weighted Jacobian.

design_check() gives the last two before fitting, at the start point: which directions in parameter space the data do not constrain, and which parameters they load on. With ~25 reflections and ADPs refined it is the first thing to look at. hamilton_test() compares nested models.

Device: the fit runs on the CPU by default (device='cpu'). With some 30 reflections it is a few hundred small float64 kernels, which a CPU does in ~15 s for the STO ladder while a consumer GPU – float64 at a small fraction of its float32 rate, plus launch overhead – takes ~50 s, to identical results. device='default' uses whatever JAX picks; the GPU matters for refining against intensities (Option 2), not here.

Classes¶

`Restraint`	A soft restraint on one physical parameter: `(p[name] - target) / sigma`.
`FitResult`	Outcome of `fit()`.
`Fitter`	A model and a data set with their residual and Jacobian compiled once.

Functions¶

`residual_fn`(model, obs, pmap[, restraints])	`x -> weighted residual vector` (real), plus `x -> calc`.
`fit`(model, obs, pmap[, restraints, x0, max_nfev, ...])	Fits the free parameters of `pmap` to `obs`.
`design_check`(model, obs, pmap[, restraints, x, top, ...])	How well the data constrain each direction in parameter space.
`hamilton_test`(fit_restricted, fit_full[, alpha])	Hamilton's R-factor ratio test of a nested pair of fits.

Module Contents¶

class pyextal.jaxMultipole.fit.Restraint¶

A soft restraint on one physical parameter: (p[name] - target) / sigma.

name: str¶

target: float¶

sigma: float¶

class pyextal.jaxMultipole.fit.FitResult¶

Outcome of fit().

Attributes are NumPy arrays or floats; names label x, esd and the rows/columns of cov/corr; p and p_esd follow the model layout.

names: list¶

x: numpy.ndarray¶

p: numpy.ndarray¶

p_esd: numpy.ndarray¶

cov: numpy.ndarray¶

esd: numpy.ndarray¶

esd_scaled: numpy.ndarray¶

corr: numpy.ndarray¶

chi2: float¶

n_obs: int¶

n_par: int¶

gof: float¶

R: float¶

wR: float¶

calc: numpy.ndarray¶

obs: object¶

singular_values: numpy.ndarray = None¶

condition: float = None¶

success: bool = True¶

message: str = ''¶

pmap: object = None¶

value(name)¶

p_esd_of(name)¶

has(name)¶

strong_correlations(threshold=0.9)¶

pyextal.jaxMultipole.fit.residual_fn(model, obs, pmap, restraints=())¶

x -> weighted residual vector (real), plus x -> calc.

Real parts always count; imaginary parts only where obs.use_imag.

class pyextal.jaxMultipole.fit.Fitter(model, obs, restraints=(), device='cpu')¶

A model and a data set with their residual and Jacobian compiled once.

The residuals are compiled as a function of the full physical vector p, so any ParameterMap – any choice of frozen and free parameters – reuses them: r(x) = r_p(p0 + M x), J(x) = J_p(p) M. Refinements that change what is refined many times over (ladders, iterative block schemes, chi^2 profiles) then cost one compilation instead of one per fit. fit keeps a small cache of these, so ordinary calls benefit too.

restraints = ()¶

fit(pmap, x0=None, max_nfev=200, verbose=0)¶: Fits the free parameters of pmap; see fit().

pyextal.jaxMultipole.fit.fit(model, obs, pmap, restraints=(), x0=None, max_nfev=200, verbose=0, device='cpu')¶

Fits the free parameters of pmap to obs.

Parameters:

model (MultipoleModel) – The model.
obs (Observations) – Data (jaxMultipole.data).
pmap (ParameterMap) – Free parameters (jaxMultipole.params).
restraints (list[Restraint]) – Soft restraints.
x0 (array, optional) – Start; pmap.x0 by default.
max_nfev (int) – Function-evaluation budget.
verbose (int) – least_squares verbosity.
device (str) – 'cpu' (default), 'gpu' or 'default'.

Returns:

FitResult

Note

The compiled residuals are cached per (model, obs, restraints, device) object, so repeated fits of the same data with different free parameters compile once (Fitter).

pyextal.jaxMultipole.fit.design_check(model, obs, pmap, restraints=(), x=None, top=4, device='cpu')¶

How well the data constrain each direction in parameter space.

Returns:: singular_values, condition, predicted esd and corr at x (the start by default), and weak – for each singular value below 1e-3 of the largest, the parameters it loads on.
Return type:: dict

pyextal.jaxMultipole.fit.hamilton_test(fit_restricted, fit_full, alpha=0.05)¶

Hamilton’s R-factor ratio test of a nested pair of fits.

Parameters:

fit_restricted (FitResult) – The model with fewer parameters.
fit_full (FitResult) – The model with b more.
alpha (float) – Significance level.

Returns:

ratio = wR_restricted / wR_full, the critical value sqrt(1 + b/(n-m) F_{b, n-m, alpha}) and significant.

Return type:

dict