hotelling_parameters

hotelling_parameters(
    x: Union[np.ndarray, pd.DataFrame],
    k: int = 2,
    pcx: int = 1,
    pcy: int = 2,
    threshold: Optional[float] = None,
    rel_tol: float = 0.001,
    abs_tol: float = sys.float_info.epsilon,
    method: Literal['f', 'beta'] = 'f',
    conf_limit: Union[float, Sequence[float]] = (0.95, 0.99),
)

Hotelling’s T-squared statistic and Hotelling’s ellipse parameters.

Computes Hotelling’s T-squared statistic of each observation, its cutoffs at one or more confidence levels and, when two components are used, the semi-axes and rotation of the corresponding Hotelling’s ellipse. The number of components is either fixed (k) or chosen from a cumulative explained variance threshold.

Parameters

x : numpy.ndarray or pandas.DataFrame

Scores from PCA, PLS, ICA, or similar methods, with one row per observation and one column per component.

k : int = 2

Number of components to use. Ignored when threshold is given.

pcx : int = 1

Component (1-based) on the x-axis of the ellipse when k=2.

pcy : int = 2

Component (1-based) on the y-axis of the ellipse when k=2. Must differ from pcx.

threshold : float = None

Cumulative explained variance threshold, in (0, 1]. When given, the smallest number of leading components that explain at least this proportion of the total variance is used (at least two).

rel_tol : float = 0.001

Minimum proportion of the total variance a component must explain to be kept.

abs_tol : float = sys.float_info.epsilon

Minimum variance a component must have to be kept. Must not exceed rel_tol.

method : str = 'f'

Distribution of the T-squared cutoffs, 'f' or 'beta'. See Notes.

conf_limit : float or sequence of float = (0.95, 0.99)

Confidence level, or levels, of the cutoffs and semi-axes, each strictly between 0 and 1.

Returns

: dict

T-squared values, cutoffs, number of components and, when k=2, ellipse parameters, under the keys listed in Returned keys.

Returned keys

  • 'Tsquared': DataFrame with a column 'value', the T-squared statistic of each observation (its squared Mahalanobis distance), on the same scale as the cutoffs. When k=2, it is computed on components pcx and pcy.
  • 'cutoff_<level>pct': T-squared cutoff at each level of conf_limit, from the highest to the lowest level: 'cutoff_99pct' and 'cutoff_95pct' by default, 'cutoff_99.9pct' for a level of 0.999.
  • 'nb_comp': number of components used.
  • 'Ellipse' (only when k=2): one-row DataFrame with the semi-axes 'a_<level>pct' and 'b_<level>pct' at each level, and the rotation 'angle' of the ellipse in radians. a is the semi-axis closest to the pcx direction. 'angle' is 0 for uncorrelated scores, such as the PCA or PLS scores of the samples a model was fitted on.

Raises

: TypeError

If x is not a NumPy array or a pandas DataFrame.

: ValueError

If an argument is invalid, if there are fewer than k + 2 observations, or if fewer than two components remain after removing near-zero variance components.

: RuntimeError

If the covariance matrix of the selected components is singular.

Warnings

A UserWarning is issued when near-zero variance components are removed, or when threshold is lower than the variance explained by the first component (two components are then used).

Notes

For observation \(i\) with scores \(\mathbf{x}_i\) on the \(k\) selected components, \(T^2_i = (\mathbf{x}_i - \bar{\mathbf{x}})^\top \mathbf{S}^{-1} (\mathbf{x}_i - \bar{\mathbf{x}})\), where \(\bar{\mathbf{x}}\) and \(\mathbf{S}\) are the sample mean and covariance matrix of the \(n\) observations.

With method='f' (default), the cutoff at confidence level \(1 - \alpha\) is \(\frac{k(n-1)}{n-k} F_{1-\alpha}(k, n-k)\), as in previous versions of the package. With method='beta', it is \(\frac{(n-1)^2}{n} B_{1-\alpha}(k/2, (n-k-1)/2)\), the exact distribution of \(T^2_i\) for the observations used to estimate the mean and covariance (Tracy, Young and Mason, 1992), e.g. the scores of the samples a PCA or PLS model was built on. For these observations, the F-based cutoff is conservative, especially for small \(n\).

The ellipse is the set of points whose T-squared value equals the cutoff, so an observation lies outside it exactly when its T-squared value exceeds the cutoff. Derivations are given on the Theory page of the documentation, https://christiangoueguel.com/pyEllipse/theory.html.

References

Tracy, N. D., Young, J. C. and Mason, R. L. (1992). Multivariate control charts for individual observations. Journal of Quality Technology, 24(2), 88-95.

Examples

>>> import numpy as np
>>> from pyEllipse import hotelling_parameters
>>> scores = np.random.default_rng(0).standard_normal((50, 3)) * [3.0, 2.0, 1.0]
>>> res = hotelling_parameters(scores, k=2)
>>> list(res)
['Tsquared', 'cutoff_99pct', 'cutoff_95pct', 'nb_comp', 'Ellipse']
>>> round(res['cutoff_95pct'], 3)
6.514
>>> outliers = res['Tsquared']['value'] > res['cutoff_95pct']

Exact cutoffs for the calibration samples, at custom confidence levels:

>>> res = hotelling_parameters(scores, k=3, method='beta', conf_limit=(0.975, 0.999))
>>> [key for key in res if key.startswith('cutoff')]
['cutoff_99.9pct', 'cutoff_97.5pct']