clophfit.prenspire.mcr ====================== .. py:module:: clophfit.prenspire.mcr .. autoapi-nested-parse:: MCR-ALS: multivariate curve resolution by alternating least squares. Factor the corrected spectrum matrix of one sample/temperature, ``D[n_wells, n_lambda]``, into non-negative concentrations and spectra:: D ~= C @ S.T C >= 0, S >= 0, rows of C closed (sum to 1) Alternating non-negative least squares (scipy ``nnls``) solves one factor with the other fixed; closure (``C`` rows normalised to 1) pins the rotational freedom that ordinary PCA leaves open, which is what makes the recovered species spectra interpretable without a mechanistic model. The pKa is read off afterwards: the component whose closed profile tracks the titrant (largest |Spearman rho|) is fitted against a Henderson-Hasselbalch fraction, ``profile ~= F0 + (F1 - F0)/(1 + 10**(K - x))``, leaving K (and a residual affine mixing) free. Contrast with :func:`clophfit.fitting.spectral.fit_spectra_global`: there the two-state model *is* the factorisation (one K, species spectra solved linearly). MCR-ALS relaxes the profile shape - it can expose a third species or a non-Hill transition - at the price of a less tight K. ``fit_spectra_global``'s ``sv_ratio`` says when that trade is worth it. Classes ------- .. autoapisummary:: clophfit.prenspire.mcr.McrFit clophfit.prenspire.mcr.StackedScans Functions --------- .. autoapisummary:: clophfit.prenspire.mcr.dead_components clophfit.prenspire.mcr.mcr_als clophfit.prenspire.mcr.model_fraction clophfit.prenspire.mcr.component_pka clophfit.prenspire.mcr.stack_scans clophfit.prenspire.mcr.scan_blocks clophfit.prenspire.mcr.species_by_scan clophfit.prenspire.mcr.mcr_titration Module Contents --------------- .. py:class:: McrFit Bases: :py:obj:`NamedTuple` A converged MCR-ALS factorisation and its reconstruction quality. .. py:attribute:: c :type: numpy.ndarray [n_wells, n_species] closed, non-negative concentration profiles. .. py:attribute:: s :type: numpy.ndarray [n_species, n_lambda] non-negative species spectra (arbitrary scale). .. py:attribute:: rel_err :type: float ||D - C S.T|| / ||D||. .. py:attribute:: dead_components :type: tuple[int, ...] :value: () Components whose profile never varies across the wells. A dead component is the signature of the MCR scale degeneracy: the species was parked at an arbitrary scale while its concentration column sat at a constant, so it holds no titration information. .. py:property:: degenerate :type: bool Whether any component is dead, i.e. the decomposition is unusable. ``rel_err`` cannot detect this: a spurious component is a free lamp that lowers the residual while contributing nothing, so a degenerate fit can report a *better* error than a sound one. .. py:function:: dead_components(c, *, absolute = _DEAD_PTP_ABS, ratio = _DEAD_PTP_RATIO) Find the concentration columns that carry no well-to-well change. A component whose profile is the same number in every well cannot be fitted for a K: the curve would be fitted to round-off, and the result is a confident, meaningless number. A column counts as dead when it varies by less than the absolute floor, or by less than a small fraction of the most mobile column. The absolute floor catches the case where every column is dead; the relative term catches a redundant column parked on a scale that is not itself tiny. :param c: Closed concentration profiles ``[n_wells, n_species]``. :type c: np.ndarray :param absolute: Peak-to-peak below which a column counts as dead outright. :type absolute: float :param ratio: A column whose peak-to-peak is below this fraction of the most mobile column's also counts as dead. :type ratio: float :returns: The indices of the dead columns, in order. :rtype: tuple[int, ...] .. py:function:: mcr_als(d, n_species = 2, *, closure = True, max_iter = _MAX_ITER, tol = _TOL, s_init = None) Resolve D [n_wells, n_lambda] into C [n_wells, k] and S [n_lambda, k]. :param d: Data matrix, one row per well, one column per wavelength (already corrected: buffer-subtracted). :type d: np.ndarray :param n_species: Number of components to resolve. :type n_species: int :param closure: Normalise each row of C to sum to 1 (fractions), pinning rotation. :type closure: bool :param max_iter: ALS sweeps per start. :type max_iter: int :param tol: Relative error change that counts as converged. :type tol: float :param s_init: Optional extra spectra start [n_lambda, k] - e.g. the two extreme wells of a titration. ALS from a single start is a non-convex problem; all starts (this one, clipped and absolute PCA loadings) are run and the lowest-error factorisation returned. :type s_init: np.ndarray | None :returns: Non-negative factors and the relative reconstruction error. :rtype: McrFit :raises ValueError: If the requested number of species exceeds the matrix rank bound. .. py:function:: model_fraction(x, k, *, is_ph = True) Return the model's bound fraction on the 0-1 scale. This is the quantity the figure draws over the closed profiles, and the one the component is compared against. Keeping it a fraction on one scale means the curve and the points it overlays share an axis, so the distance between them is a real measure of fit quality. :param x: Titrant values. :type x: np.ndarray :param k: pKa for a pH titration, Kd in the units of ``x`` for a ligand one. :type k: float :param is_ph: Henderson-Hasselbalch fraction (pH) if true, Langmuir (ligand) else. :type is_ph: bool :returns: The fraction in [0, 1] for every value of ``x``. :rtype: np.ndarray .. py:function:: component_pka(x, c, *, is_ph = True) Equilibrium constant of the component that tracks the titrant. Returns K on the natural scale of the model fitted to the component: pKa for a pH titration (Henderson-Hasselbalch fraction), Kd in the units of x for a ligand titration (Langmuir fraction - linear, so no log-space propagation happens here). The fitted curve parameters (F0, F1) are the tracking component's own plateau pair, i.e. the affine mixing that maps the model fraction onto that component's scale; they are properties of the component, not of the model, and they are free because the component's baseline need not sit at 0 and its saturation need not sit at 1. :param x: Titrant per well. :type x: np.ndarray :param c: Closed concentration profiles [n_wells, k]. :type c: np.ndarray :param is_ph: Henderson-Hasselbalch fraction (pH) if true, Langmuir (ligand) else. :type is_ph: bool :returns: (K, standard error, F0, F1, component index, signed rho of that component). With k = 2 a positively tracking component is preferred: under closure the columns are complements, and the signed choice keeps the reported component index stable instead of arbitrarily picking whichever ties first on |rho|. :rtype: tuple[float, float, float, float, int, float] .. py:class:: StackedScans Bases: :py:obj:`NamedTuple` The scans of one plate, joined into one data matrix. The columns keep one block per scan: the scans are different physical measurements, so they do not share a wavelength axis. .. py:attribute:: tags :type: numpy.ndarray Scan name of every column. .. py:attribute:: wavelengths :type: numpy.ndarray Wavelength of every column, in its own scan's units. .. py:attribute:: data :type: numpy.ndarray D[n_wells, n_columns], each scan scaled by the common maximum. .. py:function:: stack_scans(scans) Join per-scan (wavelengths, Y[n_lambda, n_wells]) blocks into one matrix. The scans of one plate are separate physical measurements - an excitation scan and an emission scan are different observables - and they do not share one wavelength axis: their ranges overlap but the values at an overlapping wavelength are not the same measurement. Stacking them therefore keeps one column block per scan, tagged with the scan it came from, rather than gluing them into a single ramp that runs backwards across the seam. :param scans: One (scan name, (wavelengths, Y[n_lambda, n_wells])) block per scan. :type scans: Sequence[tuple[str, tuple[ArrayF, ArrayF]]] :returns: The scan tag and the wavelength of every column, and the data matrix itself, each scan scaled by the common maximum so blocks stay commensurable. :rtype: StackedScans .. py:function:: scan_blocks(stacked) Map each scan name to the column block it occupies. :param stacked: The stacked scans. :type stacked: StackedScans :returns: Scan name to the slice of columns it holds, in order of appearance. :rtype: dict[str, slice] .. py:function:: species_by_scan(spec, n_species) Split the resolved species spectra into one curve set per scan. The species spectra of a scan are that scan's own measurement: an excitation scan's species spectrum is an excitation spectrum, and the emission scan's is an emission spectrum. Drawing them on one axis would put two different observables on the same ramp. :param spec: The resolved spectra, with ``scan``, ``lambda`` and one column per species. :type spec: pd.DataFrame :param n_species: How many species columns to read. :type n_species: int :returns: Scan name to (wavelengths, one array per species). :rtype: dict[str, tuple[np.ndarray, list[np.ndarray]]] .. py:function:: mcr_titration(spectra_by_temp, x_by_temp, *, n_species = 2, closure = True, is_ph = True) MCR-ALS every temperature's scans; return (K table, spectra, profiles). :param spectra_by_temp: Temperature label to its scans, each ``(name, (lambda, Y[n_lambda, n_wells]))``. :type spectra_by_temp: dict[str, Sequence[tuple[str, tuple[ArrayF, ArrayF]]]] :param x_by_temp: Titrant value per well, same temperature keys. :type x_by_temp: dict[str, ArrayF] :param n_species: Components resolved per temperature. :type n_species: int :param closure: Constrain profiles to fractions. :type closure: bool :param is_ph: Titrant axis kind: pH (Henderson-Hasselbalch) or ligand (Langmuir). :type is_ph: bool :returns: One table row per (temp): K, se, component, rho, rel_err, n; per temperature a table of the resolved species spectra, tagged with the scan each row came from; and per temperature the titrant with the closed profiles, for plotting. :rtype: tuple[pd.DataFrame, dict[str, pd.DataFrame], dict[str, tuple[ArrayF, np.ndarray]]]